From ef152a15c4aed1f50a2340f951d11e48818e552c Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Mon, 14 Aug 2023 14:11:34 -0600 Subject: [PATCH 01/60] draft program to experiment with reading table values into corresponding types --- .../modules/assimilation/type_read_table.f90 | 129 ++++++++++++++++++ 1 file changed, 129 insertions(+) create mode 100644 assimilation_code/modules/assimilation/type_read_table.f90 diff --git a/assimilation_code/modules/assimilation/type_read_table.f90 b/assimilation_code/modules/assimilation/type_read_table.f90 new file mode 100644 index 0000000000..f164221db3 --- /dev/null +++ b/assimilation_code/modules/assimilation/type_read_table.f90 @@ -0,0 +1,129 @@ +program read_table + +implicit none +type obs_error_info_type + logical :: bounded_below, bounded_above + real :: lower_bound, upper_bound +end type + +type probit_inflation_type + integer :: dist_type + logical :: bounded_below, bounded_above + real :: lower_bound, upper_bound +end type + +type probit_state_type + integer :: dist_type + logical :: bounded_below, bounded_above + real :: lower_bound, upper_bound +end type + +type probit_extended_state_type + integer :: dist_type + logical :: bounded_below, bounded_above + real :: lower_bound, upper_bound +end type + +type obs_inc_info_type + integer :: filter_kind + logical :: rectangular_quadrature, gaussian_likelihood_tails + logical :: sort_obs_inc, spread_restoration + logical :: bounded_below, bounded_above + real :: lower_bound, upper_bound +end type + +type qcf_table_data_type + type(obs_error_info_type) :: obs_error_info + type(probit_inflation_type) :: probit_inflation + type(probit_state_type) :: probit_state + type(probit_extended_state_type) :: probit_extended_state + type(obs_inc_info_type) :: obs_inc_info +end type + +! Reads in the QCEFF input options from tabular data file +!character(len=50), intent(in) :: qcf_table_filename +!real(r8), intent(out) :: qcf_table_data +!real, dimension(:, :), allocatable :: qcf_table_data_rows +type(qcf_table_data_type), allocatable :: qcf_table_data(:) +character(len=30), dimension(:), allocatable :: rowheaders !!!!! might need to change len=30 + +integer, parameter :: fileid = 10 !file identifier +character(len=30), parameter :: tester_QTY = 'QTY_GPSRO' +integer :: QTY_loc(1) + +!integer, parameter :: num_columns = 28 +integer :: nlines +integer :: io +integer :: numrows +integer :: row + +!real, dimension(1:num_columns, 1:num_rows) :: table_data +!integer :: table_data_1, table_data_2 +character(len=30), dimension(4) :: header1 +character(len=30), dimension(29) :: header2 +!variables for table values ^^^ + +open(unit=fileid, file='cam_qcf_table.dat') +nlines = 0 + +do !do loop to get number of rows (or QTY's) in the table + read(fileid,*,iostat=io) + if(io/=0) exit + nlines = nlines + 1 +end do +close(fileid) + +print*, nlines + +numrows = nlines - 2 +print *, 'numrows: ', numrows + +allocate(qcf_table_data(numrows)) +allocate(rowheaders(numrows)) +write(*,*) shape(qcf_table_data) + +open(unit=fileid, file='cam_qcf_table.dat') + +read(fileid, *) header1 +read(fileid, *) header2 !! skip the headers +Write(*, *) "header1: ", header1 +Write(*, *) "header2: ", header2 + +do row = 1, numrows + read(fileid, *) rowheaders(row), qcf_table_data(row)%obs_error_info%bounded_below, qcf_table_data(row)%obs_error_info%bounded_above, & + qcf_table_data(row)%obs_error_info%lower_bound, qcf_table_data(row)%obs_error_info%upper_bound, qcf_table_data(row)%probit_inflation%dist_type, & + qcf_table_data(row)%probit_inflation%bounded_below, qcf_table_data(row)%probit_inflation%bounded_above, & + qcf_table_data(row)%probit_inflation%lower_bound, qcf_table_data(row)%probit_inflation%upper_bound, qcf_table_data(row)%probit_state%dist_type, & + qcf_table_data(row)%probit_state%bounded_below, qcf_table_data(row)%probit_state%bounded_above, & + qcf_table_data(row)%probit_state%lower_bound, qcf_table_data(row)%probit_state%upper_bound, qcf_table_data(row)%probit_extended_state%dist_type, & + qcf_table_data(row)%probit_extended_state%bounded_below, qcf_table_data(row)%probit_extended_state%bounded_above, & + qcf_table_data(row)%probit_extended_state%lower_bound, qcf_table_data(row)%probit_extended_state%upper_bound, & + qcf_table_data(row)%obs_inc_info%filter_kind, qcf_table_data(row)%obs_inc_info%rectangular_quadrature, & + qcf_table_data(row)%obs_inc_info%gaussian_likelihood_tails, qcf_table_data(row)%obs_inc_info%sort_obs_inc, & + qcf_table_data(row)%obs_inc_info%spread_restoration, qcf_table_data(row)%obs_inc_info%bounded_below, qcf_table_data(row)%obs_inc_info%bounded_above, & + qcf_table_data(row)%obs_inc_info%lower_bound, qcf_table_data(row)%obs_inc_info%upper_bound + + write(*, *) "rowheader(", row, "): ", rowheaders(row) + write(*, *) "qcf_table_data(", row, "): " + write(*, *) qcf_table_data(row)%obs_error_info%bounded_below, qcf_table_data(row)%obs_error_info%bounded_above, & + qcf_table_data(row)%obs_error_info%lower_bound, qcf_table_data(row)%obs_error_info%upper_bound, qcf_table_data(row)%probit_inflation%dist_type, & + qcf_table_data(row)%probit_inflation%bounded_below, qcf_table_data(row)%probit_inflation%bounded_above, & + qcf_table_data(row)%probit_inflation%lower_bound, qcf_table_data(row)%probit_inflation%upper_bound, qcf_table_data(row)%probit_state%dist_type, & + qcf_table_data(row)%probit_state%bounded_below, qcf_table_data(row)%probit_state%bounded_above, & + qcf_table_data(row)%probit_state%lower_bound, qcf_table_data(row)%probit_state%upper_bound, qcf_table_data(row)%probit_extended_state%dist_type, & + qcf_table_data(row)%probit_extended_state%bounded_below, qcf_table_data(row)%probit_extended_state%bounded_above, & + qcf_table_data(row)%probit_extended_state%lower_bound, qcf_table_data(row)%probit_extended_state%upper_bound, & + qcf_table_data(row)%obs_inc_info%filter_kind, qcf_table_data(row)%obs_inc_info%rectangular_quadrature, & + qcf_table_data(row)%obs_inc_info%gaussian_likelihood_tails, qcf_table_data(row)%obs_inc_info%sort_obs_inc, & + qcf_table_data(row)%obs_inc_info%spread_restoration, qcf_table_data(row)%obs_inc_info%bounded_below, qcf_table_data(row)%obs_inc_info%bounded_above, & + qcf_table_data(row)%obs_inc_info%lower_bound, qcf_table_data(row)%obs_inc_info%upper_bound +end do + +close(fileid) + +QTY_loc = findloc(rowheaders, tester_QTY) +write(*, *) 'findloc of GPSRO: ', QTY_loc(1) + +deallocate(qcf_table_data, rowheaders) + +end program read_table From ebdbf31baa53f8ecb94c305940b6e6638a89f0c0 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Mon, 14 Aug 2023 14:18:48 -0600 Subject: [PATCH 02/60] prototype table data file that uses CAM-FV QTYs --- .../modules/assimilation/cam_qcf_table.txt | 10 ++++++++++ 1 file changed, 10 insertions(+) create mode 100644 assimilation_code/modules/assimilation/cam_qcf_table.txt diff --git a/assimilation_code/modules/assimilation/cam_qcf_table.txt b/assimilation_code/modules/assimilation/cam_qcf_table.txt new file mode 100644 index 0000000000..56e205a534 --- /dev/null +++ b/assimilation_code/modules/assimilation/cam_qcf_table.txt @@ -0,0 +1,10 @@ +QCF table version 1: +QTY bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound filter_kind rectangular_quadrature gaussian_likelihood_tails sort_obs_inc spread_restoration bounded_below bounded_above lower_bound upper_bound +QTY_U_WIND_COMPONENT, .false. .false. 3 4 5 .false. .false. 8 9 10 .false. .false. 13 14 15 .false. .false. 18 19 20 .false. .false. .false. .false. .false. .false. 27 28 +QTY_V_WIND_COMPONENT, .false. .false. 3 4 5 .false. .false. 8 9 10 .false. .false. 13 14 15 .false. .false. 18 19 20 .false. .false. .false. .false. .false. .false. 27 28 +QTY_SURFACE_PRESSURE .false. .false. 3 4 5 .false. .false. 8 9 10 .false. .false. 13 14 15 .false. .false. 18 19 20 .false. .false. .false. .false. .false. .false. 27 28 +QTY_TEMPERATURE .false. .false. 3 4 5 .false. .false. 8 9 10 .false. .false. 13 14 15 .false. .false. 18 19 20 .false. .false. .false. .false. .false. .false. 27 28 +QTY_SPECIFIC_HUMIDITY .false. .false. 3 4 5 .false. .false. 8 9 10 .false. .false. 13 14 15 .false. .false. 18 19 20 .false. .false. .false. .false. .false. .false. 27 28 +QTY_CLOUD_LIQUID_WATER .false. .false. 3 4 5 .false. .false. 8 9 10 .false. .false. 13 14 15 .false. .false. 18 19 20 .false. .false. .false. .false. .false. .false. 27 28 +QTY_CLOUD_ICE .false. .false. 3 4 5 .false. .false. 8 9 10 .false. .false. 13 14 15 .false. .false. 18 19 20 .false. .false. .false. .false. .false. .false. 27 28 +QTY_GPSRO .false. .false. 3 4 5 .false. .false. 8 9 10 .false. .false. 13 14 15 .false. .false. 18 19 20 .false. .false. .false. .false. .false. .false. 27 28 From 262965b0b100f3eecc4b65296b856072be48bd3f Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Mon, 14 Aug 2023 14:44:50 -0600 Subject: [PATCH 03/60] adding new subroutine init_qcf_table to return number of rows in table --- .../assimilation/algorithm_info_mod.f90 | 28 +++++++++++++++++++ 1 file changed, 28 insertions(+) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index 7ea474c664..78c2f190ab 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -35,6 +35,7 @@ module algorithm_info_mod integer, parameter :: BOUNDED_NORMAL_RHF = 101 public :: obs_error_info, probit_dist_info, obs_inc_info, & + init_qcf_table, & EAKF, ENKF, BOUNDED_NORMAL_RHF, UNBOUNDED_RHF, GAMMA_FILTER ! Provides routines that give information about details of algorithms for @@ -238,4 +239,31 @@ end subroutine obs_inc_info !------------------------------------------------------------------------ + +subroutine init_qcf_table(qcf_table_filename, numrows) + +character(len=50), intent(in) :: qcf_table_filename +integer, intent(out) :: numrows !return value + +integer :: nlines +integer :: io +integer, parameter :: fileid = 10 !file identifier + +open(unit=fileid, file=qcf_table_filename) +nlines = 0 + +do !do loop to get number of rows (or QTY's) in the table + read(fileid,*,iostat=io) + if(io/=0) exit + nlines = nlines + 1 +end do +close(fileid) + +numrows = nlines - 2 +print *, 'numrows: ', numrows + +end subroutine init_qcf_table + +!------------------------------------------------------------------------ + end module algorithm_info_mod From 715875db83a1b66c523eb75f7d106999e8ef4940 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Mon, 14 Aug 2023 15:25:08 -0600 Subject: [PATCH 04/60] Adding a new namelist variable to the assim_tools_nml --- .../modules/assimilation/assim_tools_mod.f90 | 11 ++++++++--- 1 file changed, 8 insertions(+), 3 deletions(-) diff --git a/assimilation_code/modules/assimilation/assim_tools_mod.f90 b/assimilation_code/modules/assimilation/assim_tools_mod.f90 index 72a01dc383..c73a02d31e 100644 --- a/assimilation_code/modules/assimilation/assim_tools_mod.f90 +++ b/assimilation_code/modules/assimilation/assim_tools_mod.f90 @@ -76,7 +76,7 @@ module assim_tools_mod use normal_distribution_mod, only : normal_cdf, inv_weighted_normal_cdf -use algorithm_info_mod, only : probit_dist_info, obs_inc_info +use algorithm_info_mod, only : probit_dist_info, obs_inc_info, init_qcf_table use gamma_distribution_mod, only : gamma_cdf, inv_gamma_cdf, gamma_mn_var_to_shape_scale, & gamma_gamma_prod @@ -143,6 +143,7 @@ module assim_tools_mod ! special_localization_obs_types -> Special treatment for the specified observation types ! special_localization_cutoffs -> Different cutoff value for each specified obs type ! +character(len = 129) :: qcf_table_filename = '' !not sure if the len should be 129 here, but it is consistent with other nml variables logical :: use_algorithm_info_mod = .true. integer :: filter_kind = 1 real(r8) :: cutoff = 0.2_r8 @@ -203,7 +204,7 @@ module assim_tools_mod ! compared to previous versions of this namelist item. logical :: distribute_mean = .false. -namelist / assim_tools_nml / use_algorithm_info_mod, & +namelist / assim_tools_nml / qcf_table_filename, use_algorithm_info_mod, & filter_kind, cutoff, sort_obs_inc, & spread_restoration, sampling_error_correction, & adaptive_localization_threshold, adaptive_cutoff_floor, & @@ -222,7 +223,7 @@ module assim_tools_mod subroutine assim_tools_init() -integer :: iunit, io, i, j +integer :: iunit, io, i, j, numrows integer :: num_special_cutoff, type_index logical :: cache_override = .false. @@ -313,6 +314,10 @@ subroutine assim_tools_init() call log_namelist_selections(num_special_cutoff, cache_override) +if(qcf_table_filename) then + call init_qcf_table(qcf_table_filename, numrows) +endif + end subroutine assim_tools_init !------------------------------------------------------------- From 1f9dcc31ce92c27cf601b7e2ec6a2c1d88e0a700 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Mon, 14 Aug 2023 15:32:48 -0600 Subject: [PATCH 05/60] Adding QCF table type definitions to algorithm_info_mod --- .../assimilation/algorithm_info_mod.f90 | 46 ++++++++++++++++++- 1 file changed, 45 insertions(+), 1 deletion(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index 78c2f190ab..fb6a95c1f0 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -35,9 +35,51 @@ module algorithm_info_mod integer, parameter :: BOUNDED_NORMAL_RHF = 101 public :: obs_error_info, probit_dist_info, obs_inc_info, & - init_qcf_table, & + init_qcf_table, read_qcf_table, & + obs_error_info_type, probit_inflation_type, probit_state_type, & + probit_extended_state_type, obs_inc_info_type, qcf_table_data_type, & EAKF, ENKF, BOUNDED_NORMAL_RHF, UNBOUNDED_RHF, GAMMA_FILTER +!Creates the type definitions for the QCF table +type obs_error_info_type + logical :: bounded_below, bounded_above + real :: lower_bound, upper_bound +end type + +type probit_inflation_type + integer :: dist_type + logical :: bounded_below, bounded_above + real :: lower_bound, upper_bound +end type + +type probit_state_type + integer :: dist_type + logical :: bounded_below, bounded_above + real :: lower_bound, upper_bound +end type + +type probit_extended_state_type + integer :: dist_type + logical :: bounded_below, bounded_above + real :: lower_bound, upper_bound +end type + +type obs_inc_info_type + integer :: filter_kind + logical :: rectangular_quadrature, gaussian_likelihood_tails + logical :: sort_obs_inc, spread_restoration + logical :: bounded_below, bounded_above + real :: lower_bound, upper_bound +end type + +type qcf_table_data_type + type(obs_error_info_type) :: obs_error_info + type(probit_inflation_type) :: probit_inflation + type(probit_state_type) :: probit_state + type(probit_extended_state_type) :: probit_extended_state + type(obs_inc_info_type) :: obs_inc_info +end type + ! Provides routines that give information about details of algorithms for ! observation error sampling, observation increments, and the transformations ! for regression and inflation in probit space. @@ -266,4 +308,6 @@ end subroutine init_qcf_table !------------------------------------------------------------------------ + + end module algorithm_info_mod From e8ff87f67e25bef25801f63e9ad4854721816a3e Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Mon, 14 Aug 2023 15:38:16 -0600 Subject: [PATCH 06/60] adding type defs to use statement for algorithm_info_mod in assim_tools_mod --- assimilation_code/modules/assimilation/assim_tools_mod.f90 | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/assimilation_code/modules/assimilation/assim_tools_mod.f90 b/assimilation_code/modules/assimilation/assim_tools_mod.f90 index c73a02d31e..1c1d6bc68e 100644 --- a/assimilation_code/modules/assimilation/assim_tools_mod.f90 +++ b/assimilation_code/modules/assimilation/assim_tools_mod.f90 @@ -76,7 +76,9 @@ module assim_tools_mod use normal_distribution_mod, only : normal_cdf, inv_weighted_normal_cdf -use algorithm_info_mod, only : probit_dist_info, obs_inc_info, init_qcf_table +use algorithm_info_mod, only : probit_dist_info, obs_inc_info, init_qcf_table, & + qcf_table_data_type, obs_error_info_type, obs_inc_info_type, & + probit_inflation_type, probit_state_type, probit_extended_state_type use gamma_distribution_mod, only : gamma_cdf, inv_gamma_cdf, gamma_mn_var_to_shape_scale, & gamma_gamma_prod From 7024e2b4e41c307cd0744f8a72ec76266bd96d30 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Mon, 14 Aug 2023 15:47:21 -0600 Subject: [PATCH 07/60] Adding allocatable variables for table data, allocating after determining size of table --- assimilation_code/modules/assimilation/assim_tools_mod.f90 | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/assimilation_code/modules/assimilation/assim_tools_mod.f90 b/assimilation_code/modules/assimilation/assim_tools_mod.f90 index 1c1d6bc68e..d29717c1f2 100644 --- a/assimilation_code/modules/assimilation/assim_tools_mod.f90 +++ b/assimilation_code/modules/assimilation/assim_tools_mod.f90 @@ -129,6 +129,9 @@ module assim_tools_mod character(len=*), parameter :: source = 'assim_tools_mod.f90' +type(qcf_table_data_type), allocatable :: qcf_table_data(:) +character(len=129), allocatable :: qcf_table_row_headers(:) !!!!! might need to change len=129 + !============================================================================ !---- namelist with default values @@ -318,6 +321,8 @@ subroutine assim_tools_init() if(qcf_table_filename) then call init_qcf_table(qcf_table_filename, numrows) + allocate(qcf_table_row_headers(numrows)) + allocate(qcf_table_data(numrows)) endif end subroutine assim_tools_init From 935bcb535850363001aeab2223b7d25b67eed572 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Mon, 14 Aug 2023 16:12:59 -0600 Subject: [PATCH 08/60] New subroutine to read through the values in the QCF table and assign them to the variables in the qcf_table_data_type --- .../assimilation/algorithm_info_mod.f90 | 60 +++++++++++++++++++ .../modules/assimilation/assim_tools_mod.f90 | 3 +- 2 files changed, 62 insertions(+), 1 deletion(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index fb6a95c1f0..8a795717e1 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -309,5 +309,65 @@ end subroutine init_qcf_table !------------------------------------------------------------------------ +subroutine read_qcf_table(qcf_table_filename, numrows, qcf_table_data, rowheaders) + +! Reads in the QCEFF input options from tabular data file + +character(len=129), intent(in) :: qcf_table_filename +integer, intent(in) :: numrows +type(qcf_table_data_type), intent(inout) :: qcf_table_data(:) +character(len=129), intent(inout) :: rowheaders(:) !!!!! might need to change len=129 + +integer, parameter :: fileid = 10 !file identifier +integer :: row + +character(len=129), dimension(4) :: header1 +character(len=129), dimension(29) :: header2 + +open(unit=fileid, file=qcf_table_filename) + +read(fileid, *) header1 +read(fileid, *) header2 !! skip the headers +write(*, *) "header1: ", header1 +write(*, *) "header2: ", header2 + +! read in table values directly to qcf_table_data type +do row = 1, numrows + read(fileid, *) rowheaders(row), qcf_table_data(row)%obs_error_info%bounded_below, qcf_table_data(row)%obs_error_info%bounded_above, & + qcf_table_data(row)%obs_error_info%lower_bound, qcf_table_data(row)%obs_error_info%upper_bound, qcf_table_data(row)%probit_inflation%dist_type, & + qcf_table_data(row)%probit_inflation%bounded_below, qcf_table_data(row)%probit_inflation%bounded_above, & + qcf_table_data(row)%probit_inflation%lower_bound, qcf_table_data(row)%probit_inflation%upper_bound, qcf_table_data(row)%probit_state%dist_type, & + qcf_table_data(row)%probit_state%bounded_below, qcf_table_data(row)%probit_state%bounded_above, & + qcf_table_data(row)%probit_state%lower_bound, qcf_table_data(row)%probit_state%upper_bound, qcf_table_data(row)%probit_extended_state%dist_type, & + qcf_table_data(row)%probit_extended_state%bounded_below, qcf_table_data(row)%probit_extended_state%bounded_above, & + qcf_table_data(row)%probit_extended_state%lower_bound, qcf_table_data(row)%probit_extended_state%upper_bound, & + qcf_table_data(row)%obs_inc_info%filter_kind, qcf_table_data(row)%obs_inc_info%rectangular_quadrature, & + qcf_table_data(row)%obs_inc_info%gaussian_likelihood_tails, qcf_table_data(row)%obs_inc_info%sort_obs_inc, & + qcf_table_data(row)%obs_inc_info%spread_restoration, qcf_table_data(row)%obs_inc_info%bounded_below, qcf_table_data(row)%obs_inc_info%bounded_above, & + qcf_table_data(row)%obs_inc_info%lower_bound, qcf_table_data(row)%obs_inc_info%upper_bound + +! write to check values were correctly assigned + write(*, *) "rowheader(", row, "): ", rowheaders(row) + write(*, *) "qcf_table_data(", row, "): " + write(*, *) qcf_table_data(row)%obs_error_info%bounded_below, qcf_table_data(row)%obs_error_info%bounded_above, & + qcf_table_data(row)%obs_error_info%lower_bound, qcf_table_data(row)%obs_error_info%upper_bound, qcf_table_data(row)%probit_inflation%dist_type, & + qcf_table_data(row)%probit_inflation%bounded_below, qcf_table_data(row)%probit_inflation%bounded_above, & + qcf_table_data(row)%probit_inflation%lower_bound, qcf_table_data(row)%probit_inflation%upper_bound, qcf_table_data(row)%probit_state%dist_type, & + qcf_table_data(row)%probit_state%bounded_below, qcf_table_data(row)%probit_state%bounded_above, & + qcf_table_data(row)%probit_state%lower_bound, qcf_table_data(row)%probit_state%upper_bound, qcf_table_data(row)%probit_extended_state%dist_type, & + qcf_table_data(row)%probit_extended_state%bounded_below, qcf_table_data(row)%probit_extended_state%bounded_above, & + qcf_table_data(row)%probit_extended_state%lower_bound, qcf_table_data(row)%probit_extended_state%upper_bound, & + qcf_table_data(row)%obs_inc_info%filter_kind, qcf_table_data(row)%obs_inc_info%rectangular_quadrature, & + qcf_table_data(row)%obs_inc_info%gaussian_likelihood_tails, qcf_table_data(row)%obs_inc_info%sort_obs_inc, & + qcf_table_data(row)%obs_inc_info%spread_restoration, qcf_table_data(row)%obs_inc_info%bounded_below, qcf_table_data(row)%obs_inc_info%bounded_above, & + qcf_table_data(row)%obs_inc_info%lower_bound, qcf_table_data(row)%obs_inc_info%upper_bound +end do + +close(fileid) + + +end subroutine read_qcf_table + +!------------------------------------------------------------------------ end module algorithm_info_mod diff --git a/assimilation_code/modules/assimilation/assim_tools_mod.f90 b/assimilation_code/modules/assimilation/assim_tools_mod.f90 index d29717c1f2..1cc8d7a909 100644 --- a/assimilation_code/modules/assimilation/assim_tools_mod.f90 +++ b/assimilation_code/modules/assimilation/assim_tools_mod.f90 @@ -76,7 +76,7 @@ module assim_tools_mod use normal_distribution_mod, only : normal_cdf, inv_weighted_normal_cdf -use algorithm_info_mod, only : probit_dist_info, obs_inc_info, init_qcf_table, & +use algorithm_info_mod, only : probit_dist_info, obs_inc_info, init_qcf_table, read_qcf_table, & qcf_table_data_type, obs_error_info_type, obs_inc_info_type, & probit_inflation_type, probit_state_type, probit_extended_state_type @@ -323,6 +323,7 @@ subroutine assim_tools_init() call init_qcf_table(qcf_table_filename, numrows) allocate(qcf_table_row_headers(numrows)) allocate(qcf_table_data(numrows)) + call read_qcf_table(qcf_table_filename, numrows, qcf_table_data, qcf_table_row_headers) endif end subroutine assim_tools_init From be4ee6395eb40b4ff6af2a440d95107faae6f45c Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Mon, 14 Aug 2023 16:47:08 -0600 Subject: [PATCH 09/60] Removing qcf table data types from assim_tools_mod and reorganizing so that these type structs are only used in algorithm_info_mod --- .../modules/assimilation/algorithm_info_mod.f90 | 17 +++++++++++++---- .../modules/assimilation/assim_tools_mod.f90 | 12 ++---------- 2 files changed, 15 insertions(+), 14 deletions(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index 8a795717e1..521d5dca4f 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -80,6 +80,9 @@ module algorithm_info_mod type(obs_inc_info_type) :: obs_inc_info end type +type(qcf_table_data_type), allocatable :: qcf_table_data(:) +character(len=129), allocatable :: qcf_table_row_headers(:) !!!!! might need to change len=129 + ! Provides routines that give information about details of algorithms for ! observation error sampling, observation increments, and the transformations ! for regression and inflation in probit space. @@ -282,11 +285,11 @@ end subroutine obs_inc_info !------------------------------------------------------------------------ -subroutine init_qcf_table(qcf_table_filename, numrows) +subroutine init_qcf_table(qcf_table_filename) character(len=50), intent(in) :: qcf_table_filename -integer, intent(out) :: numrows !return value +integer :: numrows integer :: nlines integer :: io integer, parameter :: fileid = 10 !file identifier @@ -304,6 +307,11 @@ subroutine init_qcf_table(qcf_table_filename, numrows) numrows = nlines - 2 print *, 'numrows: ', numrows +allocate(qcf_table_data(numrows)) +allocate(rowheaders(numrows)) + +call read_qcf_table(qcf_table_filename, numrows, qcf_table_data, rowheaders) + end subroutine init_qcf_table !------------------------------------------------------------------------ @@ -315,8 +323,9 @@ subroutine read_qcf_table(qcf_table_filename, numrows, qcf_table_data, rowheader character(len=129), intent(in) :: qcf_table_filename integer, intent(in) :: numrows -type(qcf_table_data_type), intent(inout) :: qcf_table_data(:) -character(len=129), intent(inout) :: rowheaders(:) !!!!! might need to change len=129 + +type(qcf_table_data_type) :: qcf_table_data(:) +character(len=129) :: rowheaders(:) !!!!! might need to change len=129 integer, parameter :: fileid = 10 !file identifier integer :: row diff --git a/assimilation_code/modules/assimilation/assim_tools_mod.f90 b/assimilation_code/modules/assimilation/assim_tools_mod.f90 index 1cc8d7a909..d37556afcc 100644 --- a/assimilation_code/modules/assimilation/assim_tools_mod.f90 +++ b/assimilation_code/modules/assimilation/assim_tools_mod.f90 @@ -76,9 +76,7 @@ module assim_tools_mod use normal_distribution_mod, only : normal_cdf, inv_weighted_normal_cdf -use algorithm_info_mod, only : probit_dist_info, obs_inc_info, init_qcf_table, read_qcf_table, & - qcf_table_data_type, obs_error_info_type, obs_inc_info_type, & - probit_inflation_type, probit_state_type, probit_extended_state_type +use algorithm_info_mod, only : probit_dist_info, obs_inc_info, init_qcf_table, read_qcf_table use gamma_distribution_mod, only : gamma_cdf, inv_gamma_cdf, gamma_mn_var_to_shape_scale, & gamma_gamma_prod @@ -129,9 +127,6 @@ module assim_tools_mod character(len=*), parameter :: source = 'assim_tools_mod.f90' -type(qcf_table_data_type), allocatable :: qcf_table_data(:) -character(len=129), allocatable :: qcf_table_row_headers(:) !!!!! might need to change len=129 - !============================================================================ !---- namelist with default values @@ -320,10 +315,7 @@ subroutine assim_tools_init() call log_namelist_selections(num_special_cutoff, cache_override) if(qcf_table_filename) then - call init_qcf_table(qcf_table_filename, numrows) - allocate(qcf_table_row_headers(numrows)) - allocate(qcf_table_data(numrows)) - call read_qcf_table(qcf_table_filename, numrows, qcf_table_data, qcf_table_row_headers) + call init_qcf_table(qcf_table_filename) endif end subroutine assim_tools_init From 66356a67c740717329a14a614ea5069397b6c4b2 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Mon, 14 Aug 2023 17:15:53 -0600 Subject: [PATCH 10/60] Fixing small inconsistencies/typos --- .../modules/assimilation/algorithm_info_mod.f90 | 6 +++--- assimilation_code/modules/assimilation/assim_tools_mod.f90 | 7 +++++-- 2 files changed, 8 insertions(+), 5 deletions(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index 521d5dca4f..77a35392c4 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -287,7 +287,7 @@ end subroutine obs_inc_info subroutine init_qcf_table(qcf_table_filename) -character(len=50), intent(in) :: qcf_table_filename +character(len=129), intent(in) :: qcf_table_filename integer :: numrows integer :: nlines @@ -308,9 +308,9 @@ subroutine init_qcf_table(qcf_table_filename) print *, 'numrows: ', numrows allocate(qcf_table_data(numrows)) -allocate(rowheaders(numrows)) +allocate(qcf_table_row_headers(numrows)) -call read_qcf_table(qcf_table_filename, numrows, qcf_table_data, rowheaders) +call read_qcf_table(qcf_table_filename, numrows, qcf_table_data, qcf_table_row_headers) end subroutine init_qcf_table diff --git a/assimilation_code/modules/assimilation/assim_tools_mod.f90 b/assimilation_code/modules/assimilation/assim_tools_mod.f90 index d37556afcc..608fe799a3 100644 --- a/assimilation_code/modules/assimilation/assim_tools_mod.f90 +++ b/assimilation_code/modules/assimilation/assim_tools_mod.f90 @@ -76,7 +76,7 @@ module assim_tools_mod use normal_distribution_mod, only : normal_cdf, inv_weighted_normal_cdf -use algorithm_info_mod, only : probit_dist_info, obs_inc_info, init_qcf_table, read_qcf_table +use algorithm_info_mod, only : probit_dist_info, obs_inc_info, init_qcf_table use gamma_distribution_mod, only : gamma_cdf, inv_gamma_cdf, gamma_mn_var_to_shape_scale, & gamma_gamma_prod @@ -314,7 +314,10 @@ subroutine assim_tools_init() call log_namelist_selections(num_special_cutoff, cache_override) -if(qcf_table_filename) then +write(*,*), "HERE!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!" +if(qcf_table_filename == '') then + write(*,*), "no qcf table in namelist" +else call init_qcf_table(qcf_table_filename) endif From b62339e92adffe247ad031a84533c9993811d1ea Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Mon, 14 Aug 2023 17:18:41 -0600 Subject: [PATCH 11/60] moving the location of draft program outside /assimilation_code/modules/assimilation --- .../modules/assimilation => qcf_table}/type_read_table.f90 | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) rename {assimilation_code/modules/assimilation => qcf_table}/type_read_table.f90 (97%) diff --git a/assimilation_code/modules/assimilation/type_read_table.f90 b/qcf_table/type_read_table.f90 similarity index 97% rename from assimilation_code/modules/assimilation/type_read_table.f90 rename to qcf_table/type_read_table.f90 index f164221db3..f7cd548acf 100644 --- a/assimilation_code/modules/assimilation/type_read_table.f90 +++ b/qcf_table/type_read_table.f90 @@ -45,7 +45,7 @@ program read_table !real(r8), intent(out) :: qcf_table_data !real, dimension(:, :), allocatable :: qcf_table_data_rows type(qcf_table_data_type), allocatable :: qcf_table_data(:) -character(len=30), dimension(:), allocatable :: rowheaders !!!!! might need to change len=30 +character(len=129), dimension(:), allocatable :: rowheaders !!!!! might need to change len=30 integer, parameter :: fileid = 10 !file identifier character(len=30), parameter :: tester_QTY = 'QTY_GPSRO' @@ -63,7 +63,7 @@ program read_table character(len=30), dimension(29) :: header2 !variables for table values ^^^ -open(unit=fileid, file='cam_qcf_table.dat') +open(unit=fileid, file='cam_qcf_table.txt') nlines = 0 do !do loop to get number of rows (or QTY's) in the table @@ -82,7 +82,7 @@ program read_table allocate(rowheaders(numrows)) write(*,*) shape(qcf_table_data) -open(unit=fileid, file='cam_qcf_table.dat') +open(unit=fileid, file='cam_qcf_table.txt') read(fileid, *) header1 read(fileid, *) header2 !! skip the headers From 7659472cf7e16c08fc910be746ba4258672ae950 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Tue, 15 Aug 2023 15:14:06 -0600 Subject: [PATCH 12/60] Adding draft subroutine write_qcf_table to test that values are being read in correctly; removed rowheaders argument from subroutines where not needed --- .../assimilation/algorithm_info_mod.f90 | 96 ++++++++++++------- 1 file changed, 63 insertions(+), 33 deletions(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index 77a35392c4..4f47c64e3f 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -155,6 +155,8 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & logical, intent(out) :: bounded_below, bounded_above real(r8), intent(out) :: lower_bound, upper_bound +integer :: QTY_loc(1) + ! Have input information about the kind of the state or observation being transformed ! along with additional logical info that indicates whether this is an observation ! or state variable and about whether the transformation is being done for inflation @@ -176,18 +178,19 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & if(is_inflation) then ! Case for inflation transformation - if(kind == QTY_STATE_VARIABLE) then - dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION - bounded_below = .false.; bounded_above = .false. - lower_bound = missing_r8; upper_bound = missing_r8 - elseif(kind == QTY_TRACER_CONCENTRATION) then - dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION - bounded_below = .true.; bounded_above = .false. - lower_bound = 0.0_r8; upper_bound = missing_r8 - elseif(kind == QTY_TRACER_SOURCE) then - dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION - bounded_below = .true.; bounded_above = .false. - lower_bound = 0.0_r8; upper_bound = missing_r8 +! if(kind == QTY_STATE_VARIABLE) then +! dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION +! bounded_below = .false.; bounded_above = .false. +! lower_bound = missing_r8; upper_bound = missing_r8 +! elseif(kind == QTY_TRACER_CONCENTRATION) then +! dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION +! bounded_below = .true.; bounded_above = .false. +! lower_bound = 0.0_r8; upper_bound = missing_r8 +! elseif(kind == QTY_TRACER_SOURCE) then +! dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION +! bounded_below = .true.; bounded_above = .false. +! lower_bound = 0.0_r8; upper_bound = missing_r8 + if findloc else write(*, *) 'Illegal kind in obs_error_info' stop @@ -310,22 +313,21 @@ subroutine init_qcf_table(qcf_table_filename) allocate(qcf_table_data(numrows)) allocate(qcf_table_row_headers(numrows)) -call read_qcf_table(qcf_table_filename, numrows, qcf_table_data, qcf_table_row_headers) +call read_qcf_table(qcf_table_filename) end subroutine init_qcf_table !------------------------------------------------------------------------ -subroutine read_qcf_table(qcf_table_filename, numrows, qcf_table_data, rowheaders) +subroutine read_qcf_table(qcf_table_filename) ! Reads in the QCEFF input options from tabular data file character(len=129), intent(in) :: qcf_table_filename -integer, intent(in) :: numrows -type(qcf_table_data_type) :: qcf_table_data(:) -character(len=129) :: rowheaders(:) !!!!! might need to change len=129 +!type(qcf_table_data_type) :: qcf_table_data(:) +!character(len=129) :: rowheaders(:) !!!!! might need to change len=129 integer, parameter :: fileid = 10 !file identifier integer :: row @@ -341,8 +343,8 @@ subroutine read_qcf_table(qcf_table_filename, numrows, qcf_table_data, rowheader write(*, *) "header2: ", header2 ! read in table values directly to qcf_table_data type -do row = 1, numrows - read(fileid, *) rowheaders(row), qcf_table_data(row)%obs_error_info%bounded_below, qcf_table_data(row)%obs_error_info%bounded_above, & +do row = 1, size(qcf_table_data) + read(fileid, *) qcf_table_row_headers(row), qcf_table_data(row)%obs_error_info%bounded_below, qcf_table_data(row)%obs_error_info%bounded_above, & qcf_table_data(row)%obs_error_info%lower_bound, qcf_table_data(row)%obs_error_info%upper_bound, qcf_table_data(row)%probit_inflation%dist_type, & qcf_table_data(row)%probit_inflation%bounded_below, qcf_table_data(row)%probit_inflation%bounded_above, & qcf_table_data(row)%probit_inflation%lower_bound, qcf_table_data(row)%probit_inflation%upper_bound, qcf_table_data(row)%probit_state%dist_type, & @@ -354,28 +356,56 @@ subroutine read_qcf_table(qcf_table_filename, numrows, qcf_table_data, rowheader qcf_table_data(row)%obs_inc_info%gaussian_likelihood_tails, qcf_table_data(row)%obs_inc_info%sort_obs_inc, & qcf_table_data(row)%obs_inc_info%spread_restoration, qcf_table_data(row)%obs_inc_info%bounded_below, qcf_table_data(row)%obs_inc_info%bounded_above, & qcf_table_data(row)%obs_inc_info%lower_bound, qcf_table_data(row)%obs_inc_info%upper_bound +end do + +close(fileid) + +call write_qcf_table() + +end subroutine read_qcf_table + +!------------------------------------------------------------------------ + + +subroutine write_qcf_table() ! write to check values were correctly assigned - write(*, *) "rowheader(", row, "): ", rowheaders(row) +! testing for findloc + +character(len=30), parameter :: tester_QTY = 'QTY_GPSRO' +integer :: QTY_loc(1) + +character(len=30), parameter :: tester_QTY0 = 'QTY_DUMMY' +integer :: QTY_loc0(1) + +integer :: row + +write(*,*), 'SIZE: ', size(qcf_table_data) + +do row = 1, size(qcf_table_data) + write(*, *) "qcf_table_row_headers(", row, "): ", qcf_table_row_headers(row) write(*, *) "qcf_table_data(", row, "): " write(*, *) qcf_table_data(row)%obs_error_info%bounded_below, qcf_table_data(row)%obs_error_info%bounded_above, & - qcf_table_data(row)%obs_error_info%lower_bound, qcf_table_data(row)%obs_error_info%upper_bound, qcf_table_data(row)%probit_inflation%dist_type, & - qcf_table_data(row)%probit_inflation%bounded_below, qcf_table_data(row)%probit_inflation%bounded_above, & - qcf_table_data(row)%probit_inflation%lower_bound, qcf_table_data(row)%probit_inflation%upper_bound, qcf_table_data(row)%probit_state%dist_type, & - qcf_table_data(row)%probit_state%bounded_below, qcf_table_data(row)%probit_state%bounded_above, & - qcf_table_data(row)%probit_state%lower_bound, qcf_table_data(row)%probit_state%upper_bound, qcf_table_data(row)%probit_extended_state%dist_type, & - qcf_table_data(row)%probit_extended_state%bounded_below, qcf_table_data(row)%probit_extended_state%bounded_above, & - qcf_table_data(row)%probit_extended_state%lower_bound, qcf_table_data(row)%probit_extended_state%upper_bound, & - qcf_table_data(row)%obs_inc_info%filter_kind, qcf_table_data(row)%obs_inc_info%rectangular_quadrature, & - qcf_table_data(row)%obs_inc_info%gaussian_likelihood_tails, qcf_table_data(row)%obs_inc_info%sort_obs_inc, & - qcf_table_data(row)%obs_inc_info%spread_restoration, qcf_table_data(row)%obs_inc_info%bounded_below, qcf_table_data(row)%obs_inc_info%bounded_above, & - qcf_table_data(row)%obs_inc_info%lower_bound, qcf_table_data(row)%obs_inc_info%upper_bound + qcf_table_data(row)%obs_error_info%lower_bound, qcf_table_data(row)%obs_error_info%upper_bound, qcf_table_data(row)%probit_inflation%dist_type, & + qcf_table_data(row)%probit_inflation%bounded_below, qcf_table_data(row)%probit_inflation%bounded_above, & + qcf_table_data(row)%probit_inflation%lower_bound, qcf_table_data(row)%probit_inflation%upper_bound, qcf_table_data(row)%probit_state%dist_type, & + qcf_table_data(row)%probit_state%bounded_below, qcf_table_data(row)%probit_state%bounded_above, & + qcf_table_data(row)%probit_state%lower_bound, qcf_table_data(row)%probit_state%upper_bound, qcf_table_data(row)%probit_extended_state%dist_type, & + qcf_table_data(row)%probit_extended_state%bounded_below, qcf_table_data(row)%probit_extended_state%bounded_above, & + qcf_table_data(row)%probit_extended_state%lower_bound, qcf_table_data(row)%probit_extended_state%upper_bound, & + qcf_table_data(row)%obs_inc_info%filter_kind, qcf_table_data(row)%obs_inc_info%rectangular_quadrature, & + qcf_table_data(row)%obs_inc_info%gaussian_likelihood_tails, qcf_table_data(row)%obs_inc_info%sort_obs_inc, & + qcf_table_data(row)%obs_inc_info%spread_restoration, qcf_table_data(row)%obs_inc_info%bounded_below, qcf_table_data(row)%obs_inc_info%bounded_above, & + qcf_table_data(row)%obs_inc_info%lower_bound, qcf_table_data(row)%obs_inc_info%upper_bound end do -close(fileid) +QTY_loc = findloc(qcf_table_row_headers, tester_QTY) +write(*, *) 'findloc of QTY_GPSRO: ', QTY_loc(1) +QTY_loc0 = findloc(qcf_table_row_headers, tester_QTY0) +write(*, *) 'findloc of invalid QTY (QTY_DUMMY): ', QTY_loc0(1) -end subroutine read_qcf_table +end subroutine write_qcf_table !------------------------------------------------------------------------ From e9955639da2e635e0b7bac2c1d8d11fb043f3ff0 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Tue, 15 Aug 2023 16:06:58 -0600 Subject: [PATCH 13/60] replaicing conditionals and hardcoded values in probit_dist_info --- .../assimilation/algorithm_info_mod.f90 | 113 ++++++++++++------ 1 file changed, 74 insertions(+), 39 deletions(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index 4f47c64e3f..c80f94ccde 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -176,8 +176,26 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & ! In the long run, may not have to have separate controls for each of the input possibilities ! However, for now these are things that need to be explored for science understanding -if(is_inflation) then +QTY_loc = findloc(qcf_table_row_headers, kind) +write(*, *) 'findloc of kind: ', QTY_loc(1) + +if (QTY_loc(1) == 0) then + write(*,*) 'QTY not in table!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!' + + !using default values here + dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION + bounded_below = .false.; bounded_above = .false. + lower_bound = missing_r8; upper_bound = missing_r8 + + elseif(is_inflation) then ! Case for inflation transformation + + dist_type = qcf_table_data(QTY_loc(1))%probit_inflation%dist_type + bounded_below = qcf_table_data(QTY_loc(1))%probit_inflation%bounded_below + bounded_above = qcf_table_data(QTY_loc(1))%probit_inflation%bounded_above + lower_bound = qcf_table_data(QTY_loc(1))%probit_inflation%lower_bound + upper_bound = qcf_table_data(QTY_loc(1))%probit_inflation%upper_bound + ! if(kind == QTY_STATE_VARIABLE) then ! dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION ! bounded_below = .false.; bounded_above = .false. @@ -190,49 +208,66 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & ! dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION ! bounded_below = .true.; bounded_above = .false. ! lower_bound = 0.0_r8; upper_bound = missing_r8 - if findloc - else - write(*, *) 'Illegal kind in obs_error_info' - stop - endif -elseif(is_state) then +! else +! write(*, *) 'Illegal kind in obs_error_info' +! stop +! endif + + elseif(is_state) then ! Case for state variable priors - if(kind == QTY_STATE_VARIABLE) then - dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION - bounded_below = .false.; bounded_above = .false. - lower_bound = missing_r8; upper_bound = missing_r8 - elseif(kind == QTY_TRACER_CONCENTRATION) then - dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION - bounded_below = .true.; bounded_above = .false. - lower_bound = 0.0_r8; upper_bound = missing_r8 - elseif(kind == QTY_TRACER_SOURCE) then - dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION - bounded_below = .true.; bounded_above = .false. - lower_bound = 0.0_r8; upper_bound = missing_r8 + + dist_type = qcf_table_data(QTY_loc(1))%probit_state%dist_type + bounded_below = qcf_table_data(QTY_loc(1))%probit_state%bounded_below + bounded_above = qcf_table_data(QTY_loc(1))%probit_state%bounded_above + lower_bound = qcf_table_data(QTY_loc(1))%probit_state%lower_bound + upper_bound = qcf_table_data(QTY_loc(1))%probit_state%upper_bound + + !if(kind == QTY_STATE_VARIABLE) then + ! dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION + ! bounded_below = .false.; bounded_above = .false. + ! lower_bound = missing_r8; upper_bound = missing_r8 + ! elseif(kind == QTY_TRACER_CONCENTRATION) then + ! dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION + ! bounded_below = .true.; bounded_above = .false. + ! lower_bound = 0.0_r8; upper_bound = missing_r8 +! elseif(kind == QTY_TRACER_SOURCE) then + ! dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION + ! bounded_below = .true.; bounded_above = .false. + ! lower_bound = 0.0_r8; upper_bound = missing_r8 +! else +! write(*, *) 'Illegal kind in obs_error_info' + ! stop + ! endif + else - write(*, *) 'Illegal kind in obs_error_info' - stop - endif -else ! This case is for observation (extended state) priors - if(kind == QTY_STATE_VARIABLE) then - dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION - bounded_below = .false.; bounded_above = .false. - lower_bound = missing_r8; upper_bound = missing_r8 - elseif(kind == QTY_TRACER_CONCENTRATION) then - dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION - bounded_below = .true.; bounded_above = .false. - lower_bound = 0.0_r8; upper_bound = missing_r8 - elseif(kind == QTY_TRACER_SOURCE) then - dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION - bounded_below = .true.; bounded_above = .false. - lower_bound = 0.0_r8; upper_bound = missing_r8 - else - write(*, *) 'Illegal kind in obs_error_info' - stop - endif + + dist_type = qcf_table_data(QTY_loc(1))%probit_extended_state%dist_type + bounded_below = qcf_table_data(QTY_loc(1))%probit_extended_state%bounded_below + bounded_above = qcf_table_data(QTY_loc(1))%probit_extended_state%bounded_above + lower_bound = qcf_table_data(QTY_loc(1))%probit_extended_state%lower_bound + upper_bound = qcf_table_data(QTY_loc(1))%probit_extended_state%upper_bound + +! if(kind == QTY_STATE_VARIABLE) then +! dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION +! bounded_below = .false.; bounded_above = .false. +! lower_bound = missing_r8; upper_bound = missing_r8 +! elseif(kind == QTY_TRACER_CONCENTRATION) then +! dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION +! bounded_below = .true.; bounded_above = .false. +! lower_bound = 0.0_r8; upper_bound = missing_r8 +! elseif(kind == QTY_TRACER_SOURCE) then +! dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION +! bounded_below = .true.; bounded_above = .false. +! lower_bound = 0.0_r8; upper_bound = missing_r8 +! else +! write(*, *) 'Illegal kind in obs_error_info' +! stop +! endif endif +write(*,*) dist_type, bounded_below, bounded_above, lower_bound, upper_bound + end subroutine probit_dist_info !------------------------------------------------------------------------ From fa514c254c340923e38b851e316ab3e68ed01443 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Wed, 16 Aug 2023 14:26:13 -0600 Subject: [PATCH 14/60] using get_name_for_quantity to get generic quantity from integer index --- .../assimilation/algorithm_info_mod.f90 | 41 +++++++++++++++---- 1 file changed, 33 insertions(+), 8 deletions(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index c80f94ccde..83937b807f 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -7,7 +7,8 @@ module algorithm_info_mod use types_mod, only : r8, i8, missing_r8 use obs_def_mod, only : obs_def_type, get_obs_def_type_of_obs, get_obs_def_error_variance -use obs_kind_mod, only : get_quantity_for_type_of_obs +use obs_kind_mod, only : get_quantity_for_type_of_obs, get_name_for_quantity + ! Get the QTY definitions that are needed (aka kind) use obs_kind_mod, only : QTY_STATE_VARIABLE, QTY_STATE_VAR_POWER, QTY_TRACER_CONCENTRATION, & @@ -156,6 +157,7 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & real(r8), intent(out) :: lower_bound, upper_bound integer :: QTY_loc(1) +character(len=129) :: kind_name ! Have input information about the kind of the state or observation being transformed ! along with additional logical info that indicates whether this is an observation @@ -176,13 +178,18 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & ! In the long run, may not have to have separate controls for each of the input possibilities ! However, for now these are things that need to be explored for science understanding -QTY_loc = findloc(qcf_table_row_headers, kind) -write(*, *) 'findloc of kind: ', QTY_loc(1) +!get actual name of QTY from integer index +kind_name = get_name_for_quantity(kind) +write(*,*) 'kind_name: ', kind_name + +!find location of QTY in qcf_table_data structure +QTY_loc = findloc(qcf_table_row_headers, kind_name) +write(*,*) 'findloc of kind: ', QTY_loc(1) if (QTY_loc(1) == 0) then write(*,*) 'QTY not in table!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!' - !using default values here + !use default values dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION bounded_below = .false.; bounded_above = .false. lower_bound = missing_r8; upper_bound = missing_r8 @@ -190,10 +197,10 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & elseif(is_inflation) then ! Case for inflation transformation - dist_type = qcf_table_data(QTY_loc(1))%probit_inflation%dist_type + dist_type = qcf_table_data(QTY_loc(1))%probit_inflation%dist_type !dist_type has checks in transform_to_probit, transform_from_probit bounded_below = qcf_table_data(QTY_loc(1))%probit_inflation%bounded_below bounded_above = qcf_table_data(QTY_loc(1))%probit_inflation%bounded_above - lower_bound = qcf_table_data(QTY_loc(1))%probit_inflation%lower_bound + lower_bound = qcf_table_data(QTY_loc(1))%probit_inflation%lower_bound !NEED TO ADD CHECKS THAT THESE ARE VALID VALUES upper_bound = qcf_table_data(QTY_loc(1))%probit_inflation%upper_bound ! if(kind == QTY_STATE_VARIABLE) then @@ -266,7 +273,7 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & ! endif endif -write(*,*) dist_type, bounded_below, bounded_above, lower_bound, upper_bound +write(*,*) 'probit_dist_info: ', dist_type, bounded_below, bounded_above, lower_bound, upper_bound end subroutine probit_dist_info @@ -291,6 +298,22 @@ subroutine obs_inc_info(obs_kind, filter_kind, rectangular_quadrature, gaussian_ ! Temporary approach for setting the details of how to assimilate this observation ! This example is designed to reproduce the squared forward operator results from paper +!get actual name of QTY from integer index +kind_name = get_name_for_quantity(kind) +write(*,*) 'kind_name: ', kind_name + +!find location of QTY in qcf_table_data structure +QTY_loc = findloc(qcf_table_row_headers, kind_name) +write(*,*) 'findloc of kind: ', QTY_loc(1) + +if (QTY_loc(1) == 0) then + write(*,*) 'QTY not in table!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!' + + !use default values + dist_type = BOUNDED_NORMAL_RHF + bounded_below = .false.; bounded_above = .false. + lower_bound = missing_r8; upper_bound = missing_r8 + ! Set the observation increment details for each type of quantity if(obs_kind == QTY_STATE_VARIABLE) then @@ -372,8 +395,9 @@ subroutine read_qcf_table(qcf_table_filename) open(unit=fileid, file=qcf_table_filename) +! skip the headers read(fileid, *) header1 -read(fileid, *) header2 !! skip the headers +read(fileid, *) header2 write(*, *) "header1: ", header1 write(*, *) "header2: ", header2 @@ -404,6 +428,7 @@ end subroutine read_qcf_table subroutine write_qcf_table() +! DRAFT SUBROUTINE ! write to check values were correctly assigned ! testing for findloc From f291025e4e363c646758fcd5d7c0a7fe7c18e768 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Wed, 16 Aug 2023 14:58:50 -0600 Subject: [PATCH 15/60] Replacing conditionals and hard coded values with qcf_table_data in obs_inc_info subroutine --- .../assimilation/algorithm_info_mod.f90 | 57 ++++++++++++------- 1 file changed, 35 insertions(+), 22 deletions(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index 83937b807f..fe9d1c3e26 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -291,6 +291,9 @@ subroutine obs_inc_info(obs_kind, filter_kind, rectangular_quadrature, gaussian_ logical, intent(inout) :: bounded_below, bounded_above real(r8), intent(inout) :: lower_bound, upper_bound +integer :: QTY_loc(1) +character(len=129) :: kind_name + ! The information arguments are all intent (inout). This means that if they are not set ! here, they retain the default values from the assim_tools_mod namelist. Bounds don't exist ! in that namelist, so default values are set in assim_tools_mod just before the call to here. @@ -299,7 +302,7 @@ subroutine obs_inc_info(obs_kind, filter_kind, rectangular_quadrature, gaussian_ ! This example is designed to reproduce the squared forward operator results from paper !get actual name of QTY from integer index -kind_name = get_name_for_quantity(kind) +kind_name = get_name_for_quantity(obs_kind) write(*,*) 'kind_name: ', kind_name !find location of QTY in qcf_table_data structure @@ -310,32 +313,42 @@ subroutine obs_inc_info(obs_kind, filter_kind, rectangular_quadrature, gaussian_ write(*,*) 'QTY not in table!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!' !use default values - dist_type = BOUNDED_NORMAL_RHF - bounded_below = .false.; bounded_above = .false. - lower_bound = missing_r8; upper_bound = missing_r8 - - -! Set the observation increment details for each type of quantity -if(obs_kind == QTY_STATE_VARIABLE) then filter_kind = BOUNDED_NORMAL_RHF bounded_below = .false.; bounded_above = .false. lower_bound = missing_r8; upper_bound = missing_r8 -elseif(obs_kind == QTY_TRACER_CONCENTRATION) then - filter_kind = BOUNDED_NORMAL_RHF - bounded_below = .true.; bounded_above = .false. - lower_bound = 0.0_r8; upper_bound = missing_r8 -elseif(obs_kind == QTY_TRACER_SOURCE) then - filter_kind = BOUNDED_NORMAL_RHF - bounded_below = .true.; bounded_above = .false. - lower_bound = 0.0_r8; upper_bound = missing_r8 -else - write(*, *) 'Illegal obs_kind in obs_error_info' - stop + sort_obs_inc = .false.; spread_restoration = .false. + ! Default settings for now for Icepack and tracer model tests (sort_obs_inc, spread_restoration) + + else + filter_kind = qcf_table_data(QTY_loc(1))%obs_inc_info%filter_kind !filter_kind has a check in obs_increment + sort_obs_inc = qcf_table_data(QTY_loc(1))%obs_inc_info%sort_obs_inc + spread_restoration = qcf_table_data(QTY_loc(1))%obs_inc_info%spread_restoration + bounded_below = qcf_table_data(QTY_loc(1))%obs_inc_info%bounded_below + bounded_above = qcf_table_data(QTY_loc(1))%obs_inc_info%bounded_above + lower_bound = qcf_table_data(QTY_loc(1))%obs_inc_info%lower_bound !NEED TO ADD CHECKS THAT THESE ARE VALID VALUES + upper_bound = qcf_table_data(QTY_loc(1))%obs_inc_info%upper_bound + endif -! Default settings for now for Icepack and tracer model tests -sort_obs_inc = .false. -spread_restoration = .false. +write(*,*) 'obs_inc_info: ', filter_kind, sort_obs_inc, spread_restoration, bounded_below, bounded_above, lower_bound, upper_bound + +! Set the observation increment details for each type of quantity +!if(obs_kind == QTY_STATE_VARIABLE) then +! filter_kind = BOUNDED_NORMAL_RHF +! bounded_below = .false.; bounded_above = .false. +! lower_bound = missing_r8; upper_bound = missing_r8 +!elseif(obs_kind == QTY_TRACER_CONCENTRATION) then +! filter_kind = BOUNDED_NORMAL_RHF +! bounded_below = .true.; bounded_above = .false. +! lower_bound = 0.0_r8; upper_bound = missing_r8 +!elseif(obs_kind == QTY_TRACER_SOURCE) then +! filter_kind = BOUNDED_NORMAL_RHF +! bounded_below = .true.; bounded_above = .false. +! lower_bound = 0.0_r8; upper_bound = missing_r8 +!else +! write(*, *) 'Illegal obs_kind in obs_error_info' +! stop +!endif ! Only need to set these two for options the original RHF implementation !!!rectangular_quadrature = .true. From 5a26816bcd0097535eb70b51a053168f689a7e40 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Wed, 16 Aug 2023 15:13:52 -0600 Subject: [PATCH 16/60] Replacing conditionals and hard coded values with qcf_table_data in obs_error_info subroutine --- .../assimilation/algorithm_info_mod.f90 | 50 +++++++++++++++---- 1 file changed, 39 insertions(+), 11 deletions(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index fe9d1c3e26..da3260042e 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -109,6 +109,9 @@ subroutine obs_error_info(obs_def, error_variance, & integer(i8) :: state_var_index type(location_type) :: temp_loc +integer :: QTY_loc(1) +character(len=129) :: kind_name + ! Get the kind of the observation obs_type = get_obs_def_type_of_obs(obs_def) ! If it is negative, it is an identity obs @@ -122,21 +125,46 @@ subroutine obs_error_info(obs_def, error_variance, & ! Get the default error variance error_variance = get_obs_def_error_variance(obs_def) -! Set the observation error details for each type of quantity -if(obs_kind == QTY_STATE_VARIABLE) then +!get actual name of QTY from integer index +kind_name = get_name_for_quantity(obs_kind) +write(*,*) 'kind_name: ', kind_name + +!find location of QTY in qcf_table_data structure +QTY_loc = findloc(qcf_table_row_headers, kind_name) +write(*,*) 'findloc of kind: ', QTY_loc(1) + +if (QTY_loc(1) == 0) then + write(*,*) 'QTY not in table!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!' + + !use default values bounded_below = .false.; bounded_above = .false. lower_bound = missing_r8; upper_bound = missing_r8 -elseif(obs_kind == QTY_TRACER_CONCENTRATION) then - bounded_below = .true.; bounded_above = .false. - lower_bound = 0.0_r8; upper_bound = missing_r8 -elseif(obs_kind == QTY_TRACER_SOURCE) then - bounded_below = .true.; bounded_above = .false. - lower_bound = 0.0_r8; upper_bound = missing_r8 -else - write(*, *) 'Illegal obs_kind in obs_error_info' - stop + + else + bounded_below = qcf_table_data(QTY_loc(1))%obs_error_info%bounded_below + bounded_above = qcf_table_data(QTY_loc(1))%obs_error_info%bounded_above + lower_bound = qcf_table_data(QTY_loc(1))%obs_error_info%lower_bound !NEED TO ADD CHECKS THAT THESE ARE VALID VALUES + upper_bound = qcf_table_data(QTY_loc(1))%obs_error_info%upper_bound + endif +write(*,*) 'obs_error_info: ', bounded_below, bounded_above, lower_bound, upper_bound + +! Set the observation error details for each type of quantity +!if(obs_kind == QTY_STATE_VARIABLE) then +! bounded_below = .false.; bounded_above = .false. +! lower_bound = missing_r8; upper_bound = missing_r8 +!elseif(obs_kind == QTY_TRACER_CONCENTRATION) then +! bounded_below = .true.; bounded_above = .false. +! lower_bound = 0.0_r8; upper_bound = missing_r8 +!elseif(obs_kind == QTY_TRACER_SOURCE) then +! bounded_below = .true.; bounded_above = .false. +! lower_bound = 0.0_r8; upper_bound = missing_r8 +!else +! write(*, *) 'Illegal obs_kind in obs_error_info' +! stop +!endif + end subroutine obs_error_info From 2376bbcebad5cd1f4d21da9b639d715a7e59be7f Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Wed, 16 Aug 2023 15:18:59 -0600 Subject: [PATCH 17/60] add subroutine to deallocate qcf table data structures --- .../modules/assimilation/algorithm_info_mod.f90 | 12 +++++++++++- 1 file changed, 11 insertions(+), 1 deletion(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index da3260042e..ca8120ec0b 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -36,7 +36,7 @@ module algorithm_info_mod integer, parameter :: BOUNDED_NORMAL_RHF = 101 public :: obs_error_info, probit_dist_info, obs_inc_info, & - init_qcf_table, read_qcf_table, & + init_qcf_table, deallocate_qcf_table, & obs_error_info_type, probit_inflation_type, probit_state_type, & probit_extended_state_type, obs_inc_info_type, qcf_table_data_type, & EAKF, ENKF, BOUNDED_NORMAL_RHF, UNBOUNDED_RHF, GAMMA_FILTER @@ -510,4 +510,14 @@ end subroutine write_qcf_table !------------------------------------------------------------------------ + +subroutine deallocate_qcf_table() + +deallocate(qcf_table_data) +deallocate(qcf_table_row_headers) + +end subroutine deallocate_qcf_table + +!---------------------------------------------------------------------- + end module algorithm_info_mod From 90c9bebb5fd533c1c92264a3485f909150243949 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Thu, 17 Aug 2023 11:25:34 -0600 Subject: [PATCH 18/60] making dealloc subroutine available to assim_tools_mod --- assimilation_code/modules/assimilation/assim_tools_mod.f90 | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/assimilation_code/modules/assimilation/assim_tools_mod.f90 b/assimilation_code/modules/assimilation/assim_tools_mod.f90 index 608fe799a3..b4f714eb93 100644 --- a/assimilation_code/modules/assimilation/assim_tools_mod.f90 +++ b/assimilation_code/modules/assimilation/assim_tools_mod.f90 @@ -76,7 +76,7 @@ module assim_tools_mod use normal_distribution_mod, only : normal_cdf, inv_weighted_normal_cdf -use algorithm_info_mod, only : probit_dist_info, obs_inc_info, init_qcf_table +use algorithm_info_mod, only : probit_dist_info, obs_inc_info, init_qcf_table, deallocate_qcf_table use gamma_distribution_mod, only : gamma_cdf, inv_gamma_cdf, gamma_mn_var_to_shape_scale, & gamma_gamma_prod @@ -314,7 +314,6 @@ subroutine assim_tools_init() call log_namelist_selections(num_special_cutoff, cache_override) -write(*,*), "HERE!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!" if(qcf_table_filename == '') then write(*,*), "no qcf table in namelist" else From 6de99ba08ad6079784d2016a4b34cccb7fdcbcfd Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Thu, 17 Aug 2023 13:28:52 -0600 Subject: [PATCH 19/60] removing comment blocks of old code --- .../assimilation/algorithm_info_mod.f90 | 95 +------------------ 1 file changed, 5 insertions(+), 90 deletions(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index ca8120ec0b..3eee46e986 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -95,6 +95,8 @@ module algorithm_info_mod contains !------------------------------------------------------------------------- + + subroutine obs_error_info(obs_def, error_variance, & bounded_below, bounded_above, lower_bound, upper_bound) @@ -134,7 +136,7 @@ subroutine obs_error_info(obs_def, error_variance, & write(*,*) 'findloc of kind: ', QTY_loc(1) if (QTY_loc(1) == 0) then - write(*,*) 'QTY not in table!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!' + write(*,*) 'QTY not in table, using default values' !use default values bounded_below = .false.; bounded_above = .false. @@ -150,24 +152,8 @@ subroutine obs_error_info(obs_def, error_variance, & write(*,*) 'obs_error_info: ', bounded_below, bounded_above, lower_bound, upper_bound -! Set the observation error details for each type of quantity -!if(obs_kind == QTY_STATE_VARIABLE) then -! bounded_below = .false.; bounded_above = .false. -! lower_bound = missing_r8; upper_bound = missing_r8 -!elseif(obs_kind == QTY_TRACER_CONCENTRATION) then -! bounded_below = .true.; bounded_above = .false. -! lower_bound = 0.0_r8; upper_bound = missing_r8 -!elseif(obs_kind == QTY_TRACER_SOURCE) then -! bounded_below = .true.; bounded_above = .false. -! lower_bound = 0.0_r8; upper_bound = missing_r8 -!else -! write(*, *) 'Illegal obs_kind in obs_error_info' -! stop -!endif - end subroutine obs_error_info - !------------------------------------------------------------------------- @@ -215,7 +201,7 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & write(*,*) 'findloc of kind: ', QTY_loc(1) if (QTY_loc(1) == 0) then - write(*,*) 'QTY not in table!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!' + write(*,*) 'QTY not in table, using default values' !use default values dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION @@ -231,23 +217,6 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & lower_bound = qcf_table_data(QTY_loc(1))%probit_inflation%lower_bound !NEED TO ADD CHECKS THAT THESE ARE VALID VALUES upper_bound = qcf_table_data(QTY_loc(1))%probit_inflation%upper_bound -! if(kind == QTY_STATE_VARIABLE) then -! dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION -! bounded_below = .false.; bounded_above = .false. -! lower_bound = missing_r8; upper_bound = missing_r8 -! elseif(kind == QTY_TRACER_CONCENTRATION) then -! dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION -! bounded_below = .true.; bounded_above = .false. -! lower_bound = 0.0_r8; upper_bound = missing_r8 -! elseif(kind == QTY_TRACER_SOURCE) then -! dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION -! bounded_below = .true.; bounded_above = .false. -! lower_bound = 0.0_r8; upper_bound = missing_r8 -! else -! write(*, *) 'Illegal kind in obs_error_info' -! stop -! endif - elseif(is_state) then ! Case for state variable priors @@ -257,23 +226,6 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & lower_bound = qcf_table_data(QTY_loc(1))%probit_state%lower_bound upper_bound = qcf_table_data(QTY_loc(1))%probit_state%upper_bound - !if(kind == QTY_STATE_VARIABLE) then - ! dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION - ! bounded_below = .false.; bounded_above = .false. - ! lower_bound = missing_r8; upper_bound = missing_r8 - ! elseif(kind == QTY_TRACER_CONCENTRATION) then - ! dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION - ! bounded_below = .true.; bounded_above = .false. - ! lower_bound = 0.0_r8; upper_bound = missing_r8 -! elseif(kind == QTY_TRACER_SOURCE) then - ! dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION - ! bounded_below = .true.; bounded_above = .false. - ! lower_bound = 0.0_r8; upper_bound = missing_r8 -! else -! write(*, *) 'Illegal kind in obs_error_info' - ! stop - ! endif - else ! This case is for observation (extended state) priors @@ -283,22 +235,6 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & lower_bound = qcf_table_data(QTY_loc(1))%probit_extended_state%lower_bound upper_bound = qcf_table_data(QTY_loc(1))%probit_extended_state%upper_bound -! if(kind == QTY_STATE_VARIABLE) then -! dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION -! bounded_below = .false.; bounded_above = .false. -! lower_bound = missing_r8; upper_bound = missing_r8 -! elseif(kind == QTY_TRACER_CONCENTRATION) then -! dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION -! bounded_below = .true.; bounded_above = .false. -! lower_bound = 0.0_r8; upper_bound = missing_r8 -! elseif(kind == QTY_TRACER_SOURCE) then -! dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION -! bounded_below = .true.; bounded_above = .false. -! lower_bound = 0.0_r8; upper_bound = missing_r8 -! else -! write(*, *) 'Illegal kind in obs_error_info' -! stop -! endif endif write(*,*) 'probit_dist_info: ', dist_type, bounded_below, bounded_above, lower_bound, upper_bound @@ -338,7 +274,7 @@ subroutine obs_inc_info(obs_kind, filter_kind, rectangular_quadrature, gaussian_ write(*,*) 'findloc of kind: ', QTY_loc(1) if (QTY_loc(1) == 0) then - write(*,*) 'QTY not in table!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!' + write(*,*) 'QTY not in table, using default values' !use default values filter_kind = BOUNDED_NORMAL_RHF @@ -360,24 +296,6 @@ subroutine obs_inc_info(obs_kind, filter_kind, rectangular_quadrature, gaussian_ write(*,*) 'obs_inc_info: ', filter_kind, sort_obs_inc, spread_restoration, bounded_below, bounded_above, lower_bound, upper_bound -! Set the observation increment details for each type of quantity -!if(obs_kind == QTY_STATE_VARIABLE) then -! filter_kind = BOUNDED_NORMAL_RHF -! bounded_below = .false.; bounded_above = .false. -! lower_bound = missing_r8; upper_bound = missing_r8 -!elseif(obs_kind == QTY_TRACER_CONCENTRATION) then -! filter_kind = BOUNDED_NORMAL_RHF -! bounded_below = .true.; bounded_above = .false. -! lower_bound = 0.0_r8; upper_bound = missing_r8 -!elseif(obs_kind == QTY_TRACER_SOURCE) then -! filter_kind = BOUNDED_NORMAL_RHF -! bounded_below = .true.; bounded_above = .false. -! lower_bound = 0.0_r8; upper_bound = missing_r8 -!else -! write(*, *) 'Illegal obs_kind in obs_error_info' -! stop -!endif - ! Only need to set these two for options the original RHF implementation !!!rectangular_quadrature = .true. !!!gaussian_likelihood_tails = .false. @@ -425,9 +343,6 @@ subroutine read_qcf_table(qcf_table_filename) character(len=129), intent(in) :: qcf_table_filename -!type(qcf_table_data_type) :: qcf_table_data(:) -!character(len=129) :: rowheaders(:) !!!!! might need to change len=129 - integer, parameter :: fileid = 10 !file identifier integer :: row From c5e0eb15341e4643a26f7321b3a11444a5b8be75 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Mon, 21 Aug 2023 15:43:57 -0600 Subject: [PATCH 20/60] Adding call to deallocate routine, removing unused var and old commented code --- .../modules/assimilation/algorithm_info_mod.f90 | 6 ------ assimilation_code/modules/assimilation/assim_tools_mod.f90 | 7 +++++-- 2 files changed, 5 insertions(+), 8 deletions(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index 3eee46e986..d8e93391f0 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -9,12 +9,6 @@ module algorithm_info_mod use obs_def_mod, only : obs_def_type, get_obs_def_type_of_obs, get_obs_def_error_variance use obs_kind_mod, only : get_quantity_for_type_of_obs, get_name_for_quantity - -! Get the QTY definitions that are needed (aka kind) -use obs_kind_mod, only : QTY_STATE_VARIABLE, QTY_STATE_VAR_POWER, QTY_TRACER_CONCENTRATION, & - QTY_TRACER_SOURCE -! NOTE: Sadly, the QTY itself is not sufficient for the POWER because there is additional metadata - use assim_model_mod, only : get_state_meta_data use location_mod, only : location_type diff --git a/assimilation_code/modules/assimilation/assim_tools_mod.f90 b/assimilation_code/modules/assimilation/assim_tools_mod.f90 index b4f714eb93..d4bfa02b6a 100644 --- a/assimilation_code/modules/assimilation/assim_tools_mod.f90 +++ b/assimilation_code/modules/assimilation/assim_tools_mod.f90 @@ -223,7 +223,7 @@ module assim_tools_mod subroutine assim_tools_init() -integer :: iunit, io, i, j, numrows +integer :: iunit, io, i, j integer :: num_special_cutoff, type_index logical :: cache_override = .false. @@ -315,7 +315,7 @@ subroutine assim_tools_init() call log_namelist_selections(num_special_cutoff, cache_override) if(qcf_table_filename == '') then - write(*,*), "no qcf table in namelist" + write(*,*), "No QCF table in namelist, using default values for all QTYs" else call init_qcf_table(qcf_table_filename) endif @@ -903,6 +903,9 @@ subroutine filter_assim(ens_handle, obs_ens_handle, obs_seq, keys, & ! get rid of mpi window call free_mean_window() +! free qcf_table_data structures +call deallocate_qcf_table() + ! deallocate space deallocate(close_obs_dist, & my_obs_indx, & From d61382ed2fba63fb8a580637280db797c63250ef Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Mon, 21 Aug 2023 16:47:06 -0600 Subject: [PATCH 21/60] Fixing typo in subroutine names --- .../assimilation/algorithm_info_mod.f90 | 18 +++++++++--------- .../modules/assimilation/assim_tools_mod.f90 | 10 +++++----- 2 files changed, 14 insertions(+), 14 deletions(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index d8e93391f0..b5fa43a560 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -30,7 +30,7 @@ module algorithm_info_mod integer, parameter :: BOUNDED_NORMAL_RHF = 101 public :: obs_error_info, probit_dist_info, obs_inc_info, & - init_qcf_table, deallocate_qcf_table, & + init_algorithm_info_mod, end_algorithm_info_mod, & obs_error_info_type, probit_inflation_type, probit_state_type, & probit_extended_state_type, obs_inc_info_type, qcf_table_data_type, & EAKF, ENKF, BOUNDED_NORMAL_RHF, UNBOUNDED_RHF, GAMMA_FILTER @@ -76,7 +76,7 @@ module algorithm_info_mod end type type(qcf_table_data_type), allocatable :: qcf_table_data(:) -character(len=129), allocatable :: qcf_table_row_headers(:) !!!!! might need to change len=129 +character(len=129), allocatable :: qcf_table_row_headers(:) ! Provides routines that give information about details of algorithms for ! observation error sampling, observation increments, and the transformations @@ -139,7 +139,7 @@ subroutine obs_error_info(obs_def, error_variance, & else bounded_below = qcf_table_data(QTY_loc(1))%obs_error_info%bounded_below bounded_above = qcf_table_data(QTY_loc(1))%obs_error_info%bounded_above - lower_bound = qcf_table_data(QTY_loc(1))%obs_error_info%lower_bound !NEED TO ADD CHECKS THAT THESE ARE VALID VALUES + lower_bound = qcf_table_data(QTY_loc(1))%obs_error_info%lower_bound upper_bound = qcf_table_data(QTY_loc(1))%obs_error_info%upper_bound endif @@ -208,7 +208,7 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & dist_type = qcf_table_data(QTY_loc(1))%probit_inflation%dist_type !dist_type has checks in transform_to_probit, transform_from_probit bounded_below = qcf_table_data(QTY_loc(1))%probit_inflation%bounded_below bounded_above = qcf_table_data(QTY_loc(1))%probit_inflation%bounded_above - lower_bound = qcf_table_data(QTY_loc(1))%probit_inflation%lower_bound !NEED TO ADD CHECKS THAT THESE ARE VALID VALUES + lower_bound = qcf_table_data(QTY_loc(1))%probit_inflation%lower_bound upper_bound = qcf_table_data(QTY_loc(1))%probit_inflation%upper_bound elseif(is_state) then @@ -283,7 +283,7 @@ subroutine obs_inc_info(obs_kind, filter_kind, rectangular_quadrature, gaussian_ spread_restoration = qcf_table_data(QTY_loc(1))%obs_inc_info%spread_restoration bounded_below = qcf_table_data(QTY_loc(1))%obs_inc_info%bounded_below bounded_above = qcf_table_data(QTY_loc(1))%obs_inc_info%bounded_above - lower_bound = qcf_table_data(QTY_loc(1))%obs_inc_info%lower_bound !NEED TO ADD CHECKS THAT THESE ARE VALID VALUES + lower_bound = qcf_table_data(QTY_loc(1))%obs_inc_info%lower_bound upper_bound = qcf_table_data(QTY_loc(1))%obs_inc_info%upper_bound endif @@ -299,7 +299,7 @@ end subroutine obs_inc_info !------------------------------------------------------------------------ -subroutine init_qcf_table(qcf_table_filename) +subroutine init_algorithm_info_mod(qcf_table_filename) character(len=129), intent(in) :: qcf_table_filename @@ -326,7 +326,7 @@ subroutine init_qcf_table(qcf_table_filename) call read_qcf_table(qcf_table_filename) -end subroutine init_qcf_table +end subroutine init_algorithm_info_mod !------------------------------------------------------------------------ @@ -420,12 +420,12 @@ end subroutine write_qcf_table !------------------------------------------------------------------------ -subroutine deallocate_qcf_table() +subroutine end_algorithm_info_mod() deallocate(qcf_table_data) deallocate(qcf_table_row_headers) -end subroutine deallocate_qcf_table +end subroutine end_algorithm_info_mod !---------------------------------------------------------------------- diff --git a/assimilation_code/modules/assimilation/assim_tools_mod.f90 b/assimilation_code/modules/assimilation/assim_tools_mod.f90 index d4bfa02b6a..393627be56 100644 --- a/assimilation_code/modules/assimilation/assim_tools_mod.f90 +++ b/assimilation_code/modules/assimilation/assim_tools_mod.f90 @@ -76,7 +76,7 @@ module assim_tools_mod use normal_distribution_mod, only : normal_cdf, inv_weighted_normal_cdf -use algorithm_info_mod, only : probit_dist_info, obs_inc_info, init_qcf_table, deallocate_qcf_table +use algorithm_info_mod, only : probit_dist_info, obs_inc_info, init_algorithm_info_mod, end_algorithm_info_mod use gamma_distribution_mod, only : gamma_cdf, inv_gamma_cdf, gamma_mn_var_to_shape_scale, & gamma_gamma_prod @@ -143,7 +143,7 @@ module assim_tools_mod ! special_localization_obs_types -> Special treatment for the specified observation types ! special_localization_cutoffs -> Different cutoff value for each specified obs type ! -character(len = 129) :: qcf_table_filename = '' !not sure if the len should be 129 here, but it is consistent with other nml variables +character(len = 129) :: qcf_table_filename = '' logical :: use_algorithm_info_mod = .true. integer :: filter_kind = 1 real(r8) :: cutoff = 0.2_r8 @@ -317,7 +317,7 @@ subroutine assim_tools_init() if(qcf_table_filename == '') then write(*,*), "No QCF table in namelist, using default values for all QTYs" else - call init_qcf_table(qcf_table_filename) + call init_algorithm_info_mod(qcf_table_filename) endif end subroutine assim_tools_init @@ -903,8 +903,8 @@ subroutine filter_assim(ens_handle, obs_ens_handle, obs_seq, keys, & ! get rid of mpi window call free_mean_window() -! free qcf_table_data structures -call deallocate_qcf_table() +! deallocate qcf_table_data structures +call end_algorithm_info_mod() ! deallocate space deallocate(close_obs_dist, & From e842a2b1a9b2790922ccbb1fd83f1b58f490c660 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Tue, 22 Aug 2023 13:58:25 -0600 Subject: [PATCH 22/60] Moving the allocation and deallocation of qcf table data from assim_tools_mod to filter_main in filter_mod --- .../modules/assimilation/assim_tools_mod.f90 | 14 ++------------ .../modules/assimilation/filter_mod.f90 | 14 +++++++++++++- 2 files changed, 15 insertions(+), 13 deletions(-) diff --git a/assimilation_code/modules/assimilation/assim_tools_mod.f90 b/assimilation_code/modules/assimilation/assim_tools_mod.f90 index 393627be56..958e373201 100644 --- a/assimilation_code/modules/assimilation/assim_tools_mod.f90 +++ b/assimilation_code/modules/assimilation/assim_tools_mod.f90 @@ -76,7 +76,7 @@ module assim_tools_mod use normal_distribution_mod, only : normal_cdf, inv_weighted_normal_cdf -use algorithm_info_mod, only : probit_dist_info, obs_inc_info, init_algorithm_info_mod, end_algorithm_info_mod +use algorithm_info_mod, only : probit_dist_info, obs_inc_info use gamma_distribution_mod, only : gamma_cdf, inv_gamma_cdf, gamma_mn_var_to_shape_scale, & gamma_gamma_prod @@ -143,7 +143,6 @@ module assim_tools_mod ! special_localization_obs_types -> Special treatment for the specified observation types ! special_localization_cutoffs -> Different cutoff value for each specified obs type ! -character(len = 129) :: qcf_table_filename = '' logical :: use_algorithm_info_mod = .true. integer :: filter_kind = 1 real(r8) :: cutoff = 0.2_r8 @@ -204,7 +203,7 @@ module assim_tools_mod ! compared to previous versions of this namelist item. logical :: distribute_mean = .false. -namelist / assim_tools_nml / qcf_table_filename, use_algorithm_info_mod, & +namelist / assim_tools_nml / use_algorithm_info_mod, & filter_kind, cutoff, sort_obs_inc, & spread_restoration, sampling_error_correction, & adaptive_localization_threshold, adaptive_cutoff_floor, & @@ -314,12 +313,6 @@ subroutine assim_tools_init() call log_namelist_selections(num_special_cutoff, cache_override) -if(qcf_table_filename == '') then - write(*,*), "No QCF table in namelist, using default values for all QTYs" -else - call init_algorithm_info_mod(qcf_table_filename) -endif - end subroutine assim_tools_init !------------------------------------------------------------- @@ -903,9 +896,6 @@ subroutine filter_assim(ens_handle, obs_ens_handle, obs_seq, keys, & ! get rid of mpi window call free_mean_window() -! deallocate qcf_table_data structures -call end_algorithm_info_mod() - ! deallocate space deallocate(close_obs_dist, & my_obs_indx, & diff --git a/assimilation_code/modules/assimilation/filter_mod.f90 b/assimilation_code/modules/assimilation/filter_mod.f90 index 11e50c3038..0c4035b9a6 100644 --- a/assimilation_code/modules/assimilation/filter_mod.f90 +++ b/assimilation_code/modules/assimilation/filter_mod.f90 @@ -90,7 +90,7 @@ module filter_mod use probit_transform_mod, only : transform_to_probit, transform_from_probit -use algorithm_info_mod, only : probit_dist_info +use algorithm_info_mod, only : probit_dist_info, init_algorithm_info_mod, end_algorithm_info_mod use distribution_params_mod, only : distribution_params_type, NORMAL_DISTRIBUTION @@ -166,6 +166,7 @@ module filter_mod !---------------------------------------------------------------- ! Namelist input with default values ! +character(len = 129) :: qcf_table_filename = '' logical :: use_algorithm_info_mod = .true. integer :: async = 0, ens_size = 20 integer :: tasks_per_model_advance = 1 @@ -261,6 +262,7 @@ module filter_mod namelist /filter_nml/ async, & + qcf_table_filename, & use_algorithm_info_mod, & adv_ens_command, & ens_size, & @@ -1150,6 +1152,9 @@ subroutine filter_main() call end_assim_model() call trace_message('After end_model call') +! deallocate qcf_table_data structures +!call end_algorithm_info_mod() + call trace_message('Before ensemble and obs memory cleanup') call end_ensemble_manager(state_ens_handle) @@ -1268,6 +1273,13 @@ subroutine filter_initialize_modules_used() ! Initialize the obs sequence module call static_init_obs_sequence() +! Initialize algorothm_info_mod and read in QCF table data +if(qcf_table_filename == '') then + write(*,*) "No QCF table in namelist, using default values for all QTYs" +else + call init_algorithm_info_mod(qcf_table_filename) +endif + ! Initialize the model class data now that obs_sequence is all set up call static_init_assim_model() call state_vector_io_init() From a78be832e9c65764eef5a7eb7b5cb660cbdd57d1 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Tue, 22 Aug 2023 14:38:47 -0600 Subject: [PATCH 23/60] uncommenting call to end_alg_info_mod --- assimilation_code/modules/assimilation/filter_mod.f90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/assimilation_code/modules/assimilation/filter_mod.f90 b/assimilation_code/modules/assimilation/filter_mod.f90 index 0c4035b9a6..d967eb71e6 100644 --- a/assimilation_code/modules/assimilation/filter_mod.f90 +++ b/assimilation_code/modules/assimilation/filter_mod.f90 @@ -1153,7 +1153,7 @@ subroutine filter_main() call trace_message('After end_model call') ! deallocate qcf_table_data structures -!call end_algorithm_info_mod() +call end_algorithm_info_mod() call trace_message('Before ensemble and obs memory cleanup') call end_ensemble_manager(state_ens_handle) From dc0c617431e7e848ccca63981182bfd2a2f7bba1 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Wed, 23 Aug 2023 15:28:06 -0600 Subject: [PATCH 24/60] moving call to init_algortihm_info_mod out of conditional --- .../modules/assimilation/filter_mod.f90 | 13 +++++-------- 1 file changed, 5 insertions(+), 8 deletions(-) diff --git a/assimilation_code/modules/assimilation/filter_mod.f90 b/assimilation_code/modules/assimilation/filter_mod.f90 index d967eb71e6..874a4cb1d5 100644 --- a/assimilation_code/modules/assimilation/filter_mod.f90 +++ b/assimilation_code/modules/assimilation/filter_mod.f90 @@ -166,7 +166,7 @@ module filter_mod !---------------------------------------------------------------- ! Namelist input with default values ! -character(len = 129) :: qcf_table_filename = '' +character(len = 129) :: qcf_table_filename = 'real_qcf_table.txt' !NEED TO REMOVE THIS LATER logical :: use_algorithm_info_mod = .true. integer :: async = 0, ens_size = 20 integer :: tasks_per_model_advance = 1 @@ -1273,17 +1273,14 @@ subroutine filter_initialize_modules_used() ! Initialize the obs sequence module call static_init_obs_sequence() -! Initialize algorothm_info_mod and read in QCF table data -if(qcf_table_filename == '') then - write(*,*) "No QCF table in namelist, using default values for all QTYs" -else - call init_algorithm_info_mod(qcf_table_filename) -endif - ! Initialize the model class data now that obs_sequence is all set up call static_init_assim_model() call state_vector_io_init() call initialize_qc() + +! Initialize algorothm_info_mod and read in QCF table data +call init_algorithm_info_mod(qcf_table_filename) + call trace_message('After filter_initialize_module_used call') end subroutine filter_initialize_modules_used From 2d6808e40f8355a0a11bb602cdde59d203d62dc4 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Wed, 23 Aug 2023 15:55:22 -0600 Subject: [PATCH 25/60] Reorganizing the subroutines so that init_algorithm_info_mod is at the top of algorithm_info_mod --- .../assimilation/algorithm_info_mod.f90 | 158 +++++++++--------- 1 file changed, 81 insertions(+), 77 deletions(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index b5fa43a560..027cf747eb 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -91,6 +91,87 @@ module algorithm_info_mod !------------------------------------------------------------------------- +subroutine init_algorithm_info_mod(qcf_table_filename) + +! Gets number of lines/QTYs in the QCF table, allocates space for the table data + +character(len=129), intent(in) :: qcf_table_filename + +integer :: numrows +integer :: nlines +integer :: io +integer, parameter :: fileid = 10 !file identifier + +write(*,*) 'filename: ', qcf_table_filename + +open(unit=fileid, file=qcf_table_filename) +nlines = 0 + +do !do loop to get number of rows (or QTY's) in the table + read(fileid,*,iostat=io) + if(io/=0) exit + nlines = nlines + 1 +end do +close(fileid) + +numrows = nlines - 2 +print *, 'numrows: ', numrows + +allocate(qcf_table_data(numrows)) +allocate(qcf_table_row_headers(numrows)) + +call read_qcf_table(qcf_table_filename) +call write_qcf_table() + +end subroutine init_algorithm_info_mod + +!------------------------------------------------------------------------ + + +subroutine read_qcf_table(qcf_table_filename) + +! Reads in the QCEFF input options from tabular data file + +character(len=129), intent(in) :: qcf_table_filename + +integer, parameter :: fileid = 10 !file identifier +integer :: row + +character(len=129), dimension(4) :: header1 +character(len=129), dimension(29) :: header2 + +write(*,*) 'filename: ', qcf_table_filename +open(unit=fileid, file=qcf_table_filename) + +! skip the headers +read(fileid, *) header1 +read(fileid, *) header2 +write(*, *) "header1: ", header1 +write(*, *) "header2: ", header2 + +! read in table values directly to qcf_table_data type +do row = 1, size(qcf_table_data) + read(fileid, *) qcf_table_row_headers(row), qcf_table_data(row)%obs_error_info%bounded_below, qcf_table_data(row)%obs_error_info%bounded_above, & + qcf_table_data(row)%obs_error_info%lower_bound, qcf_table_data(row)%obs_error_info%upper_bound, qcf_table_data(row)%probit_inflation%dist_type, & + qcf_table_data(row)%probit_inflation%bounded_below, qcf_table_data(row)%probit_inflation%bounded_above, & + qcf_table_data(row)%probit_inflation%lower_bound, qcf_table_data(row)%probit_inflation%upper_bound, qcf_table_data(row)%probit_state%dist_type, & + qcf_table_data(row)%probit_state%bounded_below, qcf_table_data(row)%probit_state%bounded_above, & + qcf_table_data(row)%probit_state%lower_bound, qcf_table_data(row)%probit_state%upper_bound, qcf_table_data(row)%probit_extended_state%dist_type, & + qcf_table_data(row)%probit_extended_state%bounded_below, qcf_table_data(row)%probit_extended_state%bounded_above, & + qcf_table_data(row)%probit_extended_state%lower_bound, qcf_table_data(row)%probit_extended_state%upper_bound, & + qcf_table_data(row)%obs_inc_info%filter_kind, qcf_table_data(row)%obs_inc_info%rectangular_quadrature, & + qcf_table_data(row)%obs_inc_info%gaussian_likelihood_tails, qcf_table_data(row)%obs_inc_info%sort_obs_inc, & + qcf_table_data(row)%obs_inc_info%spread_restoration, qcf_table_data(row)%obs_inc_info%bounded_below, qcf_table_data(row)%obs_inc_info%bounded_above, & + qcf_table_data(row)%obs_inc_info%lower_bound, qcf_table_data(row)%obs_inc_info%upper_bound +end do + +close(fileid) + +end subroutine read_qcf_table + +!------------------------------------------------------------------------ + + subroutine obs_error_info(obs_def, error_variance, & bounded_below, bounded_above, lower_bound, upper_bound) @@ -299,83 +380,6 @@ end subroutine obs_inc_info !------------------------------------------------------------------------ -subroutine init_algorithm_info_mod(qcf_table_filename) - -character(len=129), intent(in) :: qcf_table_filename - -integer :: numrows -integer :: nlines -integer :: io -integer, parameter :: fileid = 10 !file identifier - -open(unit=fileid, file=qcf_table_filename) -nlines = 0 - -do !do loop to get number of rows (or QTY's) in the table - read(fileid,*,iostat=io) - if(io/=0) exit - nlines = nlines + 1 -end do -close(fileid) - -numrows = nlines - 2 -print *, 'numrows: ', numrows - -allocate(qcf_table_data(numrows)) -allocate(qcf_table_row_headers(numrows)) - -call read_qcf_table(qcf_table_filename) - -end subroutine init_algorithm_info_mod - -!------------------------------------------------------------------------ - - -subroutine read_qcf_table(qcf_table_filename) - -! Reads in the QCEFF input options from tabular data file - -character(len=129), intent(in) :: qcf_table_filename - -integer, parameter :: fileid = 10 !file identifier -integer :: row - -character(len=129), dimension(4) :: header1 -character(len=129), dimension(29) :: header2 - -open(unit=fileid, file=qcf_table_filename) - -! skip the headers -read(fileid, *) header1 -read(fileid, *) header2 -write(*, *) "header1: ", header1 -write(*, *) "header2: ", header2 - -! read in table values directly to qcf_table_data type -do row = 1, size(qcf_table_data) - read(fileid, *) qcf_table_row_headers(row), qcf_table_data(row)%obs_error_info%bounded_below, qcf_table_data(row)%obs_error_info%bounded_above, & - qcf_table_data(row)%obs_error_info%lower_bound, qcf_table_data(row)%obs_error_info%upper_bound, qcf_table_data(row)%probit_inflation%dist_type, & - qcf_table_data(row)%probit_inflation%bounded_below, qcf_table_data(row)%probit_inflation%bounded_above, & - qcf_table_data(row)%probit_inflation%lower_bound, qcf_table_data(row)%probit_inflation%upper_bound, qcf_table_data(row)%probit_state%dist_type, & - qcf_table_data(row)%probit_state%bounded_below, qcf_table_data(row)%probit_state%bounded_above, & - qcf_table_data(row)%probit_state%lower_bound, qcf_table_data(row)%probit_state%upper_bound, qcf_table_data(row)%probit_extended_state%dist_type, & - qcf_table_data(row)%probit_extended_state%bounded_below, qcf_table_data(row)%probit_extended_state%bounded_above, & - qcf_table_data(row)%probit_extended_state%lower_bound, qcf_table_data(row)%probit_extended_state%upper_bound, & - qcf_table_data(row)%obs_inc_info%filter_kind, qcf_table_data(row)%obs_inc_info%rectangular_quadrature, & - qcf_table_data(row)%obs_inc_info%gaussian_likelihood_tails, qcf_table_data(row)%obs_inc_info%sort_obs_inc, & - qcf_table_data(row)%obs_inc_info%spread_restoration, qcf_table_data(row)%obs_inc_info%bounded_below, qcf_table_data(row)%obs_inc_info%bounded_above, & - qcf_table_data(row)%obs_inc_info%lower_bound, qcf_table_data(row)%obs_inc_info%upper_bound -end do - -close(fileid) - -call write_qcf_table() - -end subroutine read_qcf_table - -!------------------------------------------------------------------------ - - subroutine write_qcf_table() ! DRAFT SUBROUTINE From 19b6150d69b887dc5ad7bd5269a68e56db4f7eba Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Wed, 23 Aug 2023 16:32:46 -0600 Subject: [PATCH 26/60] Adding qcf_table_listed logical and module_initialized checks --- .../assimilation/algorithm_info_mod.f90 | 23 +++++++++++++++++-- 1 file changed, 21 insertions(+), 2 deletions(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index 027cf747eb..c5b6ef965e 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -19,6 +19,9 @@ module algorithm_info_mod implicit none private +logical :: module_initialized = .false. +logical :: qcf_table_listed = .false. + ! Defining parameter strings for different observation space filters ! For now, retaining backwards compatibility in assim_tools_mod requires using ! these specific integer values and there is no point in using these in assim_tools. @@ -102,8 +105,17 @@ subroutine init_algorithm_info_mod(qcf_table_filename) integer :: io integer, parameter :: fileid = 10 !file identifier +if (module_initialized) return +module_initialized = .true. + write(*,*) 'filename: ', qcf_table_filename +if (qcf_table_filename == '') then + write(*,*) 'No QCF table file listed in namelist, using default values for all QTYs' + return +endif + +qcf_table_listed = .true. open(unit=fileid, file=qcf_table_filename) nlines = 0 @@ -140,6 +152,8 @@ subroutine read_qcf_table(qcf_table_filename) character(len=129), dimension(4) :: header1 character(len=129), dimension(29) :: header2 +if (.not. module_initialized) call init_algorithm_info_mod(qcf_table_filename) + write(*,*) 'filename: ', qcf_table_filename open(unit=fileid, file=qcf_table_filename) @@ -210,9 +224,12 @@ subroutine obs_error_info(obs_def, error_variance, & QTY_loc = findloc(qcf_table_row_headers, kind_name) write(*,*) 'findloc of kind: ', QTY_loc(1) -if (QTY_loc(1) == 0) then - write(*,*) 'QTY not in table, using default values' +!use default values if qcf_table_filename is not in namelist +if (.not. qcf_table_listed) then + bounded_below = .false.; bounded_above = .false. + lower_bound = missing_r8; upper_bound = missing_r8 +if (QTY_loc(1) == 0) then !use default values bounded_below = .false.; bounded_above = .false. lower_bound = missing_r8; upper_bound = missing_r8 @@ -429,6 +446,8 @@ subroutine end_algorithm_info_mod() deallocate(qcf_table_data) deallocate(qcf_table_row_headers) +module_initialized = .false. + end subroutine end_algorithm_info_mod !---------------------------------------------------------------------- From efa929bae005eabf017fd75f2aabb98aa788a817 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Wed, 23 Aug 2023 16:35:31 -0600 Subject: [PATCH 27/60] Moving location of qcf_table_listed check to before data access from findloc --- .../modules/assimilation/algorithm_info_mod.f90 | 10 ++++++---- 1 file changed, 6 insertions(+), 4 deletions(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index c5b6ef965e..67a68d28e8 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -220,14 +220,16 @@ subroutine obs_error_info(obs_def, error_variance, & kind_name = get_name_for_quantity(obs_kind) write(*,*) 'kind_name: ', kind_name -!find location of QTY in qcf_table_data structure -QTY_loc = findloc(qcf_table_row_headers, kind_name) -write(*,*) 'findloc of kind: ', QTY_loc(1) - !use default values if qcf_table_filename is not in namelist if (.not. qcf_table_listed) then bounded_below = .false.; bounded_above = .false. lower_bound = missing_r8; upper_bound = missing_r8 + return +endif + +!find location of QTY in qcf_table_data structure +QTY_loc = findloc(qcf_table_row_headers, kind_name) +write(*,*) 'findloc of kind: ', QTY_loc(1) if (QTY_loc(1) == 0) then !use default values From 16354ab12881cda23deccc04499f203c9950d801 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Mon, 28 Aug 2023 16:17:50 -0600 Subject: [PATCH 28/60] Using error_handler from utilities_mod; adding check for correct table version --- .../assimilation/algorithm_info_mod.f90 | 44 ++++++++++++------- 1 file changed, 28 insertions(+), 16 deletions(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index 67a68d28e8..e3afcfbda6 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -9,6 +9,8 @@ module algorithm_info_mod use obs_def_mod, only : obs_def_type, get_obs_def_type_of_obs, get_obs_def_error_variance use obs_kind_mod, only : get_quantity_for_type_of_obs, get_name_for_quantity +use utilities_mod, only : error_handler, E_ERR + use assim_model_mod, only : get_state_meta_data use location_mod, only : location_type @@ -19,6 +21,9 @@ module algorithm_info_mod implicit none private +character(len=512) :: errstring +character(len=*), parameter :: source = 'algorithm_info_mod.f90' + logical :: module_initialized = .false. logical :: qcf_table_listed = .false. @@ -154,14 +159,11 @@ subroutine read_qcf_table(qcf_table_filename) if (.not. module_initialized) call init_algorithm_info_mod(qcf_table_filename) -write(*,*) 'filename: ', qcf_table_filename open(unit=fileid, file=qcf_table_filename) -! skip the headers +! skip the headers, make sure user is using the correct table version read(fileid, *) header1 read(fileid, *) header2 -write(*, *) "header1: ", header1 -write(*, *) "header2: ", header2 ! read in table values directly to qcf_table_data type do row = 1, size(qcf_table_data) @@ -181,6 +183,8 @@ subroutine read_qcf_table(qcf_table_filename) close(fileid) +call assert_qcf_table_version(header1) + end subroutine read_qcf_table !------------------------------------------------------------------------ @@ -218,7 +222,6 @@ subroutine obs_error_info(obs_def, error_variance, & !get actual name of QTY from integer index kind_name = get_name_for_quantity(obs_kind) -write(*,*) 'kind_name: ', kind_name !use default values if qcf_table_filename is not in namelist if (.not. qcf_table_listed) then @@ -229,7 +232,6 @@ subroutine obs_error_info(obs_def, error_variance, & !find location of QTY in qcf_table_data structure QTY_loc = findloc(qcf_table_row_headers, kind_name) -write(*,*) 'findloc of kind: ', QTY_loc(1) if (QTY_loc(1) == 0) then !use default values @@ -244,8 +246,6 @@ subroutine obs_error_info(obs_def, error_variance, & endif -write(*,*) 'obs_error_info: ', bounded_below, bounded_above, lower_bound, upper_bound - end subroutine obs_error_info !------------------------------------------------------------------------- @@ -288,11 +288,9 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & !get actual name of QTY from integer index kind_name = get_name_for_quantity(kind) -write(*,*) 'kind_name: ', kind_name !find location of QTY in qcf_table_data structure QTY_loc = findloc(qcf_table_row_headers, kind_name) -write(*,*) 'findloc of kind: ', QTY_loc(1) if (QTY_loc(1) == 0) then write(*,*) 'QTY not in table, using default values' @@ -331,8 +329,6 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & endif -write(*,*) 'probit_dist_info: ', dist_type, bounded_below, bounded_above, lower_bound, upper_bound - end subroutine probit_dist_info !------------------------------------------------------------------------ @@ -361,11 +357,9 @@ subroutine obs_inc_info(obs_kind, filter_kind, rectangular_quadrature, gaussian_ !get actual name of QTY from integer index kind_name = get_name_for_quantity(obs_kind) -write(*,*) 'kind_name: ', kind_name !find location of QTY in qcf_table_data structure QTY_loc = findloc(qcf_table_row_headers, kind_name) -write(*,*) 'findloc of kind: ', QTY_loc(1) if (QTY_loc(1) == 0) then write(*,*) 'QTY not in table, using default values' @@ -388,8 +382,6 @@ subroutine obs_inc_info(obs_kind, filter_kind, rectangular_quadrature, gaussian_ endif -write(*,*) 'obs_inc_info: ', filter_kind, sort_obs_inc, spread_restoration, bounded_below, bounded_above, lower_bound, upper_bound - ! Only need to set these two for options the original RHF implementation !!!rectangular_quadrature = .true. !!!gaussian_likelihood_tails = .false. @@ -443,8 +435,28 @@ end subroutine write_qcf_table !------------------------------------------------------------------------ +subroutine assert_qcf_table_version(header) + +!subroutine to ensure the correct version of the QCF table is being used + +character(len=129), dimension(4), intent(in) :: header + +write(*,*) 'version: ', header(4) + +if (header(4) /= '1:') then + write(errstring,*) "Using outdated/incorrect version of the QCF table" + call error_handler(E_ERR, 'assert_qcf_table_version', errstring, source) +endif + +end subroutine assert_qcf_table_version + +!------------------------------------------------------------------------ + + subroutine end_algorithm_info_mod() +if (.not. module_initialized) return + deallocate(qcf_table_data) deallocate(qcf_table_row_headers) From 85955230a7556d0e8928b13ac29d5d102c7b633f Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Tue, 5 Sep 2023 14:26:24 -0600 Subject: [PATCH 29/60] adding qcf_table_file_listed logical to two remaining subrountines; work with log_qcf_info --- .../assimilation/algorithm_info_mod.f90 | 133 ++++++++++++++++-- 1 file changed, 120 insertions(+), 13 deletions(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index e3afcfbda6..c713467522 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -9,7 +9,7 @@ module algorithm_info_mod use obs_def_mod, only : obs_def_type, get_obs_def_type_of_obs, get_obs_def_error_variance use obs_kind_mod, only : get_quantity_for_type_of_obs, get_name_for_quantity -use utilities_mod, only : error_handler, E_ERR +use utilities_mod, only : error_handler, E_ALLMSG, E_ERR, E_MSG, log_it use assim_model_mod, only : get_state_meta_data use location_mod, only : location_type @@ -138,7 +138,11 @@ subroutine init_algorithm_info_mod(qcf_table_filename) allocate(qcf_table_row_headers(numrows)) call read_qcf_table(qcf_table_filename) +!call verify_qcf_table_data(qcf_table_filename, nlines) call write_qcf_table() +call log_qcf_table_data() + +!stop end subroutine init_algorithm_info_mod @@ -164,6 +168,8 @@ subroutine read_qcf_table(qcf_table_filename) ! skip the headers, make sure user is using the correct table version read(fileid, *) header1 read(fileid, *) header2 +write(*,*) 'header1: ', header1 +write(*,*) 'header2: ', header2 ! read in table values directly to qcf_table_data type do row = 1, size(qcf_table_data) @@ -220,9 +226,6 @@ subroutine obs_error_info(obs_def, error_variance, & ! Get the default error variance error_variance = get_obs_def_error_variance(obs_def) -!get actual name of QTY from integer index -kind_name = get_name_for_quantity(obs_kind) - !use default values if qcf_table_filename is not in namelist if (.not. qcf_table_listed) then bounded_below = .false.; bounded_above = .false. @@ -230,11 +233,14 @@ subroutine obs_error_info(obs_def, error_variance, & return endif +!get actual name of QTY from integer index +kind_name = get_name_for_quantity(obs_kind) + !find location of QTY in qcf_table_data structure QTY_loc = findloc(qcf_table_row_headers, kind_name) if (QTY_loc(1) == 0) then - !use default values + !use default values if QTY is not in table bounded_below = .false.; bounded_above = .false. lower_bound = missing_r8; upper_bound = missing_r8 @@ -286,6 +292,14 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & ! In the long run, may not have to have separate controls for each of the input possibilities ! However, for now these are things that need to be explored for science understanding +!use default values if qcf_table_filename is not in namelist +if (.not. qcf_table_listed) then + dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION + bounded_below = .false.; bounded_above = .false. + lower_bound = missing_r8; upper_bound = missing_r8 + return +endif + !get actual name of QTY from integer index kind_name = get_name_for_quantity(kind) @@ -293,9 +307,9 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & QTY_loc = findloc(qcf_table_row_headers, kind_name) if (QTY_loc(1) == 0) then - write(*,*) 'QTY not in table, using default values' + write(*,*) 'QTY not in table, using default values' !remove these writes on PR - !use default values + !use default values if QTY is not in table dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION bounded_below = .false.; bounded_above = .false. lower_bound = missing_r8; upper_bound = missing_r8 @@ -355,16 +369,27 @@ subroutine obs_inc_info(obs_kind, filter_kind, rectangular_quadrature, gaussian_ ! Temporary approach for setting the details of how to assimilate this observation ! This example is designed to reproduce the squared forward operator results from paper +!use default values if qcf_table_filename is not in namelist +if (.not. qcf_table_listed) then + filter_kind = BOUNDED_NORMAL_RHF + bounded_below = .false.; bounded_above = .false. + lower_bound = missing_r8; upper_bound = missing_r8 + sort_obs_inc = .false.; spread_restoration = .false. + return +endif + !get actual name of QTY from integer index kind_name = get_name_for_quantity(obs_kind) +write(*,*) 'kind_name: ', kind_name !find location of QTY in qcf_table_data structure QTY_loc = findloc(qcf_table_row_headers, kind_name) +write(*,*) 'findloc of kind: ', QTY_loc(1) if (QTY_loc(1) == 0) then write(*,*) 'QTY not in table, using default values' - !use default values + !use default values if QTY is not in table filter_kind = BOUNDED_NORMAL_RHF bounded_below = .false.; bounded_above = .false. lower_bound = missing_r8; upper_bound = missing_r8 @@ -382,6 +407,8 @@ subroutine obs_inc_info(obs_kind, filter_kind, rectangular_quadrature, gaussian_ endif +write(*,*) 'obs_inc_info: ', filter_kind, sort_obs_inc, spread_restoration, bounded_below, bounded_above, lower_bound, upper_bound + ! Only need to set these two for options the original RHF implementation !!!rectangular_quadrature = .true. !!!gaussian_likelihood_tails = .false. @@ -405,12 +432,10 @@ subroutine write_qcf_table() integer :: row -write(*,*), 'SIZE: ', size(qcf_table_data) - do row = 1, size(qcf_table_data) - write(*, *) "qcf_table_row_headers(", row, "): ", qcf_table_row_headers(row) - write(*, *) "qcf_table_data(", row, "): " - write(*, *) qcf_table_data(row)%obs_error_info%bounded_below, qcf_table_data(row)%obs_error_info%bounded_above, & + write(*,*) "qcf_table_row_headers(", row, "): ", qcf_table_row_headers(row) + write(*,*) "qcf_table_data(", row, "): " + write(*,*) qcf_table_data(row)%obs_error_info%bounded_below, qcf_table_data(row)%obs_error_info%bounded_above, & qcf_table_data(row)%obs_error_info%lower_bound, qcf_table_data(row)%obs_error_info%upper_bound, qcf_table_data(row)%probit_inflation%dist_type, & qcf_table_data(row)%probit_inflation%bounded_below, qcf_table_data(row)%probit_inflation%bounded_above, & qcf_table_data(row)%probit_inflation%lower_bound, qcf_table_data(row)%probit_inflation%upper_bound, qcf_table_data(row)%probit_state%dist_type, & @@ -441,6 +466,10 @@ subroutine assert_qcf_table_version(header) character(len=129), dimension(4), intent(in) :: header +!if (.not. module_initialized) call init_algorithm_info_mod(qcf_table_filename) + +if (.not. qcf_table_listed) return + write(*,*) 'version: ', header(4) if (header(4) /= '1:') then @@ -453,6 +482,84 @@ end subroutine assert_qcf_table_version !------------------------------------------------------------------------ +subroutine verify_qcf_table_data(qcf_table_filename, nlines) + +!subroutine to ensure that the data in the QCF table is valid and in +!the correct formatthe right format and is correct size + +character(len=129), intent(in) :: qcf_table_filename +integer, intent(in) :: nlines + +character(len=500) :: table_rows(nlines) +integer, parameter :: fileid = 10 !file identifier +integer :: row + +if (.not. module_initialized) call init_algorithm_info_mod(qcf_table_filename) + +if (.not. qcf_table_listed) return + +open(unit=fileid, file=qcf_table_filename) + +do row = 1, nlines + read(fileid, '(A)') table_rows(row) + print *, 'full line:' + print *, table_rows(row) + print *, 'trimmed line:' + print *, trim(table_rows(row)) + print *, 'length', len_trim(table_rows(row)) +end do + +close(fileid) + +!if (size(qcf_table_row_headers) /= 2) then !NO, this needs to be table headers, not row +! write(errstring,*) 'Incorrect number of headers in the QCF table; ' , & +! 'ensure that the latest version of this table is ', & +! 'being used and is in the same format as the example' +! call error_handler(E_ERR, 'assert_qcf_table_version', errstring, source) +!endif + +end subroutine verify_qcf_table_data + +!------------------------------------------------------------------------ + + +subroutine log_qcf_table_data() + +!subroutine to write the data in QCF table to dart_log + +character(len=500) :: log_msg +integer :: row + +if (.not. qcf_table_listed) return + +do row = 1, size(qcf_table_data) + write(log_msg, *) "qcf_table_row_headers(", row, "): ", qcf_table_row_headers(row) + write(*,*) 'log_msg: ', log_msg + call log_it(log_msg) + write(log_msg, *) "qcf_table_data(", row, "): " + write(*,*) 'log_msg: ', log_msg + call log_it(log_msg) + write(log_msg, *) qcf_table_data(row)%obs_error_info%bounded_below, qcf_table_data(row)%obs_error_info%bounded_above, & + qcf_table_data(row)%obs_error_info%lower_bound, qcf_table_data(row)%obs_error_info%upper_bound, qcf_table_data(row)%probit_inflation%dist_type, & + qcf_table_data(row)%probit_inflation%bounded_below, qcf_table_data(row)%probit_inflation%bounded_above, & + qcf_table_data(row)%probit_inflation%lower_bound, qcf_table_data(row)%probit_inflation%upper_bound, qcf_table_data(row)%probit_state%dist_type, & + qcf_table_data(row)%probit_state%bounded_below, qcf_table_data(row)%probit_state%bounded_above, & + qcf_table_data(row)%probit_state%lower_bound, qcf_table_data(row)%probit_state%upper_bound, qcf_table_data(row)%probit_extended_state%dist_type, & + qcf_table_data(row)%probit_extended_state%bounded_below, qcf_table_data(row)%probit_extended_state%bounded_above, & + qcf_table_data(row)%probit_extended_state%lower_bound, qcf_table_data(row)%probit_extended_state%upper_bound, & + qcf_table_data(row)%obs_inc_info%filter_kind, qcf_table_data(row)%obs_inc_info%rectangular_quadrature, & + qcf_table_data(row)%obs_inc_info%gaussian_likelihood_tails, qcf_table_data(row)%obs_inc_info%sort_obs_inc, & + qcf_table_data(row)%obs_inc_info%spread_restoration, qcf_table_data(row)%obs_inc_info%bounded_below, qcf_table_data(row)%obs_inc_info%bounded_above, & + qcf_table_data(row)%obs_inc_info%lower_bound, qcf_table_data(row)%obs_inc_info%upper_bound + write(*,*) 'e_allmsg: ' + call error_handler(E_ALLMSG, 'write_qcf_table', log_msg, source) +end do + +end subroutine log_qcf_table_data + +!------------------------------------------------------------------------ + + subroutine end_algorithm_info_mod() if (.not. module_initialized) return From e3f3131e885684246f08cbeaa23187586e06623d Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Tue, 5 Sep 2023 15:16:37 -0600 Subject: [PATCH 30/60] commenting out prints; remove these lines later --- .../assimilation/algorithm_info_mod.f90 | 28 +++++++++---------- 1 file changed, 14 insertions(+), 14 deletions(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index c713467522..f09d3e1656 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -113,7 +113,7 @@ subroutine init_algorithm_info_mod(qcf_table_filename) if (module_initialized) return module_initialized = .true. -write(*,*) 'filename: ', qcf_table_filename +!write(*,*) 'filename: ', qcf_table_filename if (qcf_table_filename == '') then write(*,*) 'No QCF table file listed in namelist, using default values for all QTYs' @@ -132,14 +132,14 @@ subroutine init_algorithm_info_mod(qcf_table_filename) close(fileid) numrows = nlines - 2 -print *, 'numrows: ', numrows +!print *, 'numrows: ', numrows allocate(qcf_table_data(numrows)) allocate(qcf_table_row_headers(numrows)) call read_qcf_table(qcf_table_filename) !call verify_qcf_table_data(qcf_table_filename, nlines) -call write_qcf_table() +!call write_qcf_table() call log_qcf_table_data() !stop @@ -168,8 +168,8 @@ subroutine read_qcf_table(qcf_table_filename) ! skip the headers, make sure user is using the correct table version read(fileid, *) header1 read(fileid, *) header2 -write(*,*) 'header1: ', header1 -write(*,*) 'header2: ', header2 +!write(*,*) 'header1: ', header1 +!write(*,*) 'header2: ', header2 ! read in table values directly to qcf_table_data type do row = 1, size(qcf_table_data) @@ -307,7 +307,7 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & QTY_loc = findloc(qcf_table_row_headers, kind_name) if (QTY_loc(1) == 0) then - write(*,*) 'QTY not in table, using default values' !remove these writes on PR + ! write(*,*) 'QTY not in table, using default values' !remove these writes on PR !use default values if QTY is not in table dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION @@ -380,14 +380,14 @@ subroutine obs_inc_info(obs_kind, filter_kind, rectangular_quadrature, gaussian_ !get actual name of QTY from integer index kind_name = get_name_for_quantity(obs_kind) -write(*,*) 'kind_name: ', kind_name +!write(*,*) 'kind_name: ', kind_name !find location of QTY in qcf_table_data structure QTY_loc = findloc(qcf_table_row_headers, kind_name) -write(*,*) 'findloc of kind: ', QTY_loc(1) +!write(*,*) 'findloc of kind: ', QTY_loc(1) if (QTY_loc(1) == 0) then - write(*,*) 'QTY not in table, using default values' + !write(*,*) 'QTY not in table, using default values' !use default values if QTY is not in table filter_kind = BOUNDED_NORMAL_RHF @@ -407,7 +407,7 @@ subroutine obs_inc_info(obs_kind, filter_kind, rectangular_quadrature, gaussian_ endif -write(*,*) 'obs_inc_info: ', filter_kind, sort_obs_inc, spread_restoration, bounded_below, bounded_above, lower_bound, upper_bound +!write(*,*) 'obs_inc_info: ', filter_kind, sort_obs_inc, spread_restoration, bounded_below, bounded_above, lower_bound, upper_bound ! Only need to set these two for options the original RHF implementation !!!rectangular_quadrature = .true. @@ -470,7 +470,7 @@ subroutine assert_qcf_table_version(header) if (.not. qcf_table_listed) return -write(*,*) 'version: ', header(4) +!write(*,*) 'version: ', header(4) if (header(4) /= '1:') then write(errstring,*) "Using outdated/incorrect version of the QCF table" @@ -534,10 +534,10 @@ subroutine log_qcf_table_data() do row = 1, size(qcf_table_data) write(log_msg, *) "qcf_table_row_headers(", row, "): ", qcf_table_row_headers(row) - write(*,*) 'log_msg: ', log_msg + ! print *, 'log_msg: ', log_msg call log_it(log_msg) write(log_msg, *) "qcf_table_data(", row, "): " - write(*,*) 'log_msg: ', log_msg + ! print *, 'log_msg: ', log_msg call log_it(log_msg) write(log_msg, *) qcf_table_data(row)%obs_error_info%bounded_below, qcf_table_data(row)%obs_error_info%bounded_above, & qcf_table_data(row)%obs_error_info%lower_bound, qcf_table_data(row)%obs_error_info%upper_bound, qcf_table_data(row)%probit_inflation%dist_type, & @@ -551,7 +551,7 @@ subroutine log_qcf_table_data() qcf_table_data(row)%obs_inc_info%gaussian_likelihood_tails, qcf_table_data(row)%obs_inc_info%sort_obs_inc, & qcf_table_data(row)%obs_inc_info%spread_restoration, qcf_table_data(row)%obs_inc_info%bounded_below, qcf_table_data(row)%obs_inc_info%bounded_above, & qcf_table_data(row)%obs_inc_info%lower_bound, qcf_table_data(row)%obs_inc_info%upper_bound - write(*,*) 'e_allmsg: ' + ! print *, 'e_allmsg: ' call error_handler(E_ALLMSG, 'write_qcf_table', log_msg, source) end do From 05c5fffd97e76293271ba3c106e8bf7ae3f0efdf Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Mon, 11 Sep 2023 12:43:37 -0600 Subject: [PATCH 31/60] Moving initialize_modules call to be after read of filter_nml --- assimilation_code/modules/assimilation/filter_mod.f90 | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/assimilation_code/modules/assimilation/filter_mod.f90 b/assimilation_code/modules/assimilation/filter_mod.f90 index 874a4cb1d5..c83a433a8a 100644 --- a/assimilation_code/modules/assimilation/filter_mod.f90 +++ b/assimilation_code/modules/assimilation/filter_mod.f90 @@ -166,7 +166,7 @@ module filter_mod !---------------------------------------------------------------- ! Namelist input with default values ! -character(len = 129) :: qcf_table_filename = 'real_qcf_table.txt' !NEED TO REMOVE THIS LATER +character(len = 129) :: qcf_table_filename = '' logical :: use_algorithm_info_mod = .true. integer :: async = 0, ens_size = 20 integer :: tasks_per_model_advance = 1 @@ -363,13 +363,13 @@ subroutine filter_main() real(r8), allocatable :: prior_qc_copy(:) -call filter_initialize_modules_used() ! static_init_model called in here - ! Read the namelist entry call find_namelist_in_file("input.nml", "filter_nml", iunit) read(iunit, nml = filter_nml, iostat = io) call check_namelist_read(iunit, io, "filter_nml") +call filter_initialize_modules_used() ! static_init_model called in here + ! Record the namelist values used for the run ... if (do_nml_file()) write(nmlfileunit, nml=filter_nml) if (do_nml_term()) write( * , nml=filter_nml) From 27192266cad79b50cae6a277614bd12626af67a0 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Mon, 11 Sep 2023 16:09:51 -0600 Subject: [PATCH 32/60] Adding check that all QTYs in the table exist in DART using get_index_for_quantity; adding checks that upper bounds are greater than lower bounds --- .../assimilation/algorithm_info_mod.f90 | 145 ++++++++---------- 1 file changed, 66 insertions(+), 79 deletions(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index f09d3e1656..1555ec5fdb 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -7,9 +7,9 @@ module algorithm_info_mod use types_mod, only : r8, i8, missing_r8 use obs_def_mod, only : obs_def_type, get_obs_def_type_of_obs, get_obs_def_error_variance -use obs_kind_mod, only : get_quantity_for_type_of_obs, get_name_for_quantity +use obs_kind_mod, only : get_quantity_for_type_of_obs, get_name_for_quantity, get_index_for_quantity -use utilities_mod, only : error_handler, E_ALLMSG, E_ERR, E_MSG, log_it +use utilities_mod, only : error_handler, E_ALLMSG, E_ERR, E_MSG, log_it, logfileunit use assim_model_mod, only : get_state_meta_data use location_mod, only : location_type @@ -21,6 +21,7 @@ module algorithm_info_mod implicit none private +character(len=2000) :: log_msg character(len=512) :: errstring character(len=*), parameter :: source = 'algorithm_info_mod.f90' @@ -39,8 +40,6 @@ module algorithm_info_mod public :: obs_error_info, probit_dist_info, obs_inc_info, & init_algorithm_info_mod, end_algorithm_info_mod, & - obs_error_info_type, probit_inflation_type, probit_state_type, & - probit_extended_state_type, obs_inc_info_type, qcf_table_data_type, & EAKF, ENKF, BOUNDED_NORMAL_RHF, UNBOUNDED_RHF, GAMMA_FILTER !Creates the type definitions for the QCF table @@ -83,8 +82,11 @@ module algorithm_info_mod type(obs_inc_info_type) :: obs_inc_info end type +character(len=129), dimension(4) :: header1 +character(len=129), dimension(29) :: header2 + +character(len=129), allocatable :: qcf_table_row_headers(:) type(qcf_table_data_type), allocatable :: qcf_table_data(:) -character(len=129), allocatable :: qcf_table_row_headers(:) ! Provides routines that give information about details of algorithms for ! observation error sampling, observation increments, and the transformations @@ -113,8 +115,6 @@ subroutine init_algorithm_info_mod(qcf_table_filename) if (module_initialized) return module_initialized = .true. -!write(*,*) 'filename: ', qcf_table_filename - if (qcf_table_filename == '') then write(*,*) 'No QCF table file listed in namelist, using default values for all QTYs' return @@ -124,7 +124,8 @@ subroutine init_algorithm_info_mod(qcf_table_filename) open(unit=fileid, file=qcf_table_filename) nlines = 0 -do !do loop to get number of rows (or QTY's) in the table +!do loop to get number of rows (or QTY's) in the table +do read(fileid,*,iostat=io) if(io/=0) exit nlines = nlines + 1 @@ -132,17 +133,17 @@ subroutine init_algorithm_info_mod(qcf_table_filename) close(fileid) numrows = nlines - 2 -!print *, 'numrows: ', numrows allocate(qcf_table_data(numrows)) allocate(qcf_table_row_headers(numrows)) call read_qcf_table(qcf_table_filename) -!call verify_qcf_table_data(qcf_table_filename, nlines) !call write_qcf_table() +call assert_qcf_table_version() +call verify_qcf_table_data() call log_qcf_table_data() -!stop +stop !FOR TESTING, REMOVE LATER end subroutine init_algorithm_info_mod @@ -158,9 +159,6 @@ subroutine read_qcf_table(qcf_table_filename) integer, parameter :: fileid = 10 !file identifier integer :: row -character(len=129), dimension(4) :: header1 -character(len=129), dimension(29) :: header2 - if (.not. module_initialized) call init_algorithm_info_mod(qcf_table_filename) open(unit=fileid, file=qcf_table_filename) @@ -168,8 +166,8 @@ subroutine read_qcf_table(qcf_table_filename) ! skip the headers, make sure user is using the correct table version read(fileid, *) header1 read(fileid, *) header2 -!write(*,*) 'header1: ', header1 -!write(*,*) 'header2: ', header2 +write(*,*) 'header1: ', header1 +write(*,*) 'header2: ', header2 ! read in table values directly to qcf_table_data type do row = 1, size(qcf_table_data) @@ -189,8 +187,6 @@ subroutine read_qcf_table(qcf_table_filename) close(fileid) -call assert_qcf_table_version(header1) - end subroutine read_qcf_table !------------------------------------------------------------------------ @@ -380,11 +376,9 @@ subroutine obs_inc_info(obs_kind, filter_kind, rectangular_quadrature, gaussian_ !get actual name of QTY from integer index kind_name = get_name_for_quantity(obs_kind) -!write(*,*) 'kind_name: ', kind_name !find location of QTY in qcf_table_data structure QTY_loc = findloc(qcf_table_row_headers, kind_name) -!write(*,*) 'findloc of kind: ', QTY_loc(1) if (QTY_loc(1) == 0) then !write(*,*) 'QTY not in table, using default values' @@ -424,7 +418,7 @@ subroutine write_qcf_table() ! write to check values were correctly assigned ! testing for findloc -character(len=30), parameter :: tester_QTY = 'QTY_GPSRO' +character(len=30), parameter :: tester_QTY = 'QTY_STATE_VARIABLE' integer :: QTY_loc(1) character(len=30), parameter :: tester_QTY0 = 'QTY_DUMMY' @@ -450,7 +444,7 @@ subroutine write_qcf_table() end do QTY_loc = findloc(qcf_table_row_headers, tester_QTY) -write(*, *) 'findloc of QTY_GPSRO: ', QTY_loc(1) +write(*, *) 'findloc of QTY_STATE_VARIABLE: ', QTY_loc(1) QTY_loc0 = findloc(qcf_table_row_headers, tester_QTY0) write(*, *) 'findloc of invalid QTY (QTY_DUMMY): ', QTY_loc0(1) @@ -460,19 +454,11 @@ end subroutine write_qcf_table !------------------------------------------------------------------------ -subroutine assert_qcf_table_version(header) +subroutine assert_qcf_table_version() !subroutine to ensure the correct version of the QCF table is being used -character(len=129), dimension(4), intent(in) :: header - -!if (.not. module_initialized) call init_algorithm_info_mod(qcf_table_filename) - -if (.not. qcf_table_listed) return - -!write(*,*) 'version: ', header(4) - -if (header(4) /= '1:') then +if (header1(4) /= '1:') then write(errstring,*) "Using outdated/incorrect version of the QCF table" call error_handler(E_ERR, 'assert_qcf_table_version', errstring, source) endif @@ -482,41 +468,50 @@ end subroutine assert_qcf_table_version !------------------------------------------------------------------------ -subroutine verify_qcf_table_data(qcf_table_filename, nlines) - -!subroutine to ensure that the data in the QCF table is valid and in -!the correct formatthe right format and is correct size +subroutine verify_qcf_table_data() -character(len=129), intent(in) :: qcf_table_filename -integer, intent(in) :: nlines +!subroutine to ensure that the data in the QCF table is valid -character(len=500) :: table_rows(nlines) -integer, parameter :: fileid = 10 !file identifier +integer :: varid integer :: row -if (.not. module_initialized) call init_algorithm_info_mod(qcf_table_filename) +!if (.not. module_initialized) call init_algorithm_info_mod(qcf_table_filename) if (.not. qcf_table_listed) return -open(unit=fileid, file=qcf_table_filename) - -do row = 1, nlines - read(fileid, '(A)') table_rows(row) - print *, 'full line:' - print *, table_rows(row) - print *, 'trimmed line:' - print *, trim(table_rows(row)) - print *, 'length', len_trim(table_rows(row)) +!Checks that all bounds are valid; currently checks that the lower bound in less than the upper +!Here we could add more specific checks if we have known limits on the bounds +do row = 1, size(qcf_table_data) + if(qcf_table_data(row)%obs_error_info%lower_bound > qcf_table_data(row)%obs_error_info%upper_bound) then + write(errstring,*) "Invalid bounds in obs_error_info" + call error_handler(E_ERR, 'verify_qcf_table_data', errstring, source) + endif + if(qcf_table_data(row)%probit_inflation%lower_bound > qcf_table_data(row)%probit_inflation%upper_bound) then + write(errstring,*) "Invalid bounds in probit_inflation" + call error_handler(E_ERR, 'verify_qcf_table_data', errstring, source) + endif + if(qcf_table_data(row)%probit_state%lower_bound > qcf_table_data(row)%probit_state%upper_bound) then + write(errstring,*) "Invalid bounds in probit_state" + call error_handler(E_ERR, 'verify_qcf_table_data', errstring, source) + endif + if(qcf_table_data(row)%probit_extended_state%lower_bound > qcf_table_data(row)%probit_extended_state%upper_bound) then + write(errstring,*) "Invalid bounds in probit_extended_state" + call error_handler(E_ERR, 'verify_qcf_table_data', errstring, source) + endif + if(qcf_table_data(row)%obs_inc_info%lower_bound > qcf_table_data(row)%obs_inc_info%upper_bound) then + write(errstring,*) "Invalid bounds in obs_inc_info" + call error_handler(E_ERR, 'verify_qcf_table_data', errstring, source) + endif end do -close(fileid) - -!if (size(qcf_table_row_headers) /= 2) then !NO, this needs to be table headers, not row -! write(errstring,*) 'Incorrect number of headers in the QCF table; ' , & -! 'ensure that the latest version of this table is ', & -! 'being used and is in the same format as the example' -! call error_handler(E_ERR, 'assert_qcf_table_version', errstring, source) -!endif +!Ensures that all QTYs listed in the table exist in DART +do row = 1, size(qcf_table_data) + varid = get_index_for_quantity(qcf_table_row_headers(row)) + if(varid == -1) then + write(errstring,*) trim(qcf_table_row_headers(row)), " is not a valid DART QTY" + call error_handler(E_ERR, 'verify_qcf_table_data', errstring, source) + endif +end do end subroutine verify_qcf_table_data @@ -527,34 +522,26 @@ subroutine log_qcf_table_data() !subroutine to write the data in QCF table to dart_log -character(len=500) :: log_msg integer :: row if (.not. qcf_table_listed) return +!call error_handler(E_ALLMSG, 'log_qcf_table_data', log_msg, source) +!call log_it(log_msg) + +write(logfileunit, *) header1 +write(logfileunit, *) header2 + do row = 1, size(qcf_table_data) - write(log_msg, *) "qcf_table_row_headers(", row, "): ", qcf_table_row_headers(row) - ! print *, 'log_msg: ', log_msg - call log_it(log_msg) - write(log_msg, *) "qcf_table_data(", row, "): " - ! print *, 'log_msg: ', log_msg - call log_it(log_msg) - write(log_msg, *) qcf_table_data(row)%obs_error_info%bounded_below, qcf_table_data(row)%obs_error_info%bounded_above, & - qcf_table_data(row)%obs_error_info%lower_bound, qcf_table_data(row)%obs_error_info%upper_bound, qcf_table_data(row)%probit_inflation%dist_type, & - qcf_table_data(row)%probit_inflation%bounded_below, qcf_table_data(row)%probit_inflation%bounded_above, & - qcf_table_data(row)%probit_inflation%lower_bound, qcf_table_data(row)%probit_inflation%upper_bound, qcf_table_data(row)%probit_state%dist_type, & - qcf_table_data(row)%probit_state%bounded_below, qcf_table_data(row)%probit_state%bounded_above, & - qcf_table_data(row)%probit_state%lower_bound, qcf_table_data(row)%probit_state%upper_bound, qcf_table_data(row)%probit_extended_state%dist_type, & - qcf_table_data(row)%probit_extended_state%bounded_below, qcf_table_data(row)%probit_extended_state%bounded_above, & - qcf_table_data(row)%probit_extended_state%lower_bound, qcf_table_data(row)%probit_extended_state%upper_bound, & - qcf_table_data(row)%obs_inc_info%filter_kind, qcf_table_data(row)%obs_inc_info%rectangular_quadrature, & - qcf_table_data(row)%obs_inc_info%gaussian_likelihood_tails, qcf_table_data(row)%obs_inc_info%sort_obs_inc, & - qcf_table_data(row)%obs_inc_info%spread_restoration, qcf_table_data(row)%obs_inc_info%bounded_below, qcf_table_data(row)%obs_inc_info%bounded_above, & - qcf_table_data(row)%obs_inc_info%lower_bound, qcf_table_data(row)%obs_inc_info%upper_bound - ! print *, 'e_allmsg: ' - call error_handler(E_ALLMSG, 'write_qcf_table', log_msg, source) + write(*,*) qcf_table_row_headers(row), qcf_table_data(row) end do +!write(log_msg,*) qcf_table_data +!write(*, *) trim(log_msg) +!write(logfileunit, *) trim(log_msg) +!call log_it(trim(log_msg)) +!call error_handler(E_MSG, 'log_qcf_table_data', trim(log_msg), source) + end subroutine log_qcf_table_data !------------------------------------------------------------------------ From 8f79b3c4027344c687e6b7002d9e67075b693d0c Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Mon, 11 Sep 2023 16:56:50 -0600 Subject: [PATCH 33/60] Replacing open and close of qcf_table file with open_file and close_file of utilities_mod - this fixes the issue with writing to the dart_log --- .../assimilation/algorithm_info_mod.f90 | 19 +++++++++++-------- 1 file changed, 11 insertions(+), 8 deletions(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index 1555ec5fdb..fbebef9cfb 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -9,7 +9,7 @@ module algorithm_info_mod use obs_def_mod, only : obs_def_type, get_obs_def_type_of_obs, get_obs_def_error_variance use obs_kind_mod, only : get_quantity_for_type_of_obs, get_name_for_quantity, get_index_for_quantity -use utilities_mod, only : error_handler, E_ALLMSG, E_ERR, E_MSG, log_it, logfileunit +use utilities_mod, only : error_handler, E_ALLMSG, E_ERR, E_MSG, log_it, logfileunit, open_file, close_file use assim_model_mod, only : get_state_meta_data use location_mod, only : location_type @@ -107,10 +107,11 @@ subroutine init_algorithm_info_mod(qcf_table_filename) character(len=129), intent(in) :: qcf_table_filename +integer :: fileid +integer :: io + integer :: numrows integer :: nlines -integer :: io -integer, parameter :: fileid = 10 !file identifier if (module_initialized) return module_initialized = .true. @@ -121,16 +122,18 @@ subroutine init_algorithm_info_mod(qcf_table_filename) endif qcf_table_listed = .true. -open(unit=fileid, file=qcf_table_filename) nlines = 0 +fileid = open_file(qcf_table_filename, 'formatted', 'read') + !do loop to get number of rows (or QTY's) in the table do read(fileid,*,iostat=io) if(io/=0) exit nlines = nlines + 1 end do -close(fileid) + +call close_file(fileid) numrows = nlines - 2 @@ -156,12 +159,12 @@ subroutine read_qcf_table(qcf_table_filename) character(len=129), intent(in) :: qcf_table_filename -integer, parameter :: fileid = 10 !file identifier +integer :: fileid integer :: row if (.not. module_initialized) call init_algorithm_info_mod(qcf_table_filename) -open(unit=fileid, file=qcf_table_filename) +fileid = open_file(qcf_table_filename, 'formatted', 'read') ! skip the headers, make sure user is using the correct table version read(fileid, *) header1 @@ -185,7 +188,7 @@ subroutine read_qcf_table(qcf_table_filename) qcf_table_data(row)%obs_inc_info%lower_bound, qcf_table_data(row)%obs_inc_info%upper_bound end do -close(fileid) +call close_file(fileid) end subroutine read_qcf_table From c0b00e06cd4f105ad3910fc3124a1978806d8e0b Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Mon, 11 Sep 2023 17:03:43 -0600 Subject: [PATCH 34/60] Making quotations consistent - using single quote only --- .../assimilation/algorithm_info_mod.f90 | 18 +++++++++--------- 1 file changed, 9 insertions(+), 9 deletions(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index fbebef9cfb..650c79f7a1 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -430,8 +430,8 @@ subroutine write_qcf_table() integer :: row do row = 1, size(qcf_table_data) - write(*,*) "qcf_table_row_headers(", row, "): ", qcf_table_row_headers(row) - write(*,*) "qcf_table_data(", row, "): " + write(*,*) 'qcf_table_row_headers(', row, '): ', qcf_table_row_headers(row) + write(*,*) 'qcf_table_data(', row, '): ' write(*,*) qcf_table_data(row)%obs_error_info%bounded_below, qcf_table_data(row)%obs_error_info%bounded_above, & qcf_table_data(row)%obs_error_info%lower_bound, qcf_table_data(row)%obs_error_info%upper_bound, qcf_table_data(row)%probit_inflation%dist_type, & qcf_table_data(row)%probit_inflation%bounded_below, qcf_table_data(row)%probit_inflation%bounded_above, & @@ -462,7 +462,7 @@ subroutine assert_qcf_table_version() !subroutine to ensure the correct version of the QCF table is being used if (header1(4) /= '1:') then - write(errstring,*) "Using outdated/incorrect version of the QCF table" + write(errstring,*) 'Using outdated/incorrect version of the QCF table' call error_handler(E_ERR, 'assert_qcf_table_version', errstring, source) endif @@ -486,23 +486,23 @@ subroutine verify_qcf_table_data() !Here we could add more specific checks if we have known limits on the bounds do row = 1, size(qcf_table_data) if(qcf_table_data(row)%obs_error_info%lower_bound > qcf_table_data(row)%obs_error_info%upper_bound) then - write(errstring,*) "Invalid bounds in obs_error_info" + write(errstring,*) 'Invalid bounds in obs_error_info' call error_handler(E_ERR, 'verify_qcf_table_data', errstring, source) endif if(qcf_table_data(row)%probit_inflation%lower_bound > qcf_table_data(row)%probit_inflation%upper_bound) then - write(errstring,*) "Invalid bounds in probit_inflation" + write(errstring,*) 'Invalid bounds in probit_inflation' call error_handler(E_ERR, 'verify_qcf_table_data', errstring, source) endif if(qcf_table_data(row)%probit_state%lower_bound > qcf_table_data(row)%probit_state%upper_bound) then - write(errstring,*) "Invalid bounds in probit_state" + write(errstring,*) 'Invalid bounds in probit_state' call error_handler(E_ERR, 'verify_qcf_table_data', errstring, source) endif if(qcf_table_data(row)%probit_extended_state%lower_bound > qcf_table_data(row)%probit_extended_state%upper_bound) then - write(errstring,*) "Invalid bounds in probit_extended_state" + write(errstring,*) 'Invalid bounds in probit_extended_state' call error_handler(E_ERR, 'verify_qcf_table_data', errstring, source) endif if(qcf_table_data(row)%obs_inc_info%lower_bound > qcf_table_data(row)%obs_inc_info%upper_bound) then - write(errstring,*) "Invalid bounds in obs_inc_info" + write(errstring,*) 'Invalid bounds in obs_inc_info' call error_handler(E_ERR, 'verify_qcf_table_data', errstring, source) endif end do @@ -511,7 +511,7 @@ subroutine verify_qcf_table_data() do row = 1, size(qcf_table_data) varid = get_index_for_quantity(qcf_table_row_headers(row)) if(varid == -1) then - write(errstring,*) trim(qcf_table_row_headers(row)), " is not a valid DART QTY" + write(errstring,*) trim(qcf_table_row_headers(row)), ' is not a valid DART QTY' call error_handler(E_ERR, 'verify_qcf_table_data', errstring, source) endif end do From 0ace5e96b19313e3c1729a869b2be74539bdecca Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Tue, 12 Sep 2023 15:07:22 -0600 Subject: [PATCH 35/60] Adding checks for duplicate QTYs in the table --- .../modules/assimilation/algorithm_info_mod.f90 | 10 +++++++++- 1 file changed, 9 insertions(+), 1 deletion(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index 650c79f7a1..00b3ff651a 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -461,7 +461,7 @@ subroutine assert_qcf_table_version() !subroutine to ensure the correct version of the QCF table is being used -if (header1(4) /= '1:') then +if (header1(4) /= '1') then write(errstring,*) 'Using outdated/incorrect version of the QCF table' call error_handler(E_ERR, 'assert_qcf_table_version', errstring, source) endif @@ -516,6 +516,14 @@ subroutine verify_qcf_table_data() endif end do +!Ensures that there are no duplicate QTYs in the table +do row = 1, size(qcf_table_data) + if(count(qcf_table_row_headers==qcf_table_row_headers(row)) > 1) then + write(errstring,*) trim(qcf_table_row_headers(row)), ' has multiple entries in the table' + call error_handler(E_ERR, 'verify_qcf_table_data', errstring, source) + endif +end do + end subroutine verify_qcf_table_data !------------------------------------------------------------------------ From b32ac9307ccab9f41f5a7e7e6b82822cd33233e0 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Wed, 13 Sep 2023 12:51:03 -0600 Subject: [PATCH 36/60] Switching to real(r8); experimenting in the log_qcf_table subroutine --- .../assimilation/algorithm_info_mod.f90 | 65 ++++++++++++------- 1 file changed, 42 insertions(+), 23 deletions(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index 00b3ff651a..5d215f6e2d 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -21,7 +21,6 @@ module algorithm_info_mod implicit none private -character(len=2000) :: log_msg character(len=512) :: errstring character(len=*), parameter :: source = 'algorithm_info_mod.f90' @@ -45,41 +44,41 @@ module algorithm_info_mod !Creates the type definitions for the QCF table type obs_error_info_type logical :: bounded_below, bounded_above - real :: lower_bound, upper_bound + real(r8) :: lower_bound, upper_bound end type type probit_inflation_type integer :: dist_type logical :: bounded_below, bounded_above - real :: lower_bound, upper_bound + real(r8) :: lower_bound, upper_bound end type type probit_state_type integer :: dist_type logical :: bounded_below, bounded_above - real :: lower_bound, upper_bound + real(r8) :: lower_bound, upper_bound end type type probit_extended_state_type integer :: dist_type logical :: bounded_below, bounded_above - real :: lower_bound, upper_bound + real(r8) :: lower_bound, upper_bound end type type obs_inc_info_type - integer :: filter_kind - logical :: rectangular_quadrature, gaussian_likelihood_tails - logical :: sort_obs_inc, spread_restoration - logical :: bounded_below, bounded_above - real :: lower_bound, upper_bound + integer :: filter_kind + logical :: rectangular_quadrature, gaussian_likelihood_tails + logical :: sort_obs_inc, spread_restoration + logical :: bounded_below, bounded_above + real(r8) :: lower_bound, upper_bound end type type qcf_table_data_type - type(obs_error_info_type) :: obs_error_info - type(probit_inflation_type) :: probit_inflation - type(probit_state_type) :: probit_state + type(obs_error_info_type) :: obs_error_info + type(probit_inflation_type) :: probit_inflation + type(probit_state_type) :: probit_state type(probit_extended_state_type) :: probit_extended_state - type(obs_inc_info_type) :: obs_inc_info + type(obs_inc_info_type) :: obs_inc_info end type character(len=129), dimension(4) :: header1 @@ -117,7 +116,8 @@ subroutine init_algorithm_info_mod(qcf_table_filename) module_initialized = .true. if (qcf_table_filename == '') then - write(*,*) 'No QCF table file listed in namelist, using default values for all QTYs' + write(errstring,*) 'No QCF table file listed in namelist, using default values for all QTYs' + call error_handler(E_MSG, 'init_algorithm_info_mod', errstring, source) return endif @@ -144,7 +144,7 @@ subroutine init_algorithm_info_mod(qcf_table_filename) !call write_qcf_table() call assert_qcf_table_version() call verify_qcf_table_data() -call log_qcf_table_data() +!call log_qcf_table_data() stop !FOR TESTING, REMOVE LATER @@ -169,8 +169,8 @@ subroutine read_qcf_table(qcf_table_filename) ! skip the headers, make sure user is using the correct table version read(fileid, *) header1 read(fileid, *) header2 -write(*,*) 'header1: ', header1 -write(*,*) 'header2: ', header2 +!write(*,*) 'header1: ', header1 +!write(*,*) 'header2: ', header2 ! read in table values directly to qcf_table_data type do row = 1, size(qcf_table_data) @@ -533,21 +533,40 @@ subroutine log_qcf_table_data() !subroutine to write the data in QCF table to dart_log +character(len=512) :: log_msg integer :: row +integer :: i if (.not. qcf_table_listed) return !call error_handler(E_ALLMSG, 'log_qcf_table_data', log_msg, source) !call log_it(log_msg) -write(logfileunit, *) header1 -write(logfileunit, *) header2 +!Write the headers to the dart_log and terminal -do row = 1, size(qcf_table_data) - write(*,*) qcf_table_row_headers(row), qcf_table_data(row) +!write(log_msg, *) trim(header1(1)) +do i = 1, size(header1) + write(*,*) trim(header1(i)) + write(log_msg,*) trim(header1(i)) + !log_msg = log_msg//trim(header1(i)) + ! call error_handler(E_MSG, 'log_qcf_table_data', trim(log_msg), source) + write(*, *) trim(log_msg) end do +write(*,*) 'log_msg: ', trim(log_msg) +!write(*,*) trim(log_msg) +!write(log_msg,*) header1 +!call error_handler(E_MSG, 'log_qcf_table_data', trim(log_msg), source) +!write(log_msg, *) header2 +!call error_handler(E_MSG, 'log_qcf_table_data', trim(log_msg), source) + +!Write the data to the dart_log and terminal +!do row = 1, size(qcf_table_data) +! write(log_msg,*) qcf_table_row_headers(row), qcf_table_data(row) +! call error_handler(E_MSG, 'log_qcf_table_data', trim(log_msg), source) +!end do -!write(log_msg,*) qcf_table_data +write(log_msg,*) qcf_table_data +write(*,*) 'log_msg: ', trim(log_msg) !write(*, *) trim(log_msg) !write(logfileunit, *) trim(log_msg) !call log_it(trim(log_msg)) From 1c5cfa4e09843ae7fcec6bb305ab4673e83a2a81 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Wed, 13 Sep 2023 17:49:45 -0600 Subject: [PATCH 37/60] reading in character strings from QCF table for filter_kind instead of ints --- .../assimilation/algorithm_info_mod.f90 | 55 +++++++++++++++++-- 1 file changed, 49 insertions(+), 6 deletions(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index 5d215f6e2d..42bdca1acf 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -9,7 +9,7 @@ module algorithm_info_mod use obs_def_mod, only : obs_def_type, get_obs_def_type_of_obs, get_obs_def_error_variance use obs_kind_mod, only : get_quantity_for_type_of_obs, get_name_for_quantity, get_index_for_quantity -use utilities_mod, only : error_handler, E_ALLMSG, E_ERR, E_MSG, log_it, logfileunit, open_file, close_file +use utilities_mod, only : error_handler, E_ALLMSG, E_ERR, E_MSG, log_it, logfileunit, open_file, close_file, to_upper use assim_model_mod, only : get_state_meta_data use location_mod, only : location_type @@ -66,7 +66,7 @@ module algorithm_info_mod end type type obs_inc_info_type - integer :: filter_kind + character(len=129) :: filter_kind logical :: rectangular_quadrature, gaussian_likelihood_tails logical :: sort_obs_inc, spread_restoration logical :: bounded_below, bounded_above @@ -141,12 +141,12 @@ subroutine init_algorithm_info_mod(qcf_table_filename) allocate(qcf_table_row_headers(numrows)) call read_qcf_table(qcf_table_filename) -!call write_qcf_table() +call write_qcf_table() call assert_qcf_table_version() call verify_qcf_table_data() !call log_qcf_table_data() -stop !FOR TESTING, REMOVE LATER +!stop !FOR TESTING, REMOVE LATER end subroutine init_algorithm_info_mod @@ -270,6 +270,7 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & real(r8), intent(out) :: lower_bound, upper_bound integer :: QTY_loc(1) +character(len=129) :: dist_type_string character(len=129) :: kind_name ! Have input information about the kind of the state or observation being transformed @@ -316,6 +317,8 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & elseif(is_inflation) then ! Case for inflation transformation +! dist_type_string = to_upper(qcf_table_data(QTY_loc(1))%probit_inflation%dist_type) + dist_type = qcf_table_data(QTY_loc(1))%probit_inflation%dist_type !dist_type has checks in transform_to_probit, transform_from_probit bounded_below = qcf_table_data(QTY_loc(1))%probit_inflation%bounded_below bounded_above = qcf_table_data(QTY_loc(1))%probit_inflation%bounded_above @@ -361,6 +364,11 @@ subroutine obs_inc_info(obs_kind, filter_kind, rectangular_quadrature, gaussian_ integer :: QTY_loc(1) character(len=129) :: kind_name +integer :: filter_kind_loc(1) +character(len=129), dimension(5) :: possible_filter_kinds +integer, dimension(5) :: possible_filter_kind_ints +character(len=129) :: filter_kind_string + ! The information arguments are all intent (inout). This means that if they are not set ! here, they retain the default values from the assim_tools_mod namelist. Bounds don't exist ! in that namelist, so default values are set in assim_tools_mod just before the call to here. @@ -368,6 +376,23 @@ subroutine obs_inc_info(obs_kind, filter_kind, rectangular_quadrature, gaussian_ ! Temporary approach for setting the details of how to assimilate this observation ! This example is designed to reproduce the squared forward operator results from paper +! Fill array with possible filter_kind strings +possible_filter_kinds(1) = 'EAKF' +possible_filter_kinds(2) = 'ENKF' +possible_filter_kinds(3) = 'UNBOUNDED_RHF' +possible_filter_kinds(4) = 'GAMMA_FILTER' +possible_filter_kinds(5) = 'BOUNDED_NORMAL_RHF' +write(*,*) 'possible_filter_kinds' +write(*,*) possible_filter_kinds + +possible_filter_kind_ints(1) = 1 +possible_filter_kind_ints(2) = 2 +possible_filter_kind_ints(3) = 8 +possible_filter_kind_ints(4) = 11 +possible_filter_kind_ints(5) = 101 +write(*,*) 'possible_filter_kind_ints' +write(*,*) possible_filter_kind_ints + !use default values if qcf_table_filename is not in namelist if (.not. qcf_table_listed) then filter_kind = BOUNDED_NORMAL_RHF @@ -394,7 +419,25 @@ subroutine obs_inc_info(obs_kind, filter_kind, rectangular_quadrature, gaussian_ ! Default settings for now for Icepack and tracer model tests (sort_obs_inc, spread_restoration) else - filter_kind = qcf_table_data(QTY_loc(1))%obs_inc_info%filter_kind !filter_kind has a check in obs_increment + + ! Comparing the filter_kind in string format to list of potential filter_kinds + filter_kind_string = qcf_table_data(QTY_loc(1))%obs_inc_info%filter_kind + write(*,*) 'filter_kind_string: ', filter_kind_string + call to_upper(filter_kind_string) + filter_kind_loc = findloc(possible_filter_kinds, trim(filter_kind_string)) + write(*,*) 'filter_kind_string: ', filter_kind_string + + if (filter_kind_loc(1) == 0) then + write(errstring, *) 'Invalid filter_kind' + call error_handler(E_ERR, 'obs_inc_info', errstring, source) + + else + filter_kind = possible_filter_kind_ints(filter_kind_loc(1)) + endif + +! if (filter_kind_string == + + ! filter_kind = qcf_table_data(QTY_loc(1))%obs_inc_info%filter_kind !filter_kind has a check in obs_increment sort_obs_inc = qcf_table_data(QTY_loc(1))%obs_inc_info%sort_obs_inc spread_restoration = qcf_table_data(QTY_loc(1))%obs_inc_info%spread_restoration bounded_below = qcf_table_data(QTY_loc(1))%obs_inc_info%bounded_below @@ -404,7 +447,7 @@ subroutine obs_inc_info(obs_kind, filter_kind, rectangular_quadrature, gaussian_ endif -!write(*,*) 'obs_inc_info: ', filter_kind, sort_obs_inc, spread_restoration, bounded_below, bounded_above, lower_bound, upper_bound +write(*,*) 'obs_inc_info: ', filter_kind, sort_obs_inc, spread_restoration, bounded_below, bounded_above, lower_bound, upper_bound ! Only need to set these two for options the original RHF implementation !!!rectangular_quadrature = .true. From 844a811590ce7c4f3b8a9f9e5028f11da521251e Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Thu, 14 Sep 2023 13:34:53 -0600 Subject: [PATCH 38/60] Reading in dist_type from the table as a character string instead of int --- .../assimilation/algorithm_info_mod.f90 | 99 ++++++++++++++++--- 1 file changed, 86 insertions(+), 13 deletions(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index 42bdca1acf..43cf5a8714 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -48,19 +48,19 @@ module algorithm_info_mod end type type probit_inflation_type - integer :: dist_type + character(len=129) :: dist_type logical :: bounded_below, bounded_above real(r8) :: lower_bound, upper_bound end type type probit_state_type - integer :: dist_type + character(len=129) :: dist_type logical :: bounded_below, bounded_above real(r8) :: lower_bound, upper_bound end type type probit_extended_state_type - integer :: dist_type + character(len=129) :: dist_type logical :: bounded_below, bounded_above real(r8) :: lower_bound, upper_bound end type @@ -270,9 +270,13 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & real(r8), intent(out) :: lower_bound, upper_bound integer :: QTY_loc(1) -character(len=129) :: dist_type_string character(len=129) :: kind_name +integer :: dist_type_loc(1) +character(len=129), dimension(7) :: possible_dist_types +integer, dimension(7) :: possible_dist_type_ints +character(len=129) :: dist_type_string + ! Have input information about the kind of the state or observation being transformed ! along with additional logical info that indicates whether this is an observation ! or state variable and about whether the transformation is being done for inflation @@ -292,6 +296,27 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & ! In the long run, may not have to have separate controls for each of the input possibilities ! However, for now these are things that need to be explored for science understanding +! Fill arrays with possible dist_type strings and corresponding ints +possible_dist_types(1) = 'NORMAL_DISTRIBUTION ' +possible_dist_types(2) = 'BOUNDED_NORMAL_RH_DISTRIBUTION' +possible_dist_types(3) = 'GAMMA_DISTRIBUTION' +possible_dist_types(4) = 'BETA_DISTRIBUTION' +possible_dist_types(5) = 'LOG_NORMAL_DISTRIBUTION' +possible_dist_types(6) = 'UNIFORM_DISTRIBUTION ' +possible_dist_types(7) = 'PARTICLE_FILTER_DISTRIBUTION' +!write(*,*) 'possible_dist_types' +!write(*,*) possible_dist_types + +possible_dist_type_ints(1) = 1 +possible_dist_type_ints(2) = 2 +possible_dist_type_ints(3) = 3 +possible_dist_type_ints(4) = 4 +possible_dist_type_ints(5) = 5 +possible_dist_type_ints(6) = 6 +possible_dist_type_ints(7) = 7 +!write(*,*) 'possible_dist_type_ints' +!write(*,*) possible_dist_type_ints + !use default values if qcf_table_filename is not in namelist if (.not. qcf_table_listed) then dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION @@ -317,32 +342,80 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & elseif(is_inflation) then ! Case for inflation transformation -! dist_type_string = to_upper(qcf_table_data(QTY_loc(1))%probit_inflation%dist_type) + ! Comparing the dist_type in string format to list of potential dist_types + dist_type_string = qcf_table_data(QTY_loc(1))%probit_inflation%dist_type + ! write(*,*) 'dist_type_string: ', dist_type_string + call to_upper(dist_type_string) + dist_type_loc = findloc(possible_dist_types, trim(dist_type_string)) + ! write(*,*) 'dist_type_string: ', dist_type_string + + if (dist_type_loc(1) == 0) then + write(errstring, *) 'Invalid dist_type' + call error_handler(E_ERR, 'probit_dist_info', errstring, source) + + else + dist_type = possible_dist_type_ints(dist_type_loc(1)) + endif - dist_type = qcf_table_data(QTY_loc(1))%probit_inflation%dist_type !dist_type has checks in transform_to_probit, transform_from_probit bounded_below = qcf_table_data(QTY_loc(1))%probit_inflation%bounded_below bounded_above = qcf_table_data(QTY_loc(1))%probit_inflation%bounded_above lower_bound = qcf_table_data(QTY_loc(1))%probit_inflation%lower_bound upper_bound = qcf_table_data(QTY_loc(1))%probit_inflation%upper_bound + write(*,*) 'probit_inflation: ', dist_type, bounded_below, bounded_above, lower_bound, upper_bound + elseif(is_state) then ! Case for state variable priors - dist_type = qcf_table_data(QTY_loc(1))%probit_state%dist_type + ! Comparing the dist_type in string format to list of potential dist_types + dist_type_string = qcf_table_data(QTY_loc(1))%probit_state%dist_type + write(*,*) 'dist_type_string: ', dist_type_string + call to_upper(dist_type_string) + dist_type_loc = findloc(possible_dist_types, trim(dist_type_string)) + write(*,*) 'dist_type_string: ', dist_type_string + + if (dist_type_loc(1) == 0) then + write(errstring, *) 'Invalid dist_type' + call error_handler(E_ERR, 'probit_dist_info', errstring, source) + + else + dist_type = possible_dist_type_ints(dist_type_loc(1)) + endif + + ! dist_type = qcf_table_data(QTY_loc(1))%probit_state%dist_type bounded_below = qcf_table_data(QTY_loc(1))%probit_state%bounded_below bounded_above = qcf_table_data(QTY_loc(1))%probit_state%bounded_above lower_bound = qcf_table_data(QTY_loc(1))%probit_state%lower_bound upper_bound = qcf_table_data(QTY_loc(1))%probit_state%upper_bound + write(*,*) 'probit_state: ', dist_type, bounded_below, bounded_above, lower_bound, upper_bound + else ! This case is for observation (extended state) priors - dist_type = qcf_table_data(QTY_loc(1))%probit_extended_state%dist_type + ! Comparing the dist_type in string format to list of potential dist_types + dist_type_string = qcf_table_data(QTY_loc(1))%probit_extended_state%dist_type + write(*,*) 'dist_type_string: ', dist_type_string + call to_upper(dist_type_string) + dist_type_loc = findloc(possible_dist_types, trim(dist_type_string)) + write(*,*) 'dist_type_string: ', dist_type_string + + if (dist_type_loc(1) == 0) then + write(errstring, *) 'Invalid dist_type' + call error_handler(E_ERR, 'probit_dist_info', errstring, source) + + else + dist_type = possible_dist_type_ints(dist_type_loc(1)) + endif + + ! dist_type = qcf_table_data(QTY_loc(1))%probit_extended_state%dist_type bounded_below = qcf_table_data(QTY_loc(1))%probit_extended_state%bounded_below bounded_above = qcf_table_data(QTY_loc(1))%probit_extended_state%bounded_above lower_bound = qcf_table_data(QTY_loc(1))%probit_extended_state%lower_bound upper_bound = qcf_table_data(QTY_loc(1))%probit_extended_state%upper_bound + write(*,*) 'probit_extended_state: ', dist_type, bounded_below, bounded_above, lower_bound, upper_bound + endif end subroutine probit_dist_info @@ -376,22 +449,22 @@ subroutine obs_inc_info(obs_kind, filter_kind, rectangular_quadrature, gaussian_ ! Temporary approach for setting the details of how to assimilate this observation ! This example is designed to reproduce the squared forward operator results from paper -! Fill array with possible filter_kind strings +! Fill arrays with possible filter_kind strings and corresponding ints possible_filter_kinds(1) = 'EAKF' possible_filter_kinds(2) = 'ENKF' possible_filter_kinds(3) = 'UNBOUNDED_RHF' possible_filter_kinds(4) = 'GAMMA_FILTER' possible_filter_kinds(5) = 'BOUNDED_NORMAL_RHF' -write(*,*) 'possible_filter_kinds' -write(*,*) possible_filter_kinds +!write(*,*) 'possible_filter_kinds' +!write(*,*) possible_filter_kinds possible_filter_kind_ints(1) = 1 possible_filter_kind_ints(2) = 2 possible_filter_kind_ints(3) = 8 possible_filter_kind_ints(4) = 11 possible_filter_kind_ints(5) = 101 -write(*,*) 'possible_filter_kind_ints' -write(*,*) possible_filter_kind_ints +!write(*,*) 'possible_filter_kind_ints' +!write(*,*) possible_filter_kind_ints !use default values if qcf_table_filename is not in namelist if (.not. qcf_table_listed) then From 5ebb7936bfd414906e38c20ce54f68d2adcde490 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Tue, 19 Sep 2023 17:14:57 -0600 Subject: [PATCH 39/60] Fixing log format; cleaning code --- .../assimilation/algorithm_info_mod.f90 | 199 +++++------------- .../modules/assimilation/filter_mod.nml | 1 + 2 files changed, 56 insertions(+), 144 deletions(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index 43cf5a8714..750ecd2ba9 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -117,16 +117,15 @@ subroutine init_algorithm_info_mod(qcf_table_filename) if (qcf_table_filename == '') then write(errstring,*) 'No QCF table file listed in namelist, using default values for all QTYs' - call error_handler(E_MSG, 'init_algorithm_info_mod', errstring, source) + call error_handler(E_MSG, 'init_algorithm_info_mod:', errstring, source) return endif qcf_table_listed = .true. -nlines = 0 - -fileid = open_file(qcf_table_filename, 'formatted', 'read') +fileid = open_file(trim(qcf_table_filename), 'formatted', 'read') -!do loop to get number of rows (or QTY's) in the table +! Do loop to get number of rows (or QTY's) in the table +nlines = 0 do read(fileid,*,iostat=io) if(io/=0) exit @@ -141,12 +140,9 @@ subroutine init_algorithm_info_mod(qcf_table_filename) allocate(qcf_table_row_headers(numrows)) call read_qcf_table(qcf_table_filename) -call write_qcf_table() call assert_qcf_table_version() call verify_qcf_table_data() -!call log_qcf_table_data() - -!stop !FOR TESTING, REMOVE LATER +call log_qcf_table_data() end subroutine init_algorithm_info_mod @@ -164,13 +160,11 @@ subroutine read_qcf_table(qcf_table_filename) if (.not. module_initialized) call init_algorithm_info_mod(qcf_table_filename) -fileid = open_file(qcf_table_filename, 'formatted', 'read') +fileid = open_file(trim(qcf_table_filename), 'formatted', 'read') -! skip the headers, make sure user is using the correct table version +! skip the headers read(fileid, *) header1 read(fileid, *) header2 -!write(*,*) 'header1: ', header1 -!write(*,*) 'header2: ', header2 ! read in table values directly to qcf_table_data type do row = 1, size(qcf_table_data) @@ -304,8 +298,6 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & possible_dist_types(5) = 'LOG_NORMAL_DISTRIBUTION' possible_dist_types(6) = 'UNIFORM_DISTRIBUTION ' possible_dist_types(7) = 'PARTICLE_FILTER_DISTRIBUTION' -!write(*,*) 'possible_dist_types' -!write(*,*) possible_dist_types possible_dist_type_ints(1) = 1 possible_dist_type_ints(2) = 2 @@ -314,8 +306,6 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & possible_dist_type_ints(5) = 5 possible_dist_type_ints(6) = 6 possible_dist_type_ints(7) = 7 -!write(*,*) 'possible_dist_type_ints' -!write(*,*) possible_dist_type_ints !use default values if qcf_table_filename is not in namelist if (.not. qcf_table_listed) then @@ -332,8 +322,6 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & QTY_loc = findloc(qcf_table_row_headers, kind_name) if (QTY_loc(1) == 0) then - ! write(*,*) 'QTY not in table, using default values' !remove these writes on PR - !use default values if QTY is not in table dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION bounded_below = .false.; bounded_above = .false. @@ -344,14 +332,12 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & ! Comparing the dist_type in string format to list of potential dist_types dist_type_string = qcf_table_data(QTY_loc(1))%probit_inflation%dist_type - ! write(*,*) 'dist_type_string: ', dist_type_string call to_upper(dist_type_string) dist_type_loc = findloc(possible_dist_types, trim(dist_type_string)) - ! write(*,*) 'dist_type_string: ', dist_type_string if (dist_type_loc(1) == 0) then - write(errstring, *) 'Invalid dist_type' - call error_handler(E_ERR, 'probit_dist_info', errstring, source) + write(errstring, *) 'Invalid dist_type: ', trim(dist_type_string) + call error_handler(E_ERR, 'probit_dist_info:', errstring, source) else dist_type = possible_dist_type_ints(dist_type_loc(1)) @@ -362,60 +348,48 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & lower_bound = qcf_table_data(QTY_loc(1))%probit_inflation%lower_bound upper_bound = qcf_table_data(QTY_loc(1))%probit_inflation%upper_bound - write(*,*) 'probit_inflation: ', dist_type, bounded_below, bounded_above, lower_bound, upper_bound - elseif(is_state) then ! Case for state variable priors ! Comparing the dist_type in string format to list of potential dist_types dist_type_string = qcf_table_data(QTY_loc(1))%probit_state%dist_type - write(*,*) 'dist_type_string: ', dist_type_string call to_upper(dist_type_string) dist_type_loc = findloc(possible_dist_types, trim(dist_type_string)) - write(*,*) 'dist_type_string: ', dist_type_string if (dist_type_loc(1) == 0) then - write(errstring, *) 'Invalid dist_type' - call error_handler(E_ERR, 'probit_dist_info', errstring, source) + write(errstring, *) 'Invalid dist_type: ', trim(dist_type_string) + call error_handler(E_ERR, 'probit_dist_info:', errstring, source) else dist_type = possible_dist_type_ints(dist_type_loc(1)) endif - ! dist_type = qcf_table_data(QTY_loc(1))%probit_state%dist_type bounded_below = qcf_table_data(QTY_loc(1))%probit_state%bounded_below bounded_above = qcf_table_data(QTY_loc(1))%probit_state%bounded_above lower_bound = qcf_table_data(QTY_loc(1))%probit_state%lower_bound upper_bound = qcf_table_data(QTY_loc(1))%probit_state%upper_bound - write(*,*) 'probit_state: ', dist_type, bounded_below, bounded_above, lower_bound, upper_bound - else ! This case is for observation (extended state) priors ! Comparing the dist_type in string format to list of potential dist_types dist_type_string = qcf_table_data(QTY_loc(1))%probit_extended_state%dist_type - write(*,*) 'dist_type_string: ', dist_type_string call to_upper(dist_type_string) dist_type_loc = findloc(possible_dist_types, trim(dist_type_string)) - write(*,*) 'dist_type_string: ', dist_type_string if (dist_type_loc(1) == 0) then - write(errstring, *) 'Invalid dist_type' - call error_handler(E_ERR, 'probit_dist_info', errstring, source) + write(errstring, *) 'Invalid dist_type: ', trim(dist_type_string) + call error_handler(E_ERR, 'probit_dist_info:', errstring, source) else dist_type = possible_dist_type_ints(dist_type_loc(1)) endif - ! dist_type = qcf_table_data(QTY_loc(1))%probit_extended_state%dist_type bounded_below = qcf_table_data(QTY_loc(1))%probit_extended_state%bounded_below bounded_above = qcf_table_data(QTY_loc(1))%probit_extended_state%bounded_above lower_bound = qcf_table_data(QTY_loc(1))%probit_extended_state%lower_bound upper_bound = qcf_table_data(QTY_loc(1))%probit_extended_state%upper_bound - write(*,*) 'probit_extended_state: ', dist_type, bounded_below, bounded_above, lower_bound, upper_bound - endif end subroutine probit_dist_info @@ -455,16 +429,12 @@ subroutine obs_inc_info(obs_kind, filter_kind, rectangular_quadrature, gaussian_ possible_filter_kinds(3) = 'UNBOUNDED_RHF' possible_filter_kinds(4) = 'GAMMA_FILTER' possible_filter_kinds(5) = 'BOUNDED_NORMAL_RHF' -!write(*,*) 'possible_filter_kinds' -!write(*,*) possible_filter_kinds possible_filter_kind_ints(1) = 1 possible_filter_kind_ints(2) = 2 possible_filter_kind_ints(3) = 8 possible_filter_kind_ints(4) = 11 possible_filter_kind_ints(5) = 101 -!write(*,*) 'possible_filter_kind_ints' -!write(*,*) possible_filter_kind_ints !use default values if qcf_table_filename is not in namelist if (.not. qcf_table_listed) then @@ -482,8 +452,6 @@ subroutine obs_inc_info(obs_kind, filter_kind, rectangular_quadrature, gaussian_ QTY_loc = findloc(qcf_table_row_headers, kind_name) if (QTY_loc(1) == 0) then - !write(*,*) 'QTY not in table, using default values' - !use default values if QTY is not in table filter_kind = BOUNDED_NORMAL_RHF bounded_below = .false.; bounded_above = .false. @@ -495,22 +463,17 @@ subroutine obs_inc_info(obs_kind, filter_kind, rectangular_quadrature, gaussian_ ! Comparing the filter_kind in string format to list of potential filter_kinds filter_kind_string = qcf_table_data(QTY_loc(1))%obs_inc_info%filter_kind - write(*,*) 'filter_kind_string: ', filter_kind_string call to_upper(filter_kind_string) filter_kind_loc = findloc(possible_filter_kinds, trim(filter_kind_string)) - write(*,*) 'filter_kind_string: ', filter_kind_string if (filter_kind_loc(1) == 0) then - write(errstring, *) 'Invalid filter_kind' - call error_handler(E_ERR, 'obs_inc_info', errstring, source) + write(errstring, *) 'Invalid filter_kind: ', trim(filter_kind_string) + call error_handler(E_ERR, 'obs_inc_info:', errstring, source) else filter_kind = possible_filter_kind_ints(filter_kind_loc(1)) endif -! if (filter_kind_string == - - ! filter_kind = qcf_table_data(QTY_loc(1))%obs_inc_info%filter_kind !filter_kind has a check in obs_increment sort_obs_inc = qcf_table_data(QTY_loc(1))%obs_inc_info%sort_obs_inc spread_restoration = qcf_table_data(QTY_loc(1))%obs_inc_info%spread_restoration bounded_below = qcf_table_data(QTY_loc(1))%obs_inc_info%bounded_below @@ -520,8 +483,6 @@ subroutine obs_inc_info(obs_kind, filter_kind, rectangular_quadrature, gaussian_ endif -write(*,*) 'obs_inc_info: ', filter_kind, sort_obs_inc, spread_restoration, bounded_below, bounded_above, lower_bound, upper_bound - ! Only need to set these two for options the original RHF implementation !!!rectangular_quadrature = .true. !!!gaussian_likelihood_tails = .false. @@ -531,55 +492,13 @@ end subroutine obs_inc_info !------------------------------------------------------------------------ -subroutine write_qcf_table() - -! DRAFT SUBROUTINE -! write to check values were correctly assigned -! testing for findloc - -character(len=30), parameter :: tester_QTY = 'QTY_STATE_VARIABLE' -integer :: QTY_loc(1) - -character(len=30), parameter :: tester_QTY0 = 'QTY_DUMMY' -integer :: QTY_loc0(1) - -integer :: row - -do row = 1, size(qcf_table_data) - write(*,*) 'qcf_table_row_headers(', row, '): ', qcf_table_row_headers(row) - write(*,*) 'qcf_table_data(', row, '): ' - write(*,*) qcf_table_data(row)%obs_error_info%bounded_below, qcf_table_data(row)%obs_error_info%bounded_above, & - qcf_table_data(row)%obs_error_info%lower_bound, qcf_table_data(row)%obs_error_info%upper_bound, qcf_table_data(row)%probit_inflation%dist_type, & - qcf_table_data(row)%probit_inflation%bounded_below, qcf_table_data(row)%probit_inflation%bounded_above, & - qcf_table_data(row)%probit_inflation%lower_bound, qcf_table_data(row)%probit_inflation%upper_bound, qcf_table_data(row)%probit_state%dist_type, & - qcf_table_data(row)%probit_state%bounded_below, qcf_table_data(row)%probit_state%bounded_above, & - qcf_table_data(row)%probit_state%lower_bound, qcf_table_data(row)%probit_state%upper_bound, qcf_table_data(row)%probit_extended_state%dist_type, & - qcf_table_data(row)%probit_extended_state%bounded_below, qcf_table_data(row)%probit_extended_state%bounded_above, & - qcf_table_data(row)%probit_extended_state%lower_bound, qcf_table_data(row)%probit_extended_state%upper_bound, & - qcf_table_data(row)%obs_inc_info%filter_kind, qcf_table_data(row)%obs_inc_info%rectangular_quadrature, & - qcf_table_data(row)%obs_inc_info%gaussian_likelihood_tails, qcf_table_data(row)%obs_inc_info%sort_obs_inc, & - qcf_table_data(row)%obs_inc_info%spread_restoration, qcf_table_data(row)%obs_inc_info%bounded_below, qcf_table_data(row)%obs_inc_info%bounded_above, & - qcf_table_data(row)%obs_inc_info%lower_bound, qcf_table_data(row)%obs_inc_info%upper_bound -end do - -QTY_loc = findloc(qcf_table_row_headers, tester_QTY) -write(*, *) 'findloc of QTY_STATE_VARIABLE: ', QTY_loc(1) - -QTY_loc0 = findloc(qcf_table_row_headers, tester_QTY0) -write(*, *) 'findloc of invalid QTY (QTY_DUMMY): ', QTY_loc0(1) - -end subroutine write_qcf_table - -!------------------------------------------------------------------------ - - subroutine assert_qcf_table_version() -!subroutine to ensure the correct version of the QCF table is being used +! Subroutine to ensure the correct version of the QCF table is being used -if (header1(4) /= '1') then +if (trim(header1(4)) /= '1') then write(errstring,*) 'Using outdated/incorrect version of the QCF table' - call error_handler(E_ERR, 'assert_qcf_table_version', errstring, source) + call error_handler(E_ERR, 'assert_qcf_table_version:', errstring, source) endif end subroutine assert_qcf_table_version @@ -589,13 +508,11 @@ end subroutine assert_qcf_table_version subroutine verify_qcf_table_data() -!subroutine to ensure that the data in the QCF table is valid +! Subroutine to ensure that the data in the QCF table is valid integer :: varid integer :: row -!if (.not. module_initialized) call init_algorithm_info_mod(qcf_table_filename) - if (.not. qcf_table_listed) return !Checks that all bounds are valid; currently checks that the lower bound in less than the upper @@ -603,40 +520,40 @@ subroutine verify_qcf_table_data() do row = 1, size(qcf_table_data) if(qcf_table_data(row)%obs_error_info%lower_bound > qcf_table_data(row)%obs_error_info%upper_bound) then write(errstring,*) 'Invalid bounds in obs_error_info' - call error_handler(E_ERR, 'verify_qcf_table_data', errstring, source) + call error_handler(E_ERR, 'verify_qcf_table_data:', errstring, source) endif if(qcf_table_data(row)%probit_inflation%lower_bound > qcf_table_data(row)%probit_inflation%upper_bound) then write(errstring,*) 'Invalid bounds in probit_inflation' - call error_handler(E_ERR, 'verify_qcf_table_data', errstring, source) + call error_handler(E_ERR, 'verify_qcf_table_data:', errstring, source) endif if(qcf_table_data(row)%probit_state%lower_bound > qcf_table_data(row)%probit_state%upper_bound) then write(errstring,*) 'Invalid bounds in probit_state' - call error_handler(E_ERR, 'verify_qcf_table_data', errstring, source) + call error_handler(E_ERR, 'verify_qcf_table_data:', errstring, source) endif if(qcf_table_data(row)%probit_extended_state%lower_bound > qcf_table_data(row)%probit_extended_state%upper_bound) then write(errstring,*) 'Invalid bounds in probit_extended_state' - call error_handler(E_ERR, 'verify_qcf_table_data', errstring, source) + call error_handler(E_ERR, 'verify_qcf_table_data:', errstring, source) endif if(qcf_table_data(row)%obs_inc_info%lower_bound > qcf_table_data(row)%obs_inc_info%upper_bound) then write(errstring,*) 'Invalid bounds in obs_inc_info' - call error_handler(E_ERR, 'verify_qcf_table_data', errstring, source) + call error_handler(E_ERR, 'verify_qcf_table_data:', errstring, source) endif end do !Ensures that all QTYs listed in the table exist in DART do row = 1, size(qcf_table_data) - varid = get_index_for_quantity(qcf_table_row_headers(row)) + varid = get_index_for_quantity(trim(qcf_table_row_headers(row))) if(varid == -1) then write(errstring,*) trim(qcf_table_row_headers(row)), ' is not a valid DART QTY' - call error_handler(E_ERR, 'verify_qcf_table_data', errstring, source) + call error_handler(E_ERR, 'verify_qcf_table_data:', errstring, source) endif end do !Ensures that there are no duplicate QTYs in the table do row = 1, size(qcf_table_data) - if(count(qcf_table_row_headers==qcf_table_row_headers(row)) > 1) then + if(count(qcf_table_row_headers==trim(qcf_table_row_headers(row))) > 1) then write(errstring,*) trim(qcf_table_row_headers(row)), ' has multiple entries in the table' - call error_handler(E_ERR, 'verify_qcf_table_data', errstring, source) + call error_handler(E_ERR, 'verify_qcf_table_data:', errstring, source) endif end do @@ -647,46 +564,40 @@ end subroutine verify_qcf_table_data subroutine log_qcf_table_data() -!subroutine to write the data in QCF table to dart_log - -character(len=512) :: log_msg +! Subroutine to write the data in QCF table to dart_log +character(len=2000) :: log_msg integer :: row -integer :: i if (.not. qcf_table_listed) return -!call error_handler(E_ALLMSG, 'log_qcf_table_data', log_msg, source) -!call log_it(log_msg) +call error_handler(E_MSG, '', '', source) !Writing blank line to log +call error_handler(E_MSG, 'log_qcf_table_data:', 'Logging the data in the QCF Table', source) + +! Write the table headers to the dart_log and terminal +write(log_msg, '(A4, A6, A9, A)') header1(:) +call error_handler(E_MSG, 'log_qcf_table_data:', trim(log_msg), source) -!Write the headers to the dart_log and terminal +write(log_msg,'(3A14, 2A12, 3(A10, 2A14, 2A12), A12, A23, A26, A13, A19, 2A14, 2A12)') header2(:) +call error_handler(E_MSG, 'log_qcf_table_data:', trim(log_msg), source) -!write(log_msg, *) trim(header1(1)) -do i = 1, size(header1) - write(*,*) trim(header1(i)) - write(log_msg,*) trim(header1(i)) - !log_msg = log_msg//trim(header1(i)) - ! call error_handler(E_MSG, 'log_qcf_table_data', trim(log_msg), source) - write(*, *) trim(log_msg) +! Write the table data to the dart_log and terminal +do row = 1, size(qcf_table_data) + write(log_msg, *) trim(qcf_table_row_headers(row)), qcf_table_data(row)%obs_error_info%bounded_below, qcf_table_data(row)%obs_error_info%bounded_above, & + qcf_table_data(row)%obs_error_info%lower_bound, qcf_table_data(row)%obs_error_info%upper_bound, trim(qcf_table_data(row)%probit_inflation%dist_type), & + qcf_table_data(row)%probit_inflation%bounded_below, qcf_table_data(row)%probit_inflation%bounded_above, & + qcf_table_data(row)%probit_inflation%lower_bound, qcf_table_data(row)%probit_inflation%upper_bound, trim(qcf_table_data(row)%probit_state%dist_type), & + qcf_table_data(row)%probit_state%bounded_below, qcf_table_data(row)%probit_state%bounded_above, & + qcf_table_data(row)%probit_state%lower_bound, qcf_table_data(row)%probit_state%upper_bound, trim(qcf_table_data(row)%probit_extended_state%dist_type), & + qcf_table_data(row)%probit_extended_state%bounded_below, qcf_table_data(row)%probit_extended_state%bounded_above, & + qcf_table_data(row)%probit_extended_state%lower_bound, qcf_table_data(row)%probit_extended_state%upper_bound, & + trim(qcf_table_data(row)%obs_inc_info%filter_kind), qcf_table_data(row)%obs_inc_info%rectangular_quadrature, & + qcf_table_data(row)%obs_inc_info%gaussian_likelihood_tails, qcf_table_data(row)%obs_inc_info%sort_obs_inc, & + qcf_table_data(row)%obs_inc_info%spread_restoration, qcf_table_data(row)%obs_inc_info%bounded_below, qcf_table_data(row)%obs_inc_info%bounded_above, & + qcf_table_data(row)%obs_inc_info%lower_bound, qcf_table_data(row)%obs_inc_info%upper_bound +call error_handler(E_MSG, 'log_qcf_table_data:', trim(log_msg), source) end do -write(*,*) 'log_msg: ', trim(log_msg) -!write(*,*) trim(log_msg) -!write(log_msg,*) header1 -!call error_handler(E_MSG, 'log_qcf_table_data', trim(log_msg), source) -!write(log_msg, *) header2 -!call error_handler(E_MSG, 'log_qcf_table_data', trim(log_msg), source) - -!Write the data to the dart_log and terminal -!do row = 1, size(qcf_table_data) -! write(log_msg,*) qcf_table_row_headers(row), qcf_table_data(row) -! call error_handler(E_MSG, 'log_qcf_table_data', trim(log_msg), source) -!end do - -write(log_msg,*) qcf_table_data -write(*,*) 'log_msg: ', trim(log_msg) -!write(*, *) trim(log_msg) -!write(logfileunit, *) trim(log_msg) -!call log_it(trim(log_msg)) -!call error_handler(E_MSG, 'log_qcf_table_data', trim(log_msg), source) + +call error_handler(E_MSG, '', '', source) !Writing blank line to log end subroutine log_qcf_table_data diff --git a/assimilation_code/modules/assimilation/filter_mod.nml b/assimilation_code/modules/assimilation/filter_mod.nml index 362138fc5f..79611257bf 100644 --- a/assimilation_code/modules/assimilation/filter_mod.nml +++ b/assimilation_code/modules/assimilation/filter_mod.nml @@ -1,4 +1,5 @@ &filter_nml + qcf_table_filename = '' use_algorithm_info_mod = .true., single_file_in = .false., input_state_files = '' From 2be50efce04e061ac0bc69cd2f484919752afc69 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Mon, 25 Sep 2023 11:09:02 -0600 Subject: [PATCH 40/60] Removing unnecessary files --- .../modules/assimilation/cam_qcf_table.txt | 10 -- qcf_table/type_read_table.f90 | 129 ------------------ 2 files changed, 139 deletions(-) delete mode 100644 assimilation_code/modules/assimilation/cam_qcf_table.txt delete mode 100644 qcf_table/type_read_table.f90 diff --git a/assimilation_code/modules/assimilation/cam_qcf_table.txt b/assimilation_code/modules/assimilation/cam_qcf_table.txt deleted file mode 100644 index 56e205a534..0000000000 --- a/assimilation_code/modules/assimilation/cam_qcf_table.txt +++ /dev/null @@ -1,10 +0,0 @@ -QCF table version 1: -QTY bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound filter_kind rectangular_quadrature gaussian_likelihood_tails sort_obs_inc spread_restoration bounded_below bounded_above lower_bound upper_bound -QTY_U_WIND_COMPONENT, .false. .false. 3 4 5 .false. .false. 8 9 10 .false. .false. 13 14 15 .false. .false. 18 19 20 .false. .false. .false. .false. .false. .false. 27 28 -QTY_V_WIND_COMPONENT, .false. .false. 3 4 5 .false. .false. 8 9 10 .false. .false. 13 14 15 .false. .false. 18 19 20 .false. .false. .false. .false. .false. .false. 27 28 -QTY_SURFACE_PRESSURE .false. .false. 3 4 5 .false. .false. 8 9 10 .false. .false. 13 14 15 .false. .false. 18 19 20 .false. .false. .false. .false. .false. .false. 27 28 -QTY_TEMPERATURE .false. .false. 3 4 5 .false. .false. 8 9 10 .false. .false. 13 14 15 .false. .false. 18 19 20 .false. .false. .false. .false. .false. .false. 27 28 -QTY_SPECIFIC_HUMIDITY .false. .false. 3 4 5 .false. .false. 8 9 10 .false. .false. 13 14 15 .false. .false. 18 19 20 .false. .false. .false. .false. .false. .false. 27 28 -QTY_CLOUD_LIQUID_WATER .false. .false. 3 4 5 .false. .false. 8 9 10 .false. .false. 13 14 15 .false. .false. 18 19 20 .false. .false. .false. .false. .false. .false. 27 28 -QTY_CLOUD_ICE .false. .false. 3 4 5 .false. .false. 8 9 10 .false. .false. 13 14 15 .false. .false. 18 19 20 .false. .false. .false. .false. .false. .false. 27 28 -QTY_GPSRO .false. .false. 3 4 5 .false. .false. 8 9 10 .false. .false. 13 14 15 .false. .false. 18 19 20 .false. .false. .false. .false. .false. .false. 27 28 diff --git a/qcf_table/type_read_table.f90 b/qcf_table/type_read_table.f90 deleted file mode 100644 index f7cd548acf..0000000000 --- a/qcf_table/type_read_table.f90 +++ /dev/null @@ -1,129 +0,0 @@ -program read_table - -implicit none -type obs_error_info_type - logical :: bounded_below, bounded_above - real :: lower_bound, upper_bound -end type - -type probit_inflation_type - integer :: dist_type - logical :: bounded_below, bounded_above - real :: lower_bound, upper_bound -end type - -type probit_state_type - integer :: dist_type - logical :: bounded_below, bounded_above - real :: lower_bound, upper_bound -end type - -type probit_extended_state_type - integer :: dist_type - logical :: bounded_below, bounded_above - real :: lower_bound, upper_bound -end type - -type obs_inc_info_type - integer :: filter_kind - logical :: rectangular_quadrature, gaussian_likelihood_tails - logical :: sort_obs_inc, spread_restoration - logical :: bounded_below, bounded_above - real :: lower_bound, upper_bound -end type - -type qcf_table_data_type - type(obs_error_info_type) :: obs_error_info - type(probit_inflation_type) :: probit_inflation - type(probit_state_type) :: probit_state - type(probit_extended_state_type) :: probit_extended_state - type(obs_inc_info_type) :: obs_inc_info -end type - -! Reads in the QCEFF input options from tabular data file -!character(len=50), intent(in) :: qcf_table_filename -!real(r8), intent(out) :: qcf_table_data -!real, dimension(:, :), allocatable :: qcf_table_data_rows -type(qcf_table_data_type), allocatable :: qcf_table_data(:) -character(len=129), dimension(:), allocatable :: rowheaders !!!!! might need to change len=30 - -integer, parameter :: fileid = 10 !file identifier -character(len=30), parameter :: tester_QTY = 'QTY_GPSRO' -integer :: QTY_loc(1) - -!integer, parameter :: num_columns = 28 -integer :: nlines -integer :: io -integer :: numrows -integer :: row - -!real, dimension(1:num_columns, 1:num_rows) :: table_data -!integer :: table_data_1, table_data_2 -character(len=30), dimension(4) :: header1 -character(len=30), dimension(29) :: header2 -!variables for table values ^^^ - -open(unit=fileid, file='cam_qcf_table.txt') -nlines = 0 - -do !do loop to get number of rows (or QTY's) in the table - read(fileid,*,iostat=io) - if(io/=0) exit - nlines = nlines + 1 -end do -close(fileid) - -print*, nlines - -numrows = nlines - 2 -print *, 'numrows: ', numrows - -allocate(qcf_table_data(numrows)) -allocate(rowheaders(numrows)) -write(*,*) shape(qcf_table_data) - -open(unit=fileid, file='cam_qcf_table.txt') - -read(fileid, *) header1 -read(fileid, *) header2 !! skip the headers -Write(*, *) "header1: ", header1 -Write(*, *) "header2: ", header2 - -do row = 1, numrows - read(fileid, *) rowheaders(row), qcf_table_data(row)%obs_error_info%bounded_below, qcf_table_data(row)%obs_error_info%bounded_above, & - qcf_table_data(row)%obs_error_info%lower_bound, qcf_table_data(row)%obs_error_info%upper_bound, qcf_table_data(row)%probit_inflation%dist_type, & - qcf_table_data(row)%probit_inflation%bounded_below, qcf_table_data(row)%probit_inflation%bounded_above, & - qcf_table_data(row)%probit_inflation%lower_bound, qcf_table_data(row)%probit_inflation%upper_bound, qcf_table_data(row)%probit_state%dist_type, & - qcf_table_data(row)%probit_state%bounded_below, qcf_table_data(row)%probit_state%bounded_above, & - qcf_table_data(row)%probit_state%lower_bound, qcf_table_data(row)%probit_state%upper_bound, qcf_table_data(row)%probit_extended_state%dist_type, & - qcf_table_data(row)%probit_extended_state%bounded_below, qcf_table_data(row)%probit_extended_state%bounded_above, & - qcf_table_data(row)%probit_extended_state%lower_bound, qcf_table_data(row)%probit_extended_state%upper_bound, & - qcf_table_data(row)%obs_inc_info%filter_kind, qcf_table_data(row)%obs_inc_info%rectangular_quadrature, & - qcf_table_data(row)%obs_inc_info%gaussian_likelihood_tails, qcf_table_data(row)%obs_inc_info%sort_obs_inc, & - qcf_table_data(row)%obs_inc_info%spread_restoration, qcf_table_data(row)%obs_inc_info%bounded_below, qcf_table_data(row)%obs_inc_info%bounded_above, & - qcf_table_data(row)%obs_inc_info%lower_bound, qcf_table_data(row)%obs_inc_info%upper_bound - - write(*, *) "rowheader(", row, "): ", rowheaders(row) - write(*, *) "qcf_table_data(", row, "): " - write(*, *) qcf_table_data(row)%obs_error_info%bounded_below, qcf_table_data(row)%obs_error_info%bounded_above, & - qcf_table_data(row)%obs_error_info%lower_bound, qcf_table_data(row)%obs_error_info%upper_bound, qcf_table_data(row)%probit_inflation%dist_type, & - qcf_table_data(row)%probit_inflation%bounded_below, qcf_table_data(row)%probit_inflation%bounded_above, & - qcf_table_data(row)%probit_inflation%lower_bound, qcf_table_data(row)%probit_inflation%upper_bound, qcf_table_data(row)%probit_state%dist_type, & - qcf_table_data(row)%probit_state%bounded_below, qcf_table_data(row)%probit_state%bounded_above, & - qcf_table_data(row)%probit_state%lower_bound, qcf_table_data(row)%probit_state%upper_bound, qcf_table_data(row)%probit_extended_state%dist_type, & - qcf_table_data(row)%probit_extended_state%bounded_below, qcf_table_data(row)%probit_extended_state%bounded_above, & - qcf_table_data(row)%probit_extended_state%lower_bound, qcf_table_data(row)%probit_extended_state%upper_bound, & - qcf_table_data(row)%obs_inc_info%filter_kind, qcf_table_data(row)%obs_inc_info%rectangular_quadrature, & - qcf_table_data(row)%obs_inc_info%gaussian_likelihood_tails, qcf_table_data(row)%obs_inc_info%sort_obs_inc, & - qcf_table_data(row)%obs_inc_info%spread_restoration, qcf_table_data(row)%obs_inc_info%bounded_below, qcf_table_data(row)%obs_inc_info%bounded_above, & - qcf_table_data(row)%obs_inc_info%lower_bound, qcf_table_data(row)%obs_inc_info%upper_bound -end do - -close(fileid) - -QTY_loc = findloc(rowheaders, tester_QTY) -write(*, *) 'findloc of GPSRO: ', QTY_loc(1) - -deallocate(qcf_table_data, rowheaders) - -end program read_table From 399b62bfa5ab456ac1e83181d015df16bdd15524 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Mon, 25 Sep 2023 11:15:38 -0600 Subject: [PATCH 41/60] Adding the yaml_to_table.py script and example yaml file qcf_table_template.yaml to the DART repo at /DART/assimilation_code/programs/qcf_table/. The location of these files in the repo may change --- .../qcf_table/qcf_table_template.yaml | 103 ++++++++++++++++++ .../programs/qcf_table/yaml_to_table.py | 77 +++++++++++++ 2 files changed, 180 insertions(+) create mode 100644 assimilation_code/programs/qcf_table/qcf_table_template.yaml create mode 100644 assimilation_code/programs/qcf_table/yaml_to_table.py diff --git a/assimilation_code/programs/qcf_table/qcf_table_template.yaml b/assimilation_code/programs/qcf_table/qcf_table_template.yaml new file mode 100644 index 0000000000..c2a62d9eb5 --- /dev/null +++ b/assimilation_code/programs/qcf_table/qcf_table_template.yaml @@ -0,0 +1,103 @@ +QCF table version: 1 +QTY_TEMPLATE: + obs_error_info: + bounded_below + bounded_above + lower_bound + upper_bound + probit_inflation: + dist_type + bounded_below + bounded_above + lower_bound + upper_bound + probit_state: + dist_type + bounded_below + bounded_above + lower_bound + upper_bound + probit_extended_state: + dist_type + bounded_below + bounded_above + lower_bound + upper_bound + obs_inc_info: + filter_kind + rectangular_quadrature + gaussian_likelihood_tails + sort_obs_inc + spread_restoration + bounded_below + bounded_above + lower_bound + upper_bound +QTY_STATE_VARIABLE: + obs_error_info: + bounded_below: .false. + bounded_above: .false. + lower_bound: -888888.0 + upper_bound: -888888.0 + probit_inflation: + dist_type: BOUNDED_NORMAL_RH_DISTRIBUTION + bounded_below: .false. + bounded_above: .false. + lower_bound: -888888.0 + upper_bound: -888888.0 + probit_state: + dist_type: BOUNDED_NORMAL_RH_DISTRIBUTION + bounded_below: .false. + bounded_above: .false. + lower_bound: -888888.0 + upper_bound: -888888.0 + probit_extended_state: + dist_type: BOUNDED_NORMAL_RH_DISTRIBUTION + bounded_below: .false. + bounded_above: .false. + lower_bound: -888888.0 + upper_bound: -888888.0 + obs_inc_info: + filter_kind: BOUNDED_NORMAL_RHF + rectangular_quadrature: .false. + gaussian_likelihood_tails: .false. + sort_obs_inc: .false. + spread_restoration: .false. + bounded_below: .false. + bounded_above: .false. + lower_bound: -888888.0 + upper_bound: -888888.0 +QTY_TRACER_SOURCE: + obs_error_info: + bounded_below: .false. + bounded_above: .false. + lower_bound: -888888.0 + upper_bound: -888888.0 + probit_inflation: + dist_type: BOUNDED_NORMAL_RH_DISTRIBUTION + bounded_below: .false. + bounded_above: .false. + lower_bound: -888888.0 + upper_bound: -888888.0 + probit_state: + dist_type: BOUNDED_NORMAL_RH_DISTRIBUTION + bounded_below: .false. + bounded_above: .false. + lower_bound: -888888.0 + upper_bound: -888888.0 + probit_extended_state: + dist_type: BOUNDED_NORMAL_RH_DISTRIBUTION + bounded_below: .false. + bounded_above: .false. + lower_bound: -888888.0 + upper_bound: -888888.0 + obs_inc_info: + filter_kind: BOUNDED_NORMAL_RHF + rectangular_quadrature: .false. + gaussian_likelihood_tails: .false. + sort_obs_inc: .false. + spread_restoration: .false. + bounded_below: .false. + bounded_above: .false. + lower_bound: -888888.0 + upper_bound: -888888.0 diff --git a/assimilation_code/programs/qcf_table/yaml_to_table.py b/assimilation_code/programs/qcf_table/yaml_to_table.py new file mode 100644 index 0000000000..fc765e7c36 --- /dev/null +++ b/assimilation_code/programs/qcf_table/yaml_to_table.py @@ -0,0 +1,77 @@ +import yaml + +#Prompt user for name of input and output files +input_yaml = input('Please enter the name of your input yaml file (filename must end in ".yaml") OR press enter/return to use the default filename "qcf_table.yaml"\n') +output_txt = input('Please enter the name for the output file for the table (filename must end in ".txt") OR press enter/return to use the default filename "qcf_table.txt"\n') + +#Using deault names for input/output files if not specified +if (input_yaml == ''): + input_yaml = 'qcf_table.yaml' + +if (output_txt == ''): + output_txt = 'qcf_table.txt' + +#Open and load yaml file +with open(input_yaml) as file: + dict = yaml.safe_load(file) + + column_headers = list(dict.keys()) + column_data = list(dict.values()) + + obs_errror_info_header = dict['QTY_TEMPLATE']['obs_error_info'] + probit_inflation_header = dict['QTY_TEMPLATE']['probit_inflation'] + probit_state_header = dict['QTY_TEMPLATE']['probit_state'] + probit_extended_state_header = dict['QTY_TEMPLATE']['probit_extended_state'] + obs_inc_info_header = dict['QTY_TEMPLATE']['obs_inc_info'] + + f = open(output_txt, "w") + +#Write the table's headers to the output file + f.write(column_headers[0] + ": " + str(column_data[0]) + "\n") + + f.write(column_headers[1] + ": ") + for name in obs_errror_info_header: + f.write(name) + f.write(" ") + for name in probit_inflation_header: + f.write(name) + f.write(" ") + for name in probit_state_header: + f.write(name) + f.write(" ") + for name in probit_extended_state_header: + f.write(name) + f.write(" ") + for name in obs_inc_info_header: + f.write(name) + f.write("\n") + +#Writing table data to the output file + for i in range(2, len(column_headers)): + f.write(column_headers[i] + " ") + + obs_error_info = dict[column_headers[i]]['obs_error_info'].items() + for key, value in obs_error_info: + f.write(str(value) + " ") + + probit_inflation = dict[column_headers[i]]['probit_inflation'].items() + for key, value in probit_inflation: + f.write(str(value) + " ") + + probit_state = dict[column_headers[i]]['probit_state'].items() + for key, value in probit_state: + f.write(str(value) + " ") + + probit_extended_state = dict[column_headers[i]]['probit_extended_state'].items() + for key, value in probit_extended_state: + f.write(str(value) + " ") + + obs_inc_info = dict[column_headers[i]]['obs_inc_info'].items() + for key, value in obs_inc_info: + f.write(str(value) + " ") + + f.write("\n") + + f.close + +print('QCF table produced in ' + output_txt) From 1911ed03cd152329877c03235e2c8d76d0f2b8c4 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Mon, 25 Sep 2023 11:20:49 -0600 Subject: [PATCH 42/60] Fixing deallocation by adding conditional to check if data types were allocated --- .../modules/assimilation/algorithm_info_mod.f90 | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index 750ecd2ba9..47c8f3f844 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -607,12 +607,13 @@ end subroutine log_qcf_table_data subroutine end_algorithm_info_mod() if (.not. module_initialized) return +module_initialized = .false. + +if (.not. qcf_table_listed) return deallocate(qcf_table_data) deallocate(qcf_table_row_headers) -module_initialized = .false. - end subroutine end_algorithm_info_mod !---------------------------------------------------------------------- From 41ffe9554d68a5be0f6fcd95351e2758d75ad342 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Mon, 25 Sep 2023 15:42:22 -0600 Subject: [PATCH 43/60] Adjusting the YAML template file; removing trailing spaces from possible_dist_types strings --- .../assimilation/algorithm_info_mod.f90 | 4 +-- .../qcf_table/qcf_table_template.yaml | 34 ------------------- 2 files changed, 2 insertions(+), 36 deletions(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index 47c8f3f844..dbe9027985 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -291,12 +291,12 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & ! However, for now these are things that need to be explored for science understanding ! Fill arrays with possible dist_type strings and corresponding ints -possible_dist_types(1) = 'NORMAL_DISTRIBUTION ' +possible_dist_types(1) = 'NORMAL_DISTRIBUTION' possible_dist_types(2) = 'BOUNDED_NORMAL_RH_DISTRIBUTION' possible_dist_types(3) = 'GAMMA_DISTRIBUTION' possible_dist_types(4) = 'BETA_DISTRIBUTION' possible_dist_types(5) = 'LOG_NORMAL_DISTRIBUTION' -possible_dist_types(6) = 'UNIFORM_DISTRIBUTION ' +possible_dist_types(6) = 'UNIFORM_DISTRIBUTION' possible_dist_types(7) = 'PARTICLE_FILTER_DISTRIBUTION' possible_dist_type_ints(1) = 1 diff --git a/assimilation_code/programs/qcf_table/qcf_table_template.yaml b/assimilation_code/programs/qcf_table/qcf_table_template.yaml index c2a62d9eb5..46b46358af 100644 --- a/assimilation_code/programs/qcf_table/qcf_table_template.yaml +++ b/assimilation_code/programs/qcf_table/qcf_table_template.yaml @@ -67,37 +67,3 @@ QTY_STATE_VARIABLE: bounded_above: .false. lower_bound: -888888.0 upper_bound: -888888.0 -QTY_TRACER_SOURCE: - obs_error_info: - bounded_below: .false. - bounded_above: .false. - lower_bound: -888888.0 - upper_bound: -888888.0 - probit_inflation: - dist_type: BOUNDED_NORMAL_RH_DISTRIBUTION - bounded_below: .false. - bounded_above: .false. - lower_bound: -888888.0 - upper_bound: -888888.0 - probit_state: - dist_type: BOUNDED_NORMAL_RH_DISTRIBUTION - bounded_below: .false. - bounded_above: .false. - lower_bound: -888888.0 - upper_bound: -888888.0 - probit_extended_state: - dist_type: BOUNDED_NORMAL_RH_DISTRIBUTION - bounded_below: .false. - bounded_above: .false. - lower_bound: -888888.0 - upper_bound: -888888.0 - obs_inc_info: - filter_kind: BOUNDED_NORMAL_RHF - rectangular_quadrature: .false. - gaussian_likelihood_tails: .false. - sort_obs_inc: .false. - spread_restoration: .false. - bounded_below: .false. - bounded_above: .false. - lower_bound: -888888.0 - upper_bound: -888888.0 From f5fa29bb0c604c989a568c42adac0e09f1354906 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Mon, 25 Sep 2023 17:28:21 -0600 Subject: [PATCH 44/60] Adding documentaion file to repo --- guide/qcf_table.rst | 227 ++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 227 insertions(+) create mode 100644 guide/qcf_table.rst diff --git a/guide/qcf_table.rst b/guide/qcf_table.rst new file mode 100644 index 0000000000..113bffee9a --- /dev/null +++ b/guide/qcf_table.rst @@ -0,0 +1,227 @@ +.. _QCF Table: + +############################################ +Using the QCF Table to Control Input Options +############################################ + +This file contains instructions for using an input table to set input options with the DART quantile conserving and probit transform filtering tools. +See the following link to learn more about these tools and how to use them: +https://docs.dart.ucar.edu/en/quantile_methods/models/lorenz_96_tracer_advection/work/readme.html + +Using this input table allows the user to specify the control options for the Quantile Conserving Filter (QCF), also known as the Quantile Conserving Ensemble Filtering Framework (QCEFF). The observation, state, and inflation variables are all included in this single table. + +The new quantile options are read in from the table at runtime and then set in the module algorithm_info_mod.f90 in the DART/assimilation_code/modules/assimilation directory. This module provides routines that give information about details of algorithms for observation error sampling, observation increments, and the transformations for regression and inflation in probit space. + +For individual QTYs in DART, the user can specify the options such as the bounds, distribution type, filter kind, etc. for the obs_error_info, probit_dist_info, and obs_inc_info subroutines in algorithm_info_mod.f90 + +If the user does not use a QCF input table with the DART quantile conserving and probit transform filtering tools, then the default values for these options will be used for all QTYs. + +Table Composition +----------------- +Each QTY is specified in its own column, having 28 total control options. +These control options are divided into 3 main groups, which are the options used for the obs_error_info, probit_dist_info, and obs_inc_info. However, the user is able to specify different values for probit inflation, probit state, and probit extended state, resulting in 5 total groupings for the control options. + +The obs_error_info subroutine computes information needed to compute error sample for this observation. +For obs_error_info the input options are the two bounds (lower and upper). + +The probit_dist_info subroutine computes the details of the probit transform. +From probit_dist_info, the values needed are the bounds and the distribution type. These can be different for all three cases (inflation, state, and extended_state). + +The obs_inc_info subrotuine sets the details of how to assimilate this observation. +From obs_inc_info, the bounds, plus the filter_kind, rectangular_quadrature, gaussian_likelihood_tails, sort_obs_inc, and spread_restoration are needed. However, rectangular_quadrature and gaussian_likelihood_tails are only applicable with RHF. + +Full list of options: +Obs_error_info: bounded_below, bounded_above, lower_bound, upper_bound [4 columns] +Probit_dist_info: dist_type, bounded_below, bounded_above, lower_bound, upper_bound (x3 for inflation, state, and observation (extended state) priors) [15 columns] +Obs_inc_info: filter_kind, rectangular_quadrature, gaussian_likelihood_tails, sort_obs_inc, spread_restoration, bounded_below, bounded_above, lower_bound, upper_bound [9 columns] + +Customizing the Table +--------------------- +The table can be customized by either editing a YAML file (which is then converted to a tabular data file in .txt format by a Python script) or a Google Sheet spreadsheet (which is then downloaded in .csv format). The specifics of how to manually edit both formats will be detailed in the following sections. + +Regardless of which of these formats are used, the table consists of two headers. The first states the version # of the table being used; the most recent version of the table needs to be used to ensure compatibilty with DART. The current version # is 1. The second header lists the full set of input options, or all 28 column names in other words. + +Generally, the user will add and fill in one row for each bounded QTY. If a QTY is not listed in the table, the default values will be used for all 28 options. Therefore, the user will only need to add rows for QTYs that use non-default values for any of the input options. ** + +The majority of the input options are read in as logicals, and will need to be written in the format of either 'F' or '.false.' These include bounded_below, bounded_above, rectangular_quadrature, gaussian_likelihood_tails, sort_obs_inc, and spread_restoration. + +The actual numerical values of the bounds are read in as real_r8 types. These can be specified as reals or integers in the table. + +dist_type and filter_kind are read in as strings, which the possible values for are listed below: + +dist_type: +NORMAL_DISTRIBUTION +BOUNDED_NORMAL_RH_DISTRIBUTION +GAMMA_DISTRIBUTION +BETA_DISTRIBUTION +LOG_NORMAL_DISTRIBUTION +UNIFORM_DISTRIBUTION +PARTICLE_FILTER_DISTRIBUTION + +filter_kind: +EAKF +ENKF +UNBOUNDED_RHF +GAMMA_FILTER +BOUNDED_NORMAL_RHF + +The default values for each of the options are listed below: +bounded_below = .false. +bounded_above = .false. +lower_bound = -888888 +upper_bound = -888888 +dist_type = BOUNDED_NORMAL_RH_DISTRIBUTION +filter_kind = BOUNDED_NORMAL_RHF +rectangular_quadrature = .false. +gaussian_likelihood_tails = .false. +sort_obs_inc = .false. +spread_restoration = .false. + +Note that bounds set to -888888 are missing_r8 values. + +YAML File Usage +--------------- +This section will detail how to customize the qcf_table_template.yaml file and then utilize the yaml_to_table.py Python script to convert the YAML dictionary into a table in .txt format. + +First, the user needs to access YAML template file, located in DART/assimilation/programs/qcf_table/ +This template file is then to be copied into another file. You can name this anything, but the standard name is 'qcf_table.yaml'. + +.. code:: + cp qcf_table_template.yaml qcf_table.yaml + +The YAML file needs to match the formatting in qcf_table_template.yaml, which is as follows: + +:: + + QCF table version: 1 + QTY_TEMPLATE: + obs_error_info: + bounded_below + bounded_above + lower_bound + upper_bound + probit_inflation: + dist_type + bounded_below + bounded_above + lower_bound + upper_bound + probit_state: + dist_type + bounded_below + bounded_above + lower_bound + upper_bound + probit_extended_state: + dist_type + bounded_below + bounded_above + lower_bound + upper_bound + obs_inc_info: + filter_kind + rectangular_quadrature + gaussian_likelihood_tails + sort_obs_inc + spread_restoration + bounded_below + bounded_above + lower_bound + upper_bound + QTY_STATE_VARIABLE: + obs_error_info: + bounded_below: .false. + bounded_above: .false. + lower_bound: -888888.0 + upper_bound: -888888.0 + probit_inflation: + dist_type: BOUNDED_NORMAL_RH_DISTRIBUTION + bounded_below: .false. + bounded_above: .false. + lower_bound: -888888.0 + upper_bound: -888888.0 + probit_state: + dist_type: BOUNDED_NORMAL_RH_DISTRIBUTION + bounded_below: .false. + bounded_above: .false. + lower_bound: -888888.0 + upper_bound: -888888.0 + probit_extended_state: + dist_type: BOUNDED_NORMAL_RH_DISTRIBUTION + bounded_below: .false. + bounded_above: .false. + lower_bound: -888888.0 + upper_bound: -888888.0 + obs_inc_info: + filter_kind: BOUNDED_NORMAL_RHF + rectangular_quadrature: .false. + gaussian_likelihood_tails: .false. + sort_obs_inc: .false. + spread_restoration: .false. + bounded_below: .false. + bounded_above: .false. + lower_bound: -888888.0 + upper_bound: -888888.0 + +To customize the YAML dictionary file, the user should change the name 'QTY_STATE_VARIABLE' to the name of the first QTY to be specified with non-default values. Edit the values for the vairables wanting to be changed, and leave the rest of the variables set to the default values. + +To add additional QTYs after this, simply copy the lines pertaining to first QTY, change the name of the QTY, and set the variables accordingly. + +To remove a QTY from the YAML dictionary, simply remove the lines it consists of. + +The user will then take their customized YAML file and pass it as input into a Python script. This will convert it into a text file contaning the table data. + +This script is located in DART/assimilation/programs/qcf_table/ + +To use the Python script on Derecho or Cheyenne, the user must first load the correct modules + +:: + + module load conda + conda activate npl + +Then run the python script. + +:: + + python3 yaml_to_table.py + +The user will be prompted to enter the name of the input YAML file and the name for the output text file name. +A table will be produced at the specified output filename. + +Copy or move this file to your working directory. + +Google Sheets Usage +------------------- +This section will detail how to customize the Google Sheets spreadsheet and then download the spreadsheet into a table in .csv format. + +Folow this link https://docs.google.com/spreadsheets/d/1SI4wHBXatLAAMfiMx3mUUC7x0fqz4lniKuM4_i5j6bM/edit#gid=0 to access the template spreadsheet. + +The QTYs listed in the template file (QTY_STATE_VARIABLE, QTY_TRACER_SOURCE) correspond to the lorenz_6_tracer_advection model and have the default values set for all variables. Make sure to remove these QTYs if you are not running an analagous model. ** + +Make a copy of the table by selecting 'File > Make a copy' from the menu bar. + +To customize the spreadsheet, click on the cell you want to edit and change the value of that cell. +To add a new QTY to the spreadsheet, simply copy the row of a listed QTY, change the QTY name, and edit the cells individually to set the control options. +To remove a QTY from the spreadsheet, select the row corresponding to that QTY. Then right click and choose "Delete Row" + +Ensure that there are no empty rows in between the QTYs listed in the spreadsheet. + +Download the spreadsheet as a .csv file by selecting 'File > Download > csv' from the menu bar. + +Google Sheets will append the name of the file with " - Sheet1.csv". For example a spreadsheet named "qcf_table" wil be downloaded as "qcf_table - Sheet1.csv" +Rename this file to remove this append to ensure that there are no spaces in the filename. + +Copy or move this file to your working directory. + +Using the table in DART +----------------------- +Navigate to your working directory. + +Edit your namelist file (input.nml) +Add the item "qcf_table_filename = 'your_filename' to the &filter_nml section, replacing your_filename with the actual name of the file you want to use. +Remember that the default values will be used for all QTYs if no filename is listed here. + +Build and run filter normally. + +The data read from the QCF table used is written to the output file dart_log.out From b37cf76791c3138a796202f657c17e9ce1ff5239 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Mon, 25 Sep 2023 19:02:05 -0600 Subject: [PATCH 45/60] Removing unused routines from utilities_mod use list --- assimilation_code/modules/assimilation/algorithm_info_mod.f90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index dbe9027985..2a01884641 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -9,7 +9,7 @@ module algorithm_info_mod use obs_def_mod, only : obs_def_type, get_obs_def_type_of_obs, get_obs_def_error_variance use obs_kind_mod, only : get_quantity_for_type_of_obs, get_name_for_quantity, get_index_for_quantity -use utilities_mod, only : error_handler, E_ALLMSG, E_ERR, E_MSG, log_it, logfileunit, open_file, close_file, to_upper +use utilities_mod, only : error_handler, E_ERR, E_MSG, open_file, close_file, to_upper use assim_model_mod, only : get_state_meta_data use location_mod, only : location_type From 668187559445480c4a779d3fd634424abf976fff Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Mon, 25 Sep 2023 19:18:48 -0600 Subject: [PATCH 46/60] Adding &probit_transform_nml to input.nml files for lorenz 63 and 96 - this will allow for the build and run lorenz github actions to pass --- models/lorenz_63/work/input.nml | 13 +++++++++++-- models/lorenz_96/work/input.nml | 13 +++++++++++-- 2 files changed, 22 insertions(+), 4 deletions(-) diff --git a/models/lorenz_63/work/input.nml b/models/lorenz_63/work/input.nml index 90504677f9..ab4a045d01 100644 --- a/models/lorenz_63/work/input.nml +++ b/models/lorenz_63/work/input.nml @@ -1,3 +1,9 @@ +&probit_transform_nml + fix_bound_violations = .false., + use_logit_instead_of_probit = .false. + do_inverse_check = .true. + / + &perfect_model_obs_nml read_input_state_from_file = .true., single_file_in = .true. @@ -28,6 +34,8 @@ / &filter_nml + qcf_table_filename = '', + use_algorithm_info_mod = .true., single_file_in = .true., input_state_files = '' input_state_file_list = 'filter_input_list.txt' @@ -86,6 +94,7 @@ / &assim_tools_nml + use_algorithm_info_mod = .true., filter_kind = 1, cutoff = 1000000.0 sort_obs_inc = .false., @@ -138,9 +147,9 @@ &preprocess_nml overwrite_output = .true. input_obs_def_mod_file = '../../../observations/forward_operators/DEFAULT_obs_def_mod.F90' - output_obs_def_mod_file = '../../../observations/forward_operators/obs_def_mod.f90' +output_obs_def_mod_file = './obs_def_mod.f90' input_obs_qty_mod_file = '../../../assimilation_code/modules/observations/DEFAULT_obs_kind_mod.F90' - output_obs_qty_mod_file = '../../../assimilation_code/modules/observations/obs_kind_mod.f90' +output_obs_qty_mod_file = './obs_kind_mod.f90' obs_type_files = '../../../observations/forward_operators/obs_def_1d_state_mod.f90' quantity_files = '../../../assimilation_code/modules/observations/oned_quantities_mod.f90' / diff --git a/models/lorenz_96/work/input.nml b/models/lorenz_96/work/input.nml index 39bbef5f02..03dceaf59c 100644 --- a/models/lorenz_96/work/input.nml +++ b/models/lorenz_96/work/input.nml @@ -1,3 +1,9 @@ +&probit_transform_nml + fix_bound_violations = .false., + use_logit_instead_of_probit = .false. + do_inverse_check = .true. + / + &perfect_model_obs_nml read_input_state_from_file = .true., single_file_in = .true. @@ -28,6 +34,8 @@ / &filter_nml + qcf_table_filename = '', + use_algorithm_info_mod = .true., single_file_in = .true., input_state_files = '' input_state_file_list = 'filter_input_list.txt' @@ -86,6 +94,7 @@ / &assim_tools_nml + use_algorithm_info_mod = .true., filter_kind = 1, cutoff = 0.02, sort_obs_inc = .false., @@ -144,9 +153,9 @@ &preprocess_nml overwrite_output = .true. input_obs_def_mod_file = '../../../observations/forward_operators/DEFAULT_obs_def_mod.F90' - output_obs_def_mod_file = '../../../observations/forward_operators/obs_def_mod.f90' +output_obs_def_mod_file = './obs_def_mod.f90' input_obs_qty_mod_file = '../../../assimilation_code/modules/observations/DEFAULT_obs_kind_mod.F90' - output_obs_qty_mod_file = '../../../assimilation_code/modules/observations/obs_kind_mod.f90' +output_obs_qty_mod_file = './obs_kind_mod.f90' obs_type_files = '../../../observations/forward_operators/obs_def_1d_state_mod.f90' quantity_files = '../../../assimilation_code/modules/observations/oned_quantities_mod.f90' / From 3e39cf226564458c34823a3981b7ca7cfc892726 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Tue, 26 Sep 2023 12:43:53 -0600 Subject: [PATCH 47/60] Removing the use_algorithm_info_mod logical from assim_tools_mod and filter_mod --- .../modules/assimilation/assim_tools_mod.f90 | 12 +++++------- .../modules/assimilation/filter_mod.f90 | 15 +++------------ 2 files changed, 8 insertions(+), 19 deletions(-) diff --git a/assimilation_code/modules/assimilation/assim_tools_mod.f90 b/assimilation_code/modules/assimilation/assim_tools_mod.f90 index 958e373201..b604807bac 100644 --- a/assimilation_code/modules/assimilation/assim_tools_mod.f90 +++ b/assimilation_code/modules/assimilation/assim_tools_mod.f90 @@ -143,7 +143,6 @@ module assim_tools_mod ! special_localization_obs_types -> Special treatment for the specified observation types ! special_localization_cutoffs -> Different cutoff value for each specified obs type ! -logical :: use_algorithm_info_mod = .true. integer :: filter_kind = 1 real(r8) :: cutoff = 0.2_r8 logical :: sort_obs_inc = .false. @@ -203,8 +202,7 @@ module assim_tools_mod ! compared to previous versions of this namelist item. logical :: distribute_mean = .false. -namelist / assim_tools_nml / use_algorithm_info_mod, & - filter_kind, cutoff, sort_obs_inc, & +namelist / assim_tools_nml / filter_kind, cutoff, sort_obs_inc, & spread_restoration, sampling_error_correction, & adaptive_localization_threshold, adaptive_cutoff_floor, & print_every_nth_obs, rectangular_quadrature, gaussian_likelihood_tails, & @@ -976,7 +974,8 @@ subroutine obs_increment(ens_in, ens_size, obs, obs_var, obs_kind, obs_inc, & !--------------------------begin algorithm_info control block----------------- ! More flexible abilities to control the observation space increments are ! available with this code block. It gets information about the increment method -! for the current observation is use_algorithm_info_mod is set to true in the namelist. +! for the current observation. + ! This is not an extensible mechanism for doing this as the number of ! obs increments distributions and associated information goes up ! Implications for sorting increments and for spread restoration need to be examined @@ -988,9 +987,8 @@ subroutine obs_increment(ens_in, ens_size, obs, obs_var, obs_kind, obs_inc, & bounded_below = .false.; lower_bound = 0.0_r8 bounded_above = .false.; upper_bound = 0.0_r8 -if(use_algorithm_info_mod) & - call obs_inc_info(obs_kind, filter_kind, rectangular_quadrature, gaussian_likelihood_tails, & - sort_obs_inc, spread_restoration, bounded_below, bounded_above, lower_bound, upper_bound) +call obs_inc_info(obs_kind, filter_kind, rectangular_quadrature, gaussian_likelihood_tails, & + sort_obs_inc, spread_restoration, bounded_below, bounded_above, lower_bound, upper_bound) ! Could add logic to check on sort being true when not needed. ! Could also add logic to limit the use of spread_restoration to EAKF. It will fail diff --git a/assimilation_code/modules/assimilation/filter_mod.f90 b/assimilation_code/modules/assimilation/filter_mod.f90 index c83a433a8a..1a4dbfe585 100644 --- a/assimilation_code/modules/assimilation/filter_mod.f90 +++ b/assimilation_code/modules/assimilation/filter_mod.f90 @@ -167,7 +167,6 @@ module filter_mod ! Namelist input with default values ! character(len = 129) :: qcf_table_filename = '' -logical :: use_algorithm_info_mod = .true. integer :: async = 0, ens_size = 20 integer :: tasks_per_model_advance = 1 ! if init_time_days and seconds are negative initial time is 0, 0 @@ -263,7 +262,6 @@ module filter_mod namelist /filter_nml/ async, & qcf_table_filename, & - use_algorithm_info_mod, & adv_ens_command, & ens_size, & tasks_per_model_advance, & @@ -1608,16 +1606,9 @@ subroutine filter_ensemble_inflate(ens_handle, inflate_copy, inflate, ENS_MEAN_C call get_state_meta_data(ens_handle%my_vars(j), my_state_loc, my_state_kind) ! Need to specify what kind of prior to use for each - ! Use default of untransformed if use_algorithm_info_mod is not true - if(use_algorithm_info_mod) then - call probit_dist_info(my_state_kind, .true., .true., dist_type, & - bounded_below, bounded_above, lower_bound, upper_bound) - else - ! Default is just a normal which does nothing - dist_type = NORMAL_DISTRIBUTION - bounded_below = .false. ; bounded_above = .false. - lower_bound = 0.0_r8; upper_bound = 0.0_r8 - endif + call probit_dist_info(my_state_kind, .true., .true., dist_type, & + bounded_below, bounded_above, lower_bound, upper_bound) + call transform_to_probit(grp_size, ens_handle%copies(grp_bot:grp_top, j), & dist_type, dist_params, probit_ens(1:grp_size), .false., & bounded_below, bounded_above, lower_bound, upper_bound) From 80305c74db8e24a50a902b7e83a1cc3ec1b40bf8 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Tue, 26 Sep 2023 12:53:09 -0600 Subject: [PATCH 48/60] Fixing unintentional additions to lorenx 96 and 63 input.nml files; removing use_algorithm_info_mod nml option --- models/lorenz_63/work/input.nml | 6 ++---- models/lorenz_96/work/input.nml | 6 ++---- 2 files changed, 4 insertions(+), 8 deletions(-) diff --git a/models/lorenz_63/work/input.nml b/models/lorenz_63/work/input.nml index ab4a045d01..3ed5062322 100644 --- a/models/lorenz_63/work/input.nml +++ b/models/lorenz_63/work/input.nml @@ -35,7 +35,6 @@ &filter_nml qcf_table_filename = '', - use_algorithm_info_mod = .true., single_file_in = .true., input_state_files = '' input_state_file_list = 'filter_input_list.txt' @@ -94,7 +93,6 @@ / &assim_tools_nml - use_algorithm_info_mod = .true., filter_kind = 1, cutoff = 1000000.0 sort_obs_inc = .false., @@ -147,9 +145,9 @@ &preprocess_nml overwrite_output = .true. input_obs_def_mod_file = '../../../observations/forward_operators/DEFAULT_obs_def_mod.F90' -output_obs_def_mod_file = './obs_def_mod.f90' +output_obs_def_mod_file = '../../../observations/forward_operators/obs_def_mod.f90' input_obs_qty_mod_file = '../../../assimilation_code/modules/observations/DEFAULT_obs_kind_mod.F90' -output_obs_qty_mod_file = './obs_kind_mod.f90' +output_obs_qty_mod_file = '../../../assimilation_code/modules/observations/obs_kind_mod.f90' obs_type_files = '../../../observations/forward_operators/obs_def_1d_state_mod.f90' quantity_files = '../../../assimilation_code/modules/observations/oned_quantities_mod.f90' / diff --git a/models/lorenz_96/work/input.nml b/models/lorenz_96/work/input.nml index 03dceaf59c..8aa1ce99ad 100644 --- a/models/lorenz_96/work/input.nml +++ b/models/lorenz_96/work/input.nml @@ -35,7 +35,6 @@ &filter_nml qcf_table_filename = '', - use_algorithm_info_mod = .true., single_file_in = .true., input_state_files = '' input_state_file_list = 'filter_input_list.txt' @@ -94,7 +93,6 @@ / &assim_tools_nml - use_algorithm_info_mod = .true., filter_kind = 1, cutoff = 0.02, sort_obs_inc = .false., @@ -153,9 +151,9 @@ &preprocess_nml overwrite_output = .true. input_obs_def_mod_file = '../../../observations/forward_operators/DEFAULT_obs_def_mod.F90' -output_obs_def_mod_file = './obs_def_mod.f90' +output_obs_def_mod_file = '../../../observations/forward_operators/obs_def_mod.f90' input_obs_qty_mod_file = '../../../assimilation_code/modules/observations/DEFAULT_obs_kind_mod.F90' -output_obs_qty_mod_file = './obs_kind_mod.f90' +output_obs_qty_mod_file = '../../../assimilation_code/modules/observations/obs_kind_mod.f90' obs_type_files = '../../../observations/forward_operators/obs_def_1d_state_mod.f90' quantity_files = '../../../assimilation_code/modules/observations/oned_quantities_mod.f90' / From eb38d48d99a96aebb4c51404346232167b34b813 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Tue, 26 Sep 2023 13:14:02 -0600 Subject: [PATCH 49/60] Adding spaces back in from unintentional removal --- models/lorenz_63/work/input.nml | 4 ++-- models/lorenz_96/work/input.nml | 4 ++-- 2 files changed, 4 insertions(+), 4 deletions(-) diff --git a/models/lorenz_63/work/input.nml b/models/lorenz_63/work/input.nml index 3ed5062322..86a7e2e569 100644 --- a/models/lorenz_63/work/input.nml +++ b/models/lorenz_63/work/input.nml @@ -145,9 +145,9 @@ &preprocess_nml overwrite_output = .true. input_obs_def_mod_file = '../../../observations/forward_operators/DEFAULT_obs_def_mod.F90' -output_obs_def_mod_file = '../../../observations/forward_operators/obs_def_mod.f90' + output_obs_def_mod_file = '../../../observations/forward_operators/obs_def_mod.f90' input_obs_qty_mod_file = '../../../assimilation_code/modules/observations/DEFAULT_obs_kind_mod.F90' -output_obs_qty_mod_file = '../../../assimilation_code/modules/observations/obs_kind_mod.f90' + output_obs_qty_mod_file = '../../../assimilation_code/modules/observations/obs_kind_mod.f90' obs_type_files = '../../../observations/forward_operators/obs_def_1d_state_mod.f90' quantity_files = '../../../assimilation_code/modules/observations/oned_quantities_mod.f90' / diff --git a/models/lorenz_96/work/input.nml b/models/lorenz_96/work/input.nml index 8aa1ce99ad..615ef8df5d 100644 --- a/models/lorenz_96/work/input.nml +++ b/models/lorenz_96/work/input.nml @@ -151,9 +151,9 @@ &preprocess_nml overwrite_output = .true. input_obs_def_mod_file = '../../../observations/forward_operators/DEFAULT_obs_def_mod.F90' -output_obs_def_mod_file = '../../../observations/forward_operators/obs_def_mod.f90' + output_obs_def_mod_file = '../../../observations/forward_operators/obs_def_mod.f90' input_obs_qty_mod_file = '../../../assimilation_code/modules/observations/DEFAULT_obs_kind_mod.F90' -output_obs_qty_mod_file = '../../../assimilation_code/modules/observations/obs_kind_mod.f90' + output_obs_qty_mod_file = '../../../assimilation_code/modules/observations/obs_kind_mod.f90' obs_type_files = '../../../observations/forward_operators/obs_def_1d_state_mod.f90' quantity_files = '../../../assimilation_code/modules/observations/oned_quantities_mod.f90' / From 528e1629d3a325067f7eef146784177f4b825a7f Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Tue, 26 Sep 2023 15:54:31 -0600 Subject: [PATCH 50/60] Adding &probit_transform_nml and qcf_table_filename option to model's input.nml files (9var, am2, bgrid_solo) --- models/9var/work/input.nml | 7 +++++++ models/am2/work/input.nml | 7 +++++++ models/bgrid_solo/work/input.nml | 7 +++++++ 3 files changed, 21 insertions(+) diff --git a/models/9var/work/input.nml b/models/9var/work/input.nml index cfe617f0ce..168834b85a 100644 --- a/models/9var/work/input.nml +++ b/models/9var/work/input.nml @@ -1,3 +1,9 @@ +&probit_transform_nml + fix_bound_violations = .false., + use_logit_instead_of_probit = .false. + do_inverse_check = .true. + / + &perfect_model_obs_nml read_input_state_from_file = .true., single_file_in = .true. @@ -32,6 +38,7 @@ # output_state_files = 'filter_output.nc' &filter_nml + qcf_table_filename = '' single_file_in = .true., input_state_files = '' input_state_file_list = 'filter_input_list.txt' diff --git a/models/am2/work/input.nml b/models/am2/work/input.nml index 5f453d01cd..171974666a 100644 --- a/models/am2/work/input.nml +++ b/models/am2/work/input.nml @@ -1,3 +1,9 @@ +&probit_transform_nml + fix_bound_violations = .false., + use_logit_instead_of_probit = .false. + do_inverse_check = .true. + / + &perfect_model_obs_nml start_from_restart = .true., output_restart = .true., @@ -22,6 +28,7 @@ / &filter_nml + qcf_table_filename = '' async = 2, adv_ens_command = "./advance_model.csh", ens_size = 10, diff --git a/models/bgrid_solo/work/input.nml b/models/bgrid_solo/work/input.nml index 9d65871b00..54e29e3f0b 100644 --- a/models/bgrid_solo/work/input.nml +++ b/models/bgrid_solo/work/input.nml @@ -1,3 +1,9 @@ +&probit_transform_nml + fix_bound_violations = .false., + use_logit_instead_of_probit = .false. + do_inverse_check = .true. + / + &perfect_model_obs_nml read_input_state_from_file = .true., single_file_in = .true. @@ -28,6 +34,7 @@ / &filter_nml + qcf_table_filename = '' single_file_in = .true., input_state_files = 'filter_input.nc' input_state_file_list = '' From caf5dac9ecae6e821dbe9fa391ea035eac81d4b8 Mon Sep 17 00:00:00 2001 From: Marlee Smith Date: Wed, 27 Sep 2023 15:51:34 -0600 Subject: [PATCH 51/60] Fixing missing commas in nml files --- models/9var/work/input.nml | 6 +++--- models/am2/work/input.nml | 6 +++--- models/bgrid_solo/work/input.nml | 6 +++--- 3 files changed, 9 insertions(+), 9 deletions(-) diff --git a/models/9var/work/input.nml b/models/9var/work/input.nml index 168834b85a..419aec1d3a 100644 --- a/models/9var/work/input.nml +++ b/models/9var/work/input.nml @@ -1,7 +1,7 @@ &probit_transform_nml fix_bound_violations = .false., - use_logit_instead_of_probit = .false. - do_inverse_check = .true. + use_logit_instead_of_probit = .false., + do_inverse_check = .true., / &perfect_model_obs_nml @@ -38,7 +38,7 @@ # output_state_files = 'filter_output.nc' &filter_nml - qcf_table_filename = '' + qcf_table_filename = '', single_file_in = .true., input_state_files = '' input_state_file_list = 'filter_input_list.txt' diff --git a/models/am2/work/input.nml b/models/am2/work/input.nml index 171974666a..83825dab00 100644 --- a/models/am2/work/input.nml +++ b/models/am2/work/input.nml @@ -1,7 +1,7 @@ &probit_transform_nml fix_bound_violations = .false., - use_logit_instead_of_probit = .false. - do_inverse_check = .true. + use_logit_instead_of_probit = .false., + do_inverse_check = .true., / &perfect_model_obs_nml @@ -28,7 +28,7 @@ / &filter_nml - qcf_table_filename = '' + qcf_table_filename = '', async = 2, adv_ens_command = "./advance_model.csh", ens_size = 10, diff --git a/models/bgrid_solo/work/input.nml b/models/bgrid_solo/work/input.nml index 54e29e3f0b..7311b0464e 100644 --- a/models/bgrid_solo/work/input.nml +++ b/models/bgrid_solo/work/input.nml @@ -1,7 +1,7 @@ &probit_transform_nml fix_bound_violations = .false., - use_logit_instead_of_probit = .false. - do_inverse_check = .true. + use_logit_instead_of_probit = .false., + do_inverse_check = .true., / &perfect_model_obs_nml @@ -34,7 +34,7 @@ / &filter_nml - qcf_table_filename = '' + qcf_table_filename = '', single_file_in = .true., input_state_files = 'filter_input.nc' input_state_file_list = '' From c184e907b969cf2057e755eab6825bc49dbc5d43 Mon Sep 17 00:00:00 2001 From: Helen Kershaw Date: Thu, 28 Sep 2023 15:14:43 -0400 Subject: [PATCH 52/60] add test cases for reading qcf table --- .gitignore | 1 + developer_tests/qceff/test_table_read.f90 | 23 +++++++++++ developer_tests/qceff/work/input.nml | 28 +++++++++++++ developer_tests/qceff/work/qcf_table.txt | 3 ++ .../qceff/work/qcf_table_bad_qty.txt | 3 ++ .../qceff/work/qcf_table_broke.txt | 3 ++ .../qceff/work/qcf_table_extra_columns.txt | 3 ++ developer_tests/qceff/work/qcf_table_v2.txt | 3 ++ developer_tests/qceff/work/quickbuild.sh | 40 +++++++++++++++++++ 9 files changed, 107 insertions(+) create mode 100644 developer_tests/qceff/test_table_read.f90 create mode 100644 developer_tests/qceff/work/input.nml create mode 100644 developer_tests/qceff/work/qcf_table.txt create mode 100644 developer_tests/qceff/work/qcf_table_bad_qty.txt create mode 100644 developer_tests/qceff/work/qcf_table_broke.txt create mode 100644 developer_tests/qceff/work/qcf_table_extra_columns.txt create mode 100644 developer_tests/qceff/work/qcf_table_v2.txt create mode 100755 developer_tests/qceff/work/quickbuild.sh diff --git a/.gitignore b/.gitignore index cfbf4153d3..bf6ca1dcbd 100644 --- a/.gitignore +++ b/.gitignore @@ -189,6 +189,7 @@ stacktest obs_rwtest test_quad_irreg_interp test_quad_reg_interp +test_table_read # Directories to NOT IGNORE ... same as executable names # as far as I know, these must be listed after the executables diff --git a/developer_tests/qceff/test_table_read.f90 b/developer_tests/qceff/test_table_read.f90 new file mode 100644 index 0000000000..c9e97f4d1c --- /dev/null +++ b/developer_tests/qceff/test_table_read.f90 @@ -0,0 +1,23 @@ +program test_table_read + +use algorithm_info_mod, only : init_algorithm_info_mod, end_algorithm_info_mod +use utilities_mod, only : initialize_utilities + +implicit none + +character(len=129) :: qcf_table_filename + +call initialize_utilities('test_table_read') + +!n = command_argument_count() +call get_command_argument(1,qcf_table_filename) + + +!qcf_table_filename = 'qcf_table_v2.txt' + +call init_algorithm_info_mod(qcf_table_filename) + +call end_algorithm_info_mod() + + +end program test_table_read diff --git a/developer_tests/qceff/work/input.nml b/developer_tests/qceff/work/input.nml new file mode 100644 index 0000000000..ffa7155436 --- /dev/null +++ b/developer_tests/qceff/work/input.nml @@ -0,0 +1,28 @@ +&utilities_nml + TERMLEVEL = 1, + module_details = .false. + logfilename = 'dart_log.out' + / + +# pick a random set of inputs +&preprocess_nml + overwrite_output = .true. + input_obs_qty_mod_file = '../../../assimilation_code/modules/observations/DEFAULT_obs_kind_mod.F90' + output_obs_qty_mod_file = '../../../assimilation_code/modules/observations/obs_kind_mod.f90' + input_obs_def_mod_file = '../../../observations/forward_operators/DEFAULT_obs_def_mod.F90' + output_obs_def_mod_file = '../../../observations/forward_operators/obs_def_mod.f90' + obs_type_files = '../../../observations/forward_operators/obs_def_reanalysis_bufr_mod.f90', + '../../../observations/forward_operators/obs_def_radar_mod.f90', + '../../../observations/forward_operators/obs_def_metar_mod.f90', + '../../../observations/forward_operators/obs_def_dew_point_mod.f90', + '../../../observations/forward_operators/obs_def_rel_humidity_mod.f90', + '../../../observations/forward_operators/obs_def_altimeter_mod.f90', + '../../../observations/forward_operators/obs_def_gps_mod.f90', + '../../../observations/forward_operators/obs_def_vortex_mod.f90', + '../../../observations/forward_operators/obs_def_gts_mod.f90', + '../../../observations/forward_operators/obs_def_QuikSCAT_mod.f90' + quantity_files = '../../../assimilation_code/modules/observations/default_quantities_mod.f90', + / + +&obs_kind_nml +/ diff --git a/developer_tests/qceff/work/qcf_table.txt b/developer_tests/qceff/work/qcf_table.txt new file mode 100644 index 0000000000..7d4f146540 --- /dev/null +++ b/developer_tests/qceff/work/qcf_table.txt @@ -0,0 +1,3 @@ +QCF table version: 1 +QTY_TEMPLATE: bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound filter_kind rectangular_quadrature gaussian_likelihood_tails sort_obs_inc spread_restoration bounded_below bounded_above lower_bound upper_bound +QTY_STATE_VARIABLE .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RHF .false. .false. .false. .false. .false. .false. -888888.0 -888888.0 diff --git a/developer_tests/qceff/work/qcf_table_bad_qty.txt b/developer_tests/qceff/work/qcf_table_bad_qty.txt new file mode 100644 index 0000000000..428e5fd6c5 --- /dev/null +++ b/developer_tests/qceff/work/qcf_table_bad_qty.txt @@ -0,0 +1,3 @@ +QCF table version: 1 +QTY_TEMPLATE: bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound filter_kind rectangular_quadrature gaussian_likelihood_tails sort_obs_inc spread_restoration bounded_below bounded_above lower_bound upper_bound +QTY_HAIRCUT .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RHF .false. .false. .false. .false. .false. .false. -888888.0 -888888.0 diff --git a/developer_tests/qceff/work/qcf_table_broke.txt b/developer_tests/qceff/work/qcf_table_broke.txt new file mode 100644 index 0000000000..cb78e95e49 --- /dev/null +++ b/developer_tests/qceff/work/qcf_table_broke.txt @@ -0,0 +1,3 @@ +QCF table version: 1 +QTY_TEMPLATE: bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound filter_kind rectangular_quadrature gaussian_likelihood_tails sort_obs_inc spread_restoration bounded_below bounded_above lower_bound upper_bound +QTY_STATE_VARIABLE .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RHF .false. .false. .false. .false. .false. .false. -888888.0 -888888.0 diff --git a/developer_tests/qceff/work/qcf_table_extra_columns.txt b/developer_tests/qceff/work/qcf_table_extra_columns.txt new file mode 100644 index 0000000000..d298573349 --- /dev/null +++ b/developer_tests/qceff/work/qcf_table_extra_columns.txt @@ -0,0 +1,3 @@ +QCF table version: 1 +QTY_TEMPLATE: bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound filter_kind rectangular_quadrature gaussian_likelihood_tails sort_obs_inc spread_restoration bounded_below bounded_above lower_bound upper_bound frog +QTY_STATE_VARIABLE .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RHF .false. .false. .false. .false. .false. .false. -888888.0 -888888.0 toad newt diff --git a/developer_tests/qceff/work/qcf_table_v2.txt b/developer_tests/qceff/work/qcf_table_v2.txt new file mode 100644 index 0000000000..0d1b5a46df --- /dev/null +++ b/developer_tests/qceff/work/qcf_table_v2.txt @@ -0,0 +1,3 @@ +QCF table version: 2 +QTY_TEMPLATE: bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound filter_kind rectangular_quadrature gaussian_likelihood_tails sort_obs_inc spread_restoration bounded_below bounded_above lower_bound upper_bound +QTY_STATE_VARIABLE .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RHF .false. .false. .false. .false. .false. .false. -888888.0 -888888.0 diff --git a/developer_tests/qceff/work/quickbuild.sh b/developer_tests/qceff/work/quickbuild.sh new file mode 100755 index 0000000000..81b1308494 --- /dev/null +++ b/developer_tests/qceff/work/quickbuild.sh @@ -0,0 +1,40 @@ +#!/usr/bin/env bash + +# DART software - Copyright UCAR. This open source software is provided +# by UCAR, "as is", without charge, subject to all terms of use at +# http://www.image.ucar.edu/DAReS/DART/DART_download + +main() { + + +export DART=$(git rev-parse --show-toplevel) +source "$DART"/build_templates/buildfunctions.sh + +MODEL="none" +EXTRA="$DART"/models/template/threed_model_mod.f90 +dev_test=1 +TEST="qceff" +LOCATION="threed_sphere" + +serial_programs=( +test_table_read +) + +# quickbuild arguments +arguments "$@" + +# clean the directory +\rm -f -- *.o *.mod Makefile .cppdefs + +# build and run preprocess before making any other DART executables +buildpreprocess + +# build +buildit + +# clean up +\rm -f -- *.o *.mod + +} + +main "$@" From 2e0700873d97f7d01ff32bab585274030e1087da Mon Sep 17 00:00:00 2001 From: Helen Kershaw Date: Thu, 28 Sep 2023 16:35:50 -0400 Subject: [PATCH 53/60] test script to check return codes for various qcf_table cases --- developer_tests/qceff/test_table_read.f90 | 9 ++--- developer_tests/qceff/work/runall.sh | 40 +++++++++++++++++++++++ 2 files changed, 45 insertions(+), 4 deletions(-) create mode 100755 developer_tests/qceff/work/runall.sh diff --git a/developer_tests/qceff/test_table_read.f90 b/developer_tests/qceff/test_table_read.f90 index c9e97f4d1c..aaaa279c2a 100644 --- a/developer_tests/qceff/test_table_read.f90 +++ b/developer_tests/qceff/test_table_read.f90 @@ -1,3 +1,8 @@ +! DART software - Copyright UCAR. This open source software is provided +! by UCAR, "as is", without charge, subject to all terms of use at +! http://www.image.ucar.edu/DAReS/DART/DART_download + +! qcf_table_filename expected as command line arguement program test_table_read use algorithm_info_mod, only : init_algorithm_info_mod, end_algorithm_info_mod @@ -9,12 +14,8 @@ program test_table_read call initialize_utilities('test_table_read') -!n = command_argument_count() call get_command_argument(1,qcf_table_filename) - -!qcf_table_filename = 'qcf_table_v2.txt' - call init_algorithm_info_mod(qcf_table_filename) call end_algorithm_info_mod() diff --git a/developer_tests/qceff/work/runall.sh b/developer_tests/qceff/work/runall.sh new file mode 100755 index 0000000000..b9a4cf24ed --- /dev/null +++ b/developer_tests/qceff/work/runall.sh @@ -0,0 +1,40 @@ +#!/bin/bash + +# DART software - Copyright UCAR. This open source software is provided +# by UCAR, "as is", without charge, subject to all terms of use at +# http://www.image.ucar.edu/DAReS/DART/DART_download + +# Usage: +# ./runall.sh +# ./runall.sh | grep FAIL +# ./runall.sh | grep PASS + + +should_pass () { +if [[ $? -ne 0 ]]; then + echo $1: "FAIL" +else + echo $1: "PASS" +fi +} + +should_fail () { +if [[ $? -eq 0 ]]; then + echo $1: "FAIL" +else + echo $1: "PASS" +fi +} + +./test_table_read ; should_pass "no table" + +./test_table_read qcf_table.txt /dev/null ; should_pass "correct v1 table" + +./test_table_read qcf_table_v2.txt /dev/null ; should_fail "detect wrong version" + +./test_table_read qcf_table_extra_columns.txt /dev/null ; should_pass "extra colums" + +./test_table_read qcf_table_bad_qty.txt /dev/null ; should_fail "bad qty" + +./test_table_read qcf_table_broke.txt /dev/null ; should_fail "bad value" + From 200a8f4e3fe6dbc7b82cfacb551c8b8f2a6ba373 Mon Sep 17 00:00:00 2001 From: Helen Kershaw Date: Thu, 28 Sep 2023 16:38:55 -0400 Subject: [PATCH 54/60] remove stray /dev/null left in accidentally --- developer_tests/qceff/work/runall.sh | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/developer_tests/qceff/work/runall.sh b/developer_tests/qceff/work/runall.sh index b9a4cf24ed..d47b7a2545 100755 --- a/developer_tests/qceff/work/runall.sh +++ b/developer_tests/qceff/work/runall.sh @@ -28,13 +28,13 @@ fi ./test_table_read ; should_pass "no table" -./test_table_read qcf_table.txt /dev/null ; should_pass "correct v1 table" +./test_table_read qcf_table.txt ; should_pass "correct v1 table" -./test_table_read qcf_table_v2.txt /dev/null ; should_fail "detect wrong version" +./test_table_read qcf_table_v2.txt ; should_fail "detect wrong version" -./test_table_read qcf_table_extra_columns.txt /dev/null ; should_pass "extra colums" +./test_table_read qcf_table_extra_columns.txt ; should_pass "extra colums" -./test_table_read qcf_table_bad_qty.txt /dev/null ; should_fail "bad qty" +./test_table_read qcf_table_bad_qty.txt ; should_fail "bad qty" -./test_table_read qcf_table_broke.txt /dev/null ; should_fail "bad value" +./test_table_read qcf_table_broke.txt ; should_fail "bad value" From 46f238fe36b94a9daa7db3585f82e0d23eb0c0f2 Mon Sep 17 00:00:00 2001 From: Helen Kershaw Date: Thu, 28 Sep 2023 17:09:05 -0400 Subject: [PATCH 55/60] add tests for various bounds options currently the "lower bound only" test is failing because the upper < lower check happens always rather then only when you have two bounds --- developer_tests/qceff/work/qcf_table_duplicates.txt | 6 ++++++ developer_tests/qceff/work/qcf_table_lower_bound_only.txt | 6 ++++++ developer_tests/qceff/work/qcf_table_lower_gt_upper.txt | 6 ++++++ developer_tests/qceff/work/qcf_table_no_header.txt | 4 ++++ developer_tests/qceff/work/runall.sh | 5 +++++ 5 files changed, 27 insertions(+) create mode 100644 developer_tests/qceff/work/qcf_table_duplicates.txt create mode 100644 developer_tests/qceff/work/qcf_table_lower_bound_only.txt create mode 100644 developer_tests/qceff/work/qcf_table_lower_gt_upper.txt create mode 100644 developer_tests/qceff/work/qcf_table_no_header.txt diff --git a/developer_tests/qceff/work/qcf_table_duplicates.txt b/developer_tests/qceff/work/qcf_table_duplicates.txt new file mode 100644 index 0000000000..7ffddff61f --- /dev/null +++ b/developer_tests/qceff/work/qcf_table_duplicates.txt @@ -0,0 +1,6 @@ +QCF table version: 1 +QTY_TEMPLATE: bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound filter_kind rectangular_quadrature gaussian_likelihood_tails sort_obs_inc spread_restoration bounded_below bounded_above lower_bound upper_bound +QTY_VEGETATED_AREA_FRACTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RHF .false. .false. .false. .false. .false. .false. -888888.0 -888888.0 +QTY_AQUIFER_WATER .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RHF .false. .false. .false. .false. .false. .false. -888888.0 -888888.0 +QTY_SEAICE_SALINITY008 .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RHF .false. .false. .false. .false. .false. .false. -888888.0 -888888.0 +QTY_VEGETATED_AREA_FRACTION .true. .false. 0 1000 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RHF .false. .false. .false. .false. .false. .false. -888888.0 -888888.0 diff --git a/developer_tests/qceff/work/qcf_table_lower_bound_only.txt b/developer_tests/qceff/work/qcf_table_lower_bound_only.txt new file mode 100644 index 0000000000..6f0fca4ee4 --- /dev/null +++ b/developer_tests/qceff/work/qcf_table_lower_bound_only.txt @@ -0,0 +1,6 @@ +QCF table version: 1 +QTY_TEMPLATE: bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound filter_kind rectangular_quadrature gaussian_likelihood_tails sort_obs_inc spread_restoration bounded_below bounded_above lower_bound upper_bound +QTY_STATE_VARIABLE .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RHF .false. .false. .false. .false. .false. .false. -888888.0 -888888.0 +QTY_AQUIFER_WATER .true. .false. 0.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RHF .false. .false. .false. .false. .false. .false. -888888.0 -888888.0 +QTY_SEAICE_SALINITY008 .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RHF .false. .false. .false. .false. .false. .false. -888888.0 -888888.0 +QTY_VEGETATED_AREA_FRACTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RHF .false. .false. .false. .false. .false. .false. -888888.0 -888888.0 diff --git a/developer_tests/qceff/work/qcf_table_lower_gt_upper.txt b/developer_tests/qceff/work/qcf_table_lower_gt_upper.txt new file mode 100644 index 0000000000..6370c2cdd7 --- /dev/null +++ b/developer_tests/qceff/work/qcf_table_lower_gt_upper.txt @@ -0,0 +1,6 @@ +QCF table version: 1 +QTY_TEMPLATE: bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound filter_kind rectangular_quadrature gaussian_likelihood_tails sort_obs_inc spread_restoration bounded_below bounded_above lower_bound upper_bound +QTY_STATE_VARIABLE .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RHF .false. .false. .false. .false. .false. .false. -888888.0 -888888.0 +QTY_AQUIFER_WATER .true. .true. 10.0 0.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RHF .false. .false. .false. .false. .false. .false. -888888.0 -888888.0 +QTY_SEAICE_SALINITY008 .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RHF .false. .false. .false. .false. .false. .false. -888888.0 -888888.0 +QTY_VEGETATED_AREA_FRACTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RHF .false. .false. .false. .false. .false. .false. -888888.0 -888888.0 diff --git a/developer_tests/qceff/work/qcf_table_no_header.txt b/developer_tests/qceff/work/qcf_table_no_header.txt new file mode 100644 index 0000000000..617379ffe0 --- /dev/null +++ b/developer_tests/qceff/work/qcf_table_no_header.txt @@ -0,0 +1,4 @@ +QTY_STATE_VARIABLE .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RHF .false. .false. .false. .false. .false. .false. -888888.0 -888888.0 +QTY_AQUIFER_WATER .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RHF .false. .false. .false. .false. .false. .false. -888888.0 -888888.0 +QTY_SEAICE_SALINITY008 .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RHF .false. .false. .false. .false. .false. .false. -888888.0 -888888.0 +QTY_VEGETATED_AREA_FRACTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RHF .false. .false. .false. .false. .false. .false. -888888.0 -888888.0 diff --git a/developer_tests/qceff/work/runall.sh b/developer_tests/qceff/work/runall.sh index d47b7a2545..80b3fc87b1 100755 --- a/developer_tests/qceff/work/runall.sh +++ b/developer_tests/qceff/work/runall.sh @@ -38,3 +38,8 @@ fi ./test_table_read qcf_table_broke.txt ; should_fail "bad value" +./test_table_read qcf_table_no_header.txt ; should_fail "no header" + +./test_table_read qcf_table_lower_gt_upper.txt ; should_fail "upper bound less than lower" + +./test_table_read ./test_table_read qcf_table_lower_bound_only.txt ; should_pass "lower bound only" From ca99e8bcbcb7c08226fc71be436a9ed8bf7e71d2 Mon Sep 17 00:00:00 2001 From: Helen Kershaw Date: Thu, 28 Sep 2023 17:14:26 -0400 Subject: [PATCH 56/60] test for bounds set to false, but bounds values in the table --- .../qceff/work/qcf_table_no_bounds_with_values.txt | 6 ++++++ developer_tests/qceff/work/runall.sh | 2 ++ 2 files changed, 8 insertions(+) create mode 100644 developer_tests/qceff/work/qcf_table_no_bounds_with_values.txt diff --git a/developer_tests/qceff/work/qcf_table_no_bounds_with_values.txt b/developer_tests/qceff/work/qcf_table_no_bounds_with_values.txt new file mode 100644 index 0000000000..32c8d4f8e9 --- /dev/null +++ b/developer_tests/qceff/work/qcf_table_no_bounds_with_values.txt @@ -0,0 +1,6 @@ +QCF table version: 1 +QTY_TEMPLATE: bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound filter_kind rectangular_quadrature gaussian_likelihood_tails sort_obs_inc spread_restoration bounded_below bounded_above lower_bound upper_bound +QTY_STATE_VARIABLE .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RHF .false. .false. .false. .false. .false. .false. -888888.0 -888888.0 +QTY_AQUIFER_WATER .false. .false. 10.0 0.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RHF .false. .false. .false. .false. .false. .false. -888888.0 -888888.0 +QTY_SEAICE_SALINITY008 .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RHF .false. .false. .false. .false. .false. .false. -888888.0 -888888.0 +QTY_VEGETATED_AREA_FRACTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RHF .false. .false. .false. .false. .false. .false. -888888.0 -888888.0 diff --git a/developer_tests/qceff/work/runall.sh b/developer_tests/qceff/work/runall.sh index 80b3fc87b1..0588ec6d69 100755 --- a/developer_tests/qceff/work/runall.sh +++ b/developer_tests/qceff/work/runall.sh @@ -43,3 +43,5 @@ fi ./test_table_read qcf_table_lower_gt_upper.txt ; should_fail "upper bound less than lower" ./test_table_read ./test_table_read qcf_table_lower_bound_only.txt ; should_pass "lower bound only" + +./test_table_read qcf_table_no_bounds_with_values.txt ; should_pass "bounds false, values for bounds" From 22f86c418e9e6f036b18ec8813f39335ce4cb675 Mon Sep 17 00:00:00 2001 From: Helen Kershaw Date: Fri, 29 Sep 2023 12:33:15 -0400 Subject: [PATCH 57/60] fix: need to check that a qty is bounded above and below before checking invalid bounds otherwise missing_r8 -88888 value for the upper bound is "less than" the lower bound --- .../assimilation/algorithm_info_mod.f90 | 43 ++++++++++++------- 1 file changed, 27 insertions(+), 16 deletions(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index 2a01884641..d235c289b4 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -518,26 +518,37 @@ subroutine verify_qcf_table_data() !Checks that all bounds are valid; currently checks that the lower bound in less than the upper !Here we could add more specific checks if we have known limits on the bounds do row = 1, size(qcf_table_data) - if(qcf_table_data(row)%obs_error_info%lower_bound > qcf_table_data(row)%obs_error_info%upper_bound) then - write(errstring,*) 'Invalid bounds in obs_error_info' - call error_handler(E_ERR, 'verify_qcf_table_data:', errstring, source) - endif - if(qcf_table_data(row)%probit_inflation%lower_bound > qcf_table_data(row)%probit_inflation%upper_bound) then - write(errstring,*) 'Invalid bounds in probit_inflation' - call error_handler(E_ERR, 'verify_qcf_table_data:', errstring, source) + + if (qcf_table_data(row)%obs_error_info%bounded_below .and. qcf_table_data(row)%obs_error_info%bounded_above) then + if(qcf_table_data(row)%obs_error_info%lower_bound > qcf_table_data(row)%obs_error_info%upper_bound) then + write(errstring,*) 'Invalid bounds in obs_error_info' + call error_handler(E_ERR, 'verify_qcf_table_data:', errstring, source) + endif endif - if(qcf_table_data(row)%probit_state%lower_bound > qcf_table_data(row)%probit_state%upper_bound) then - write(errstring,*) 'Invalid bounds in probit_state' - call error_handler(E_ERR, 'verify_qcf_table_data:', errstring, source) + if (qcf_table_data(row)%probit_inflation%bounded_below .and. qcf_table_data(row)%probit_inflation%bounded_above) then + if(qcf_table_data(row)%probit_inflation%lower_bound > qcf_table_data(row)%probit_inflation%upper_bound) then + write(errstring,*) 'Invalid bounds in probit_inflation' + call error_handler(E_ERR, 'verify_qcf_table_data:', errstring, source) + endif endif - if(qcf_table_data(row)%probit_extended_state%lower_bound > qcf_table_data(row)%probit_extended_state%upper_bound) then - write(errstring,*) 'Invalid bounds in probit_extended_state' - call error_handler(E_ERR, 'verify_qcf_table_data:', errstring, source) + if(qcf_table_data(row)%probit_state%bounded_below .and. qcf_table_data(row)%probit_state%bounded_above) then + if(qcf_table_data(row)%probit_state%lower_bound > qcf_table_data(row)%probit_state%upper_bound) then + write(errstring,*) 'Invalid bounds in probit_state' + call error_handler(E_ERR, 'verify_qcf_table_data:', errstring, source) + endif endif - if(qcf_table_data(row)%obs_inc_info%lower_bound > qcf_table_data(row)%obs_inc_info%upper_bound) then - write(errstring,*) 'Invalid bounds in obs_inc_info' - call error_handler(E_ERR, 'verify_qcf_table_data:', errstring, source) + if(qcf_table_data(row)%probit_extended_state%bounded_below .and. qcf_table_data(row)%probit_extended_state%bounded_above) then + if(qcf_table_data(row)%probit_extended_state%lower_bound > qcf_table_data(row)%probit_extended_state%upper_bound) then + write(errstring,*) 'Invalid bounds in probit_extended_state' + call error_handler(E_ERR, 'verify_qcf_table_data:', errstring, source) + endif endif + if(qcf_table_data(row)%obs_inc_info%bounded_below .and. qcf_table_data(row)%obs_inc_info%bounded_above) then + if(qcf_table_data(row)%obs_inc_info%lower_bound > qcf_table_data(row)%obs_inc_info%upper_bound) then + write(errstring,*) 'Invalid bounds in obs_inc_info' + call error_handler(E_ERR, 'verify_qcf_table_data:', errstring, source) + endif + endif end do !Ensures that all QTYs listed in the table exist in DART From 38353bbac9b036243a9a671631e05a7fdc891430 Mon Sep 17 00:00:00 2001 From: Helen Kershaw Date: Fri, 29 Sep 2023 13:45:07 -0400 Subject: [PATCH 58/60] fix: remove extra call to test_table read from runall.sh --- developer_tests/qceff/work/runall.sh | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/developer_tests/qceff/work/runall.sh b/developer_tests/qceff/work/runall.sh index 0588ec6d69..4202947a4c 100755 --- a/developer_tests/qceff/work/runall.sh +++ b/developer_tests/qceff/work/runall.sh @@ -42,6 +42,6 @@ fi ./test_table_read qcf_table_lower_gt_upper.txt ; should_fail "upper bound less than lower" -./test_table_read ./test_table_read qcf_table_lower_bound_only.txt ; should_pass "lower bound only" +./test_table_read qcf_table_lower_bound_only.txt ; should_pass "lower bound only" ./test_table_read qcf_table_no_bounds_with_values.txt ; should_pass "bounds false, values for bounds" From a4723e7b706dd4be054e0e103de9703514ebda8a Mon Sep 17 00:00:00 2001 From: Helen Kershaw Date: Mon, 2 Oct 2023 12:59:41 -0400 Subject: [PATCH 59/60] replace kind with qty kind is outdated terminolgy for quantity https://github.com/NCAR/DART/pull/545#issuecomment-1743365998 --- .../assimilation/algorithm_info_mod.f90 | 28 +++++++++---------- 1 file changed, 14 insertions(+), 14 deletions(-) diff --git a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 index d235c289b4..766ec4d6f8 100644 --- a/assimilation_code/modules/assimilation/algorithm_info_mod.f90 +++ b/assimilation_code/modules/assimilation/algorithm_info_mod.f90 @@ -199,21 +199,21 @@ subroutine obs_error_info(obs_def, error_variance, & logical, intent(out) :: bounded_below, bounded_above real(r8), intent(out) :: lower_bound, upper_bound -integer :: obs_type, obs_kind +integer :: obs_type, obs_qty integer(i8) :: state_var_index type(location_type) :: temp_loc integer :: QTY_loc(1) -character(len=129) :: kind_name +character(len=129) :: qty_name ! Get the kind of the observation obs_type = get_obs_def_type_of_obs(obs_def) ! If it is negative, it is an identity obs if(obs_type < 0) then state_var_index = -1 * obs_type - call get_state_meta_data(state_var_index, temp_loc, obs_kind) + call get_state_meta_data(state_var_index, temp_loc, obs_qty) else - obs_kind = get_quantity_for_type_of_obs(obs_type) + obs_qty = get_quantity_for_type_of_obs(obs_type) endif ! Get the default error variance @@ -227,10 +227,10 @@ subroutine obs_error_info(obs_def, error_variance, & endif !get actual name of QTY from integer index -kind_name = get_name_for_quantity(obs_kind) +qty_name = get_name_for_quantity(obs_qty) !find location of QTY in qcf_table_data structure -QTY_loc = findloc(qcf_table_row_headers, kind_name) +QTY_loc = findloc(qcf_table_row_headers, qty_name) if (QTY_loc(1) == 0) then !use default values if QTY is not in table @@ -264,7 +264,7 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & real(r8), intent(out) :: lower_bound, upper_bound integer :: QTY_loc(1) -character(len=129) :: kind_name +character(len=129) :: qty_name integer :: dist_type_loc(1) character(len=129), dimension(7) :: possible_dist_types @@ -316,10 +316,10 @@ subroutine probit_dist_info(kind, is_state, is_inflation, dist_type, & endif !get actual name of QTY from integer index -kind_name = get_name_for_quantity(kind) +qty_name = get_name_for_quantity(kind) !find location of QTY in qcf_table_data structure -QTY_loc = findloc(qcf_table_row_headers, kind_name) +QTY_loc = findloc(qcf_table_row_headers, qty_name) if (QTY_loc(1) == 0) then !use default values if QTY is not in table @@ -397,10 +397,10 @@ end subroutine probit_dist_info !------------------------------------------------------------------------ -subroutine obs_inc_info(obs_kind, filter_kind, rectangular_quadrature, gaussian_likelihood_tails, & +subroutine obs_inc_info(obs_qty, filter_kind, rectangular_quadrature, gaussian_likelihood_tails, & sort_obs_inc, spread_restoration, bounded_below, bounded_above, lower_bound, upper_bound) -integer, intent(in) :: obs_kind +integer, intent(in) :: obs_qty integer, intent(inout) :: filter_kind logical, intent(inout) :: rectangular_quadrature, gaussian_likelihood_tails logical, intent(inout) :: sort_obs_inc @@ -409,7 +409,7 @@ subroutine obs_inc_info(obs_kind, filter_kind, rectangular_quadrature, gaussian_ real(r8), intent(inout) :: lower_bound, upper_bound integer :: QTY_loc(1) -character(len=129) :: kind_name +character(len=129) :: qty_name integer :: filter_kind_loc(1) character(len=129), dimension(5) :: possible_filter_kinds @@ -446,10 +446,10 @@ subroutine obs_inc_info(obs_kind, filter_kind, rectangular_quadrature, gaussian_ endif !get actual name of QTY from integer index -kind_name = get_name_for_quantity(obs_kind) +qty_name = get_name_for_quantity(obs_qty) !find location of QTY in qcf_table_data structure -QTY_loc = findloc(qcf_table_row_headers, kind_name) +QTY_loc = findloc(qcf_table_row_headers, qty_name) if (QTY_loc(1) == 0) then !use default values if QTY is not in table From bd4dbab73a9e0cf6a20ebf2cd03924e5a31c4b43 Mon Sep 17 00:00:00 2001 From: Helen Kershaw Date: Mon, 2 Oct 2023 14:12:53 -0400 Subject: [PATCH 60/60] add two tests for table read: incorrect filter kind, incorrect distribution Currently the algorithm_info_mod does not catch these on table read --- .../qceff/work/qcf_table_incorrect_distribution.txt | 3 +++ .../qceff/work/qcf_table_incorrect_filter_kind.txt | 3 +++ developer_tests/qceff/work/runall.sh | 4 ++++ 3 files changed, 10 insertions(+) create mode 100644 developer_tests/qceff/work/qcf_table_incorrect_distribution.txt create mode 100644 developer_tests/qceff/work/qcf_table_incorrect_filter_kind.txt diff --git a/developer_tests/qceff/work/qcf_table_incorrect_distribution.txt b/developer_tests/qceff/work/qcf_table_incorrect_distribution.txt new file mode 100644 index 0000000000..37decd57bf --- /dev/null +++ b/developer_tests/qceff/work/qcf_table_incorrect_distribution.txt @@ -0,0 +1,3 @@ +QCF table version: 1 +QTY_TEMPLATE: bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound filter_kind rectangular_quadrature gaussian_likelihood_tails sort_obs_inc spread_restoration bounded_below bounded_above lower_bound upper_bound +QTY_STATE_VARIABLE .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 POLAR_BEAR_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RHF .false. .false. .false. .false. .false. .false. -888888.0 -888888.0 diff --git a/developer_tests/qceff/work/qcf_table_incorrect_filter_kind.txt b/developer_tests/qceff/work/qcf_table_incorrect_filter_kind.txt new file mode 100644 index 0000000000..c1125c3360 --- /dev/null +++ b/developer_tests/qceff/work/qcf_table_incorrect_filter_kind.txt @@ -0,0 +1,3 @@ +QCF table version: 1 +QTY_TEMPLATE: bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound dist_type bounded_below bounded_above lower_bound upper_bound filter_kind rectangular_quadrature gaussian_likelihood_tails sort_obs_inc spread_restoration bounded_below bounded_above lower_bound upper_bound +QTY_STATE_VARIABLE .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 BOUNDED_NORMAL_RH_DISTRIBUTION .false. .false. -888888.0 -888888.0 PENGUIN_FILTER .false. .false. .false. .false. .false. .false. -888888.0 -888888.0 diff --git a/developer_tests/qceff/work/runall.sh b/developer_tests/qceff/work/runall.sh index 4202947a4c..62d42c3fc9 100755 --- a/developer_tests/qceff/work/runall.sh +++ b/developer_tests/qceff/work/runall.sh @@ -45,3 +45,7 @@ fi ./test_table_read qcf_table_lower_bound_only.txt ; should_pass "lower bound only" ./test_table_read qcf_table_no_bounds_with_values.txt ; should_pass "bounds false, values for bounds" + +./test_table_read qcf_table_incorrect_filter_kind.txt ; should_fail "incorrect filter_kind" + +./test_table_read qcf_table_incorrect_distribution.txt ; should_fail "incorrect distribution"