diff --git a/CHANGELOG.md b/CHANGELOG.md index c93341f1..703f7313 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -11,6 +11,8 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ### Added +- Added optional NetCDF4 output mode for ObsFcstAna, including NetCDF metadata and runtime context. Changed namelist variable "out_ObsFcstAna" from logical to integer. + ### Changed - Revised and cleaned up RESTART options: diff --git a/GEOSlandassim_GridComp/GEOS_LandAssimGridComp.F90 b/GEOSlandassim_GridComp/GEOS_LandAssimGridComp.F90 index c7904fca..5fe2d892 100644 --- a/GEOSlandassim_GridComp/GEOS_LandAssimGridComp.F90 +++ b/GEOSlandassim_GridComp/GEOS_LandAssimGridComp.F90 @@ -98,14 +98,14 @@ module GEOS_LandAssimGridCompMod real :: fcsterr_inflation_fac integer :: N_obs_param logical :: out_obslog - logical :: out_ObsFcstAna + integer :: out_ObsFcstAna logical :: out_smapL4SMaup integer :: N_obsbias_max integer, dimension(:), pointer :: N_catl_vec,low_ind integer :: N_catf - !reordered tile_coord_rf and mapping l2rf + ! reordered tile_coord_rf and mapping l2rf integer, dimension(:), pointer :: l2rf, rf2l,rf2g, rf2f type(tile_coord_type), dimension(:), pointer :: tile_coord_rf => null() @@ -1407,7 +1407,7 @@ subroutine Initialize(gc, import, export, clock, rc) call MPI_BCAST(fcsterr_inflation_fac, 1, MPI_REAL, 0,MPICOMM,mpierr) call MPI_BCAST(N_obs_param, 1, MPI_INTEGER, 0,MPICOMM,mpierr) call MPI_BCAST(out_obslog, 1, MPI_LOGICAL, 0,MPICOMM,mpierr) - call MPI_BCAST(out_ObsFcstAna, 1, MPI_LOGICAL, 0,MPICOMM,mpierr) + call MPI_BCAST(out_ObsFcstAna, 1, MPI_INTEGER, 0,MPICOMM,mpierr) call MPI_BCAST(out_smapL4SMaup, 1, MPI_LOGICAL, 0,MPICOMM,mpierr) call MPI_BCAST(N_obsbias_max, 1, MPI_INTEGER, 0,MPICOMM,mpierr) @@ -1953,6 +1953,7 @@ subroutine RUN ( GC, IMPORT, EXPORT, CLOCK, RC ) N_catl, tile_coord_l, & N_catf, tile_coord_rf, tcinternal%pgrid_g, & N_catl_vec, low_ind, rf2l, & + update_type, LandAssimDTstep, & obs_param, & met_force, lai, & cat_param, cat_progn, mwRTM_param, & diff --git a/GEOSlandassim_GridComp/clsm_ensupd_enkf_update.F90 b/GEOSlandassim_GridComp/clsm_ensupd_enkf_update.F90 index b304ac01..d9cb8278 100644 --- a/GEOSlandassim_GridComp/clsm_ensupd_enkf_update.F90 +++ b/GEOSlandassim_GridComp/clsm_ensupd_enkf_update.F90 @@ -9,6 +9,9 @@ module clsm_ensupd_enkf_update + use MAPL_BaseMod, ONLY: & + MAPL_UNDEF + use MAPL_SortMod, ONLY: & MAPL_Sort @@ -30,7 +33,8 @@ module clsm_ensupd_enkf_update logit, & logunit, & nodata_generic, & - nodata_tolfrac_generic + nodata_tolfrac_generic, & + LDAS_is_nodata use LDAS_DateTimeMod, ONLY: & date_time_type, & @@ -1512,7 +1516,8 @@ end subroutine apply_enkf_increments ! ******************************************************************** subroutine output_ObsFcstAna(date_time, exp_id, & - N_obsl, Observations_l, N_obs_param, rf2f) + N_obsl, Observations_l, N_obs_param, obs_param, out_ObsFcstAna, & + update_type, dtstep_assim, rf2f) ! obs space output: observations, obs space forecast, obs space analysis, and ! associated error variances @@ -1526,11 +1531,14 @@ subroutine output_ObsFcstAna(date_time, exp_id, & character(*), intent(in) :: exp_id integer, intent(in) :: N_obsl, N_obs_param - + type(obs_param_type), dimension(N_obs_param), intent(in) :: obs_param + integer, intent(in) :: out_ObsFcstAna + integer, intent(in) :: update_type + integer, intent(in) :: dtstep_assim type(obs_type), dimension(N_obsl), intent(in) :: Observations_l - integer, dimension(:), optional, intent(in) :: rf2f + integer, dimension(:), optional, intent(in) :: rf2f ! re-ordered to LDASsa ! --------------------- @@ -1547,6 +1555,9 @@ subroutine output_ObsFcstAna(date_time, exp_id, & integer, dimension(numprocs) :: N_obsl_vec, tmp_low_ind character(300) :: fname + character( 40) :: file_ext + character(len=*), parameter :: Iam = 'output_ObsFcstAna' + character(len=400) :: err_msg #ifdef LDAS_MPI @@ -1689,60 +1700,366 @@ subroutine output_ObsFcstAna(date_time, exp_id, & ! write to file - fname = get_io_filename( './', exp_id, file_tag, date_time=date_time, & - dir_name=dir_name, ens_id=-1, no_subdirs=.true. ) + if (out_ObsFcstAna == 1 .or. out_ObsFcstAna == 3) then + + fname = get_io_filename( './', exp_id, file_tag, date_time=date_time, & + dir_name=dir_name, ens_id=-1, no_subdirs=.true. ) - open( 10, file=fname, form='unformatted', action='write') + open( 10, file=fname, form='unformatted', action='write') - ! write header + ! write header - write (10) N_obsf, date_time%year, date_time%month, & - date_time%day, date_time%hour, date_time%min, date_time%sec, & - date_time%dofyr, date_time%pentad + write (10) N_obsf, date_time%year, date_time%month, & + date_time%day, date_time%hour, date_time%min, date_time%sec, & + date_time%dofyr, date_time%pentad - ! write data + ! write data - ! Assuming a linear model and uncorrelated obs/model errors, - ! - ! the expected var of OminusF is HPHt + R, and - ! the expected var of OminusA is R - HAHt, where - ! - ! P = prior state error covariance - ! A = posterior state error covariance - ! H = observation operator + ! Assuming a linear model and uncorrelated obs/model errors, + ! + ! the expected var of OminusF is HPHt + R, and + ! the expected var of OminusA is R - HAHt, where + ! + ! P = prior state error covariance + ! A = posterior state error covariance + ! H = observation operator - write (10) (Observations_f(n)%assim, n=1,N_obsf) - write (10) (Observations_f(n)%species, n=1,N_obsf) + write (10) (Observations_f(n)%assim, n=1,N_obsf) + write (10) (Observations_f(n)%species, n=1,N_obsf) - write (10) (Observations_f(n)%tilenum, n=1,N_obsf) + write (10) (Observations_f(n)%tilenum, n=1,N_obsf) - write (10) (Observations_f(n)%lon, n=1,N_obsf) - write (10) (Observations_f(n)%lat, n=1,N_obsf) + write (10) (Observations_f(n)%lon, n=1,N_obsf) + write (10) (Observations_f(n)%lat, n=1,N_obsf) - write (10) (Observations_f(n)%obs, n=1,N_obsf) - write (10) (Observations_f(n)%obsvar, n=1,N_obsf) ! R + write (10) (Observations_f(n)%obs, n=1,N_obsf) + write (10) (Observations_f(n)%obsvar, n=1,N_obsf) ! R - write (10) (Observations_f(n)%fcst, n=1,N_obsf) - write (10) (Observations_f(n)%fcstvar, n=1,N_obsf) ! HPHt + write (10) (Observations_f(n)%fcst, n=1,N_obsf) + write (10) (Observations_f(n)%fcstvar, n=1,N_obsf) ! HPHt - write (10) (Observations_f(n)%ana, n=1,N_obsf) - write (10) (Observations_f(n)%anavar, n=1,N_obsf) ! HAHt + write (10) (Observations_f(n)%ana, n=1,N_obsf) + write (10) (Observations_f(n)%anavar, n=1,N_obsf) ! HAHt - close(10,status='keep') + close(10,status='keep') + + end if + + if (out_ObsFcstAna == 2 .or. out_ObsFcstAna == 3) then + + file_ext = '.nc4' + + fname = get_io_filename( './', exp_id, file_tag, date_time=date_time, & + dir_name=dir_name, ens_id=-1, file_ext=file_ext, no_subdirs=.true. ) + + call write_ObsFcstAna_nc4( fname, exp_id, N_obsf, Observations_f, & + N_obs_param, obs_param, update_type, dtstep_assim ) + + end if end if + if (allocated(Observations_f)) deallocate(Observations_f) - + end subroutine output_ObsFcstAna - ! ********************************************************************** + ! ******************************************************************** + + subroutine write_ObsFcstAna_nc4(fname, exp_id, N_obsf, Observations_f, & + N_obs_param, obs_param, update_type, dtstep_assim ) + + use netcdf + use pfio_NetCDF_Supplement, only: pfio_nf90_put_var_string + + implicit none + + character(*), intent(in) :: fname, exp_id + integer, intent(in) :: N_obsf, N_obs_param, update_type, dtstep_assim + type(obs_type), dimension(N_obsf), intent(in) :: Observations_f + type(obs_param_type), dimension(N_obs_param), intent(in) :: obs_param + ! ---------------- + + integer :: ncid + integer :: nobs_dimid + integer :: nspecies_dimid + integer :: assim_flag_varid, species_varid, tilenum_varid + integer :: lon_varid, lat_varid + integer :: obs_varid, obsvar_varid + integer :: fcst_varid, fcstvar_varid + integer :: ana_varid, anavar_varid + integer :: obsparam_species_id_varid, obsparam_assim_varid, obsparam_scale_varid + integer :: obsparam_errstd_varid + integer :: obsparam_varname_varid, obsparam_units_varid, obsparam_descr_varid + integer :: obsparam_fcstvarname_varid, obsparam_fcstunits_varid + + integer, dimension(:), allocatable :: species_assim_int, species_scale_int + real, dimension(:), allocatable :: species_errstd_r4 + character(len=40), dimension(:), allocatable :: species_varname_s, species_units_s + character(len=40), dimension(:), allocatable :: species_fcstvarname_s, species_fcstunits_s + character(len=40), dimension(:), allocatable :: species_descr_s + character(len=40) :: attr_name + character(len=8) :: write_date_yyyymmdd + character(len=10) :: write_time_hhmmss + character(len=5) :: write_zone + character(len=24) :: write_datetime_iso + character(len=64) :: user_name + character(len=128) :: created_by + integer :: user_len, user_status + integer :: i + + integer, dimension(N_obsf) :: tmpvecint + real, dimension(N_obsf) :: tmpvecreal + + character(len=*), parameter :: Iam = 'write_ObsFcstAna_nc4' + character(len=400) :: err_msg + + call nc4_check( nf90_create( trim(fname), nf90_clobber + NF90_NETCDF4, ncid ) ) + + call nc4_check( nf90_def_dim(ncid, 'n_obs', N_obsf, nobs_dimid) ) + call nc4_check( nf90_def_dim(ncid, 'n_species', N_obs_param, nspecies_dimid) ) + + call nc4_check( nf90_def_var(ncid, 'assim_flag', NF90_INT, [nobs_dimid], assim_flag_varid) ) + call nc4_check( nf90_def_var(ncid, 'species', NF90_INT, [nobs_dimid], species_varid) ) + call nc4_check( nf90_def_var(ncid, 'tilenum', NF90_INT, [nobs_dimid], tilenum_varid) ) + call nc4_check( nf90_def_var(ncid, 'lon', NF90_FLOAT, [nobs_dimid], lon_varid) ) + call nc4_check( nf90_def_var(ncid, 'lat', NF90_FLOAT, [nobs_dimid], lat_varid) ) + call nc4_check( nf90_def_var(ncid, 'obs', NF90_FLOAT, [nobs_dimid], obs_varid) ) + call nc4_check( nf90_def_var(ncid, 'obsvar', NF90_FLOAT, [nobs_dimid], obsvar_varid) ) + call nc4_check( nf90_def_var(ncid, 'fcst', NF90_FLOAT, [nobs_dimid], fcst_varid) ) + call nc4_check( nf90_def_var(ncid, 'fcstvar', NF90_FLOAT, [nobs_dimid], fcstvar_varid) ) + call nc4_check( nf90_def_var(ncid, 'ana', NF90_FLOAT, [nobs_dimid], ana_varid) ) + call nc4_check( nf90_def_var(ncid, 'anavar', NF90_FLOAT, [nobs_dimid], anavar_varid) ) + + call nc4_check( nf90_def_var(ncid, 'obsparam_species_id', NF90_INT, [nspecies_dimid], obsparam_species_id_varid) ) + call nc4_check( nf90_def_var(ncid, 'obsparam_assim', NF90_INT, [nspecies_dimid], obsparam_assim_varid) ) + call nc4_check( nf90_def_var(ncid, 'obsparam_scale', NF90_INT, [nspecies_dimid], obsparam_scale_varid) ) + call nc4_check( nf90_def_var(ncid, 'obsparam_errstd', NF90_FLOAT, [nspecies_dimid], obsparam_errstd_varid) ) + call nc4_check( nf90_def_var(ncid, 'obsparam_varname', NF90_STRING, [nspecies_dimid], obsparam_varname_varid) ) + call nc4_check( nf90_def_var(ncid, 'obsparam_units', NF90_STRING, [nspecies_dimid], obsparam_units_varid) ) + call nc4_check( nf90_def_var(ncid, 'obsparam_fcstvarname', NF90_STRING, [nspecies_dimid], obsparam_fcstvarname_varid) ) + call nc4_check( nf90_def_var(ncid, 'obsparam_fcstunits', NF90_STRING, [nspecies_dimid], obsparam_fcstunits_varid) ) + call nc4_check( nf90_def_var(ncid, 'obsparam_descr', NF90_STRING, [nspecies_dimid], obsparam_descr_varid) ) + + ! get date and time when file is written + call date_and_time(write_date_yyyymmdd, write_time_hhmmss, write_zone) + write(write_datetime_iso, '(A4,A1,A2,A1,A2,A1,A2,A1,A2,A1,A2,A5)') & + write_date_yyyymmdd(1:4), '-', write_date_yyyymmdd(5:6), '-', & + write_date_yyyymmdd(7:8), 'T', write_time_hhmmss(1:2), ':', & + write_time_hhmmss(3:4), ':', write_time_hhmmss(5:6), write_zone + + ! get user info (if available) + call get_environment_variable('USER', user_name, length=user_len, status=user_status) + if (user_status == 0 .and. user_len > 0) then + created_by = trim(user_name(1:user_len)) // ' via write_ObsFcstAna_nc4()' + else + created_by = 'GEOSldas write_ObsFcstAna_nc4' + end if + + ! write attributes + + call nc4_check( nf90_put_att(ncid, NF90_GLOBAL, 'Comment', 'NetCDF-4') ) + call nc4_check( nf90_put_att(ncid, NF90_GLOBAL, 'Contact', '') ) + call nc4_check( nf90_put_att(ncid, NF90_GLOBAL, 'Conventions', 'CF') ) + call nc4_check( nf90_put_att(ncid, NF90_GLOBAL, 'Filename', trim(fname)) ) + call nc4_check( nf90_put_att(ncid, NF90_GLOBAL, 'Institution', 'NASA Global Modeling and Assimilation Office') ) + call nc4_check( nf90_put_att(ncid, NF90_GLOBAL, 'Source', 'Experiment_ID: ' // trim(exp_id)) ) + call nc4_check( nf90_put_att(ncid, NF90_GLOBAL, 'Title' , 'Observation-space,Single-Level,Assimilation,Land Surface Diagnostics') ) + + call nc4_check( nf90_put_att(ncid, NF90_GLOBAL, 'GEOSldas_update_type', update_type) ) + call nc4_check( nf90_put_att(ncid, NF90_GLOBAL, 'GEOSldas_assimilation_time_step_in_seconds', dtstep_assim) ) + + do i=1,N_obs_param + write(attr_name, '(A,I0.3,A)') 'GEOSldas_observation_species_', obs_param(i)%species, '_descr' + call nc4_check( nf90_put_att(ncid, NF90_GLOBAL, trim(attr_name), trim(obs_param(i)%descr)) ) + end do + call nc4_check( nf90_put_att(ncid, NF90_GLOBAL, 'schema_version', 'ObsFcstAna_nc4_v1') ) + call nc4_check( nf90_put_att(ncid, NF90_GLOBAL, 'DateCreated', trim(write_datetime_iso)) ) + call nc4_check( nf90_put_att(ncid, NF90_GLOBAL, 'CreatedBy', trim(created_by)) ) + + call nc4_check( nf90_put_att(ncid, assim_flag_varid, 'long_name', 'observation assimilation flag (0=not assimilated, 1=assimilated)') ) + call nc4_check( nf90_put_att(ncid, assim_flag_varid, 'units', '1') ) + call nc4_check( nf90_put_att(ncid, species_varid, 'long_name', 'observation species identifier') ) + call nc4_check( nf90_put_att(ncid, species_varid, 'units', '1') ) + call nc4_check( nf90_put_att(ncid, tilenum_varid, 'long_name', 'model tile number associated with observation location') ) + call nc4_check( nf90_put_att(ncid, tilenum_varid, 'units', '1') ) + + call nc4_check( nf90_put_att(ncid, lon_varid, 'long_name', 'observation longitude') ) + call nc4_check( nf90_put_att(ncid, lon_varid, 'standard_name', 'longitude') ) + call nc4_check( nf90_put_att(ncid, lon_varid, 'units', 'degrees_east') ) + call nc4_check( nf90_put_att(ncid, lat_varid, 'long_name', 'observation latitude') ) + call nc4_check( nf90_put_att(ncid, lat_varid, 'standard_name', 'latitude') ) + call nc4_check( nf90_put_att(ncid, lat_varid, 'units', 'degrees_north') ) + + call nc4_check( nf90_put_att(ncid, obs_varid, 'long_name', 'observation value (after scaling)') ) + call nc4_check( nf90_put_att(ncid, obs_varid, 'units', 'species-dependent') ) + call nc4_check( nf90_put_att(ncid, obs_varid, '_FillValue', MAPL_UNDEF) ) + call nc4_check( nf90_put_att(ncid, obs_varid, 'missing_value', MAPL_UNDEF) ) + call nc4_check( nf90_put_att(ncid, obsvar_varid, 'long_name', 'observation error variance (after scaling)') ) + call nc4_check( nf90_put_att(ncid, obsvar_varid, 'units', 'species-dependent') ) + call nc4_check( nf90_put_att(ncid, obsvar_varid, '_FillValue', MAPL_UNDEF) ) + call nc4_check( nf90_put_att(ncid, obsvar_varid, 'missing_value', MAPL_UNDEF) ) + call nc4_check( nf90_put_att(ncid, fcst_varid, 'long_name', 'observation-equivalent model forecast value') ) + call nc4_check( nf90_put_att(ncid, fcst_varid, 'units', 'species-dependent') ) + call nc4_check( nf90_put_att(ncid, fcst_varid, '_FillValue', MAPL_UNDEF) ) + call nc4_check( nf90_put_att(ncid, fcst_varid, 'missing_value', MAPL_UNDEF) ) + call nc4_check( nf90_put_att(ncid, fcstvar_varid, 'long_name', 'observation-equivalent model forecast error variance') ) + call nc4_check( nf90_put_att(ncid, fcstvar_varid, 'units', 'species-dependent') ) + call nc4_check( nf90_put_att(ncid, fcstvar_varid, '_FillValue', MAPL_UNDEF) ) + call nc4_check( nf90_put_att(ncid, fcstvar_varid, 'missing_value', MAPL_UNDEF) ) + call nc4_check( nf90_put_att(ncid, ana_varid, 'long_name', 'observation-equivalent model analysis value') ) + call nc4_check( nf90_put_att(ncid, ana_varid, 'units', 'species-dependent') ) + call nc4_check( nf90_put_att(ncid, ana_varid, '_FillValue', MAPL_UNDEF) ) + call nc4_check( nf90_put_att(ncid, ana_varid, 'missing_value', MAPL_UNDEF) ) + call nc4_check( nf90_put_att(ncid, anavar_varid, 'long_name', 'observation-equivalent model analysis error variance') ) + call nc4_check( nf90_put_att(ncid, anavar_varid, 'units', 'species-dependent') ) + call nc4_check( nf90_put_att(ncid, anavar_varid, '_FillValue', MAPL_UNDEF) ) + call nc4_check( nf90_put_att(ncid, anavar_varid, 'missing_value', MAPL_UNDEF) ) + + call nc4_check( nf90_put_att(ncid, obsparam_species_id_varid, 'long_name', 'species-level observation parameter: observation species identifier') ) + call nc4_check( nf90_put_att(ncid, obsparam_species_id_varid, 'units', '1') ) + call nc4_check( nf90_put_att(ncid, obsparam_assim_varid, 'long_name', 'species-level observation parameter: observation assimilation flag') ) + call nc4_check( nf90_put_att(ncid, obsparam_assim_varid, 'units', '1') ) + call nc4_check( nf90_put_att(ncid, obsparam_scale_varid, 'long_name', 'species-level observation parameter: observation scaling flag') ) + call nc4_check( nf90_put_att(ncid, obsparam_scale_varid, 'units', '1') ) + call nc4_check( nf90_put_att(ncid, obsparam_errstd_varid, 'long_name', 'species-level observation parameter: default observation error standard deviation (before scaling)') ) + call nc4_check( nf90_put_att(ncid, obsparam_errstd_varid, 'units', 'species-dependent') ) + + call nc4_check( nf90_put_att(ncid, obsparam_varname_varid, 'long_name', 'species-level observation parameter: observation native variable name (before scaling)') ) + call nc4_check( nf90_put_att(ncid, obsparam_units_varid, 'long_name', 'species-level observation parameter: observation native units (before scaling)') ) + call nc4_check( nf90_put_att(ncid, obsparam_fcstvarname_varid, 'long_name', 'species-level observation parameter: observation-equivalent model variable name') ) + call nc4_check( nf90_put_att(ncid, obsparam_fcstunits_varid, 'long_name', 'species-level observation parameter: observation-equivalent model units') ) + call nc4_check( nf90_put_att(ncid, obsparam_descr_varid, 'long_name', 'species-level observation parameter: observation description') ) + + call nc4_check( nf90_enddef(ncid) ) + + ! write variables + + if (N_obsf > 0) then + + call nc4_check( nf90_put_var(ncid, species_varid, Observations_f(1:N_obsf)%species) ) + call nc4_check( nf90_put_var(ncid, tilenum_varid, Observations_f(1:N_obsf)%tilenum) ) + call nc4_check( nf90_put_var(ncid, lon_varid, Observations_f(1:N_obsf)%lon) ) + call nc4_check( nf90_put_var(ncid, lat_varid, Observations_f(1:N_obsf)%lat) ) + + ! for assim flag, convert logical to integer + tmpvecint = 0 + where (Observations_f(1:N_obsf)%assim) + tmpvecint = 1 + end where + call nc4_check( nf90_put_var(ncid, assim_flag_varid, tmpvecint) ) + + ! for data fields, replace LDAS no-data-value with MAPL_UNDEF for consistency with MAPL HISTORY output + tmpvecreal = Observations_f(1:N_obsf)%obs + do i=1,N_obsf + if (LDAS_is_nodata(tmpvecreal(i))) tmpvecreal(i) = MAPL_UNDEF + end do + call nc4_check( nf90_put_var(ncid, obs_varid, tmpvecreal) ) + + tmpvecreal = Observations_f(1:N_obsf)%obsvar + do i=1,N_obsf + if (LDAS_is_nodata(tmpvecreal(i))) tmpvecreal(i) = MAPL_UNDEF + end do + call nc4_check( nf90_put_var(ncid, obsvar_varid, tmpvecreal) ) + + tmpvecreal = Observations_f(1:N_obsf)%fcst + do i=1,N_obsf + if (LDAS_is_nodata(tmpvecreal(i))) tmpvecreal(i) = MAPL_UNDEF + end do + call nc4_check( nf90_put_var(ncid, fcst_varid, tmpvecreal) ) + + tmpvecreal = Observations_f(1:N_obsf)%fcstvar + do i=1,N_obsf + if (LDAS_is_nodata(tmpvecreal(i))) tmpvecreal(i) = MAPL_UNDEF + end do + call nc4_check( nf90_put_var(ncid, fcstvar_varid, tmpvecreal) ) + + tmpvecreal = Observations_f(1:N_obsf)%ana + do i=1,N_obsf + if (LDAS_is_nodata(tmpvecreal(i))) tmpvecreal(i) = MAPL_UNDEF + end do + call nc4_check( nf90_put_var(ncid, ana_varid, tmpvecreal) ) + + tmpvecreal = Observations_f(1:N_obsf)%anavar + do i=1,N_obsf + if (LDAS_is_nodata(tmpvecreal(i))) tmpvecreal(i) = MAPL_UNDEF + end do + call nc4_check( nf90_put_var(ncid, anavar_varid, tmpvecreal) ) + + end if + + if (N_obs_param > 0) then + + allocate(species_assim_int( N_obs_param)) + allocate(species_scale_int( N_obs_param)) + allocate(species_errstd_r4( N_obs_param)) + allocate(species_varname_s( N_obs_param)) + allocate(species_units_s( N_obs_param)) + allocate(species_fcstvarname_s(N_obs_param)) + allocate(species_fcstunits_s( N_obs_param)) + allocate(species_descr_s( N_obs_param)) + + species_assim_int = 0 + species_scale_int = 0 + + ! convert logicals to integer + where (obs_param(1:N_obs_param)%assim) species_assim_int = 1 + where (obs_param(1:N_obs_param)%scale) species_scale_int = 1 + + do i=1,N_obs_param + species_errstd_r4(i) = obs_param(i)%errstd + species_varname_s(i) = trim(obs_param(i)%varname) + species_units_s(i) = trim(obs_param(i)%units ) + species_fcstvarname_s(i) = trim(obs_param(i)%fcstvarname) + species_fcstunits_s(i) = trim(obs_param(i)%fcstunits ) + species_descr_s(i) = trim(obs_param(i)%descr ) + end do + + call nc4_check( nf90_put_var( ncid, obsparam_species_id_varid, obs_param(1:N_obs_param)%species) ) + call nc4_check( nf90_put_var( ncid, obsparam_assim_varid, species_assim_int) ) + call nc4_check( nf90_put_var( ncid, obsparam_scale_varid, species_scale_int) ) + call nc4_check( nf90_put_var( ncid, obsparam_errstd_varid, species_errstd_r4) ) + call nc4_check( pfio_nf90_put_var_string(ncid, obsparam_varname_varid, species_varname_s) ) + call nc4_check( pfio_nf90_put_var_string(ncid, obsparam_units_varid, species_units_s) ) + call nc4_check( pfio_nf90_put_var_string(ncid, obsparam_fcstvarname_varid, species_fcstvarname_s) ) + call nc4_check( pfio_nf90_put_var_string(ncid, obsparam_fcstunits_varid, species_fcstunits_s) ) + call nc4_check( pfio_nf90_put_var_string(ncid, obsparam_descr_varid, species_descr_s) ) + + deallocate(species_assim_int) + deallocate(species_scale_int) + deallocate(species_errstd_r4) + deallocate(species_varname_s) + deallocate(species_units_s) + deallocate(species_fcstvarname_s) + deallocate(species_fcstunits_s) + deallocate(species_descr_s) + + end if + + call nc4_check( nf90_close(ncid) ) + + contains + + subroutine nc4_check(status) + integer, intent(in) :: status + + if (status /= nf90_noerr) then + err_msg = 'NetCDF error in ' // Iam // ': ' // trim(nf90_strerror(status)) + call ldas_abort(LDAS_GENERIC_ERROR, Iam, err_msg) + end if + end subroutine nc4_check + + end subroutine write_ObsFcstAna_nc4 + + ! ********************************************************************** + subroutine output_ObsFcstAna_wrapper( out_ObsFcstAna, & date_time, exp_id, & N_obsl, N_obs_param, N_ens, & N_catl, tile_coord_l, & N_catf, tile_coord_f, pert_grid_g, & N_catl_vec, low_ind, f2l, & + update_type, dtstep_assim, & obs_param, & met_force, lai, cat_param, cat_progn, mwRTM_param, & Observations_l, rf2f ) @@ -1756,8 +2073,7 @@ subroutine output_ObsFcstAna_wrapper( out_ObsFcstAna, & ! major revisions for new obs handling and MPI - logical, intent(in) :: out_ObsFcstAna - + integer, intent(in) :: out_ObsFcstAna type(date_time_type), intent(in) :: date_time @@ -1765,6 +2081,8 @@ subroutine output_ObsFcstAna_wrapper( out_ObsFcstAna, & integer, intent(in) :: N_obsl, N_obs_param, N_ens, N_catl, N_catf + integer, intent(in) :: update_type, dtstep_assim + type(tile_coord_type), dimension(:), pointer :: tile_coord_l ! input type(tile_coord_type), dimension(:), pointer :: tile_coord_f ! input @@ -1787,9 +2105,9 @@ subroutine output_ObsFcstAna_wrapper( out_ObsFcstAna, & type(mwRTM_param_type), dimension(N_catl), intent(in) :: mwRTM_param - type(obs_type), dimension(:), pointer :: Observations_l ! inout + type(obs_type), dimension(:), pointer :: Observations_l ! inout - integer, dimension(N_catf), optional, intent(in) :: rf2f ! re-ordered to LDASsa + integer, dimension(N_catf), optional, intent(in) :: rf2f ! re-ordered to LDASsa ! local variables @@ -1806,7 +2124,7 @@ subroutine output_ObsFcstAna_wrapper( out_ObsFcstAna, & ! output "O-A" (obs - analysis) whenever innovations are output - if (out_ObsFcstAna) then + if (out_ObsFcstAna > 0) then ! compute model forecast of observations @@ -1828,8 +2146,9 @@ subroutine output_ObsFcstAna_wrapper( out_ObsFcstAna, & ! write out model, observations, and "OminusA" information - call output_ObsFcstAna( date_time, exp_id, N_obsl, & - Observations_l(1:N_obsl), N_obs_param, rf2f=rf2f ) + call output_ObsFcstAna( date_time, exp_id, N_obsl, & + Observations_l(1:N_obsl), N_obs_param, obs_param, out_ObsFcstAna, & + update_type, dtstep_assim, rf2f=rf2f ) end if diff --git a/GEOSlandassim_GridComp/clsm_ensupd_upd_routines.F90 b/GEOSlandassim_GridComp/clsm_ensupd_upd_routines.F90 index 795346f9..15223c2f 100644 --- a/GEOSlandassim_GridComp/clsm_ensupd_upd_routines.F90 +++ b/GEOSlandassim_GridComp/clsm_ensupd_upd_routines.F90 @@ -235,7 +235,7 @@ subroutine read_ens_upd_inputs( & type(obs_param_type), dimension(:), pointer :: obs_param ! output logical, intent(out) :: out_obslog - logical, intent(out) :: out_ObsFcstAna + integer, intent(out) :: out_ObsFcstAna logical, intent(out) :: out_smapL4SMaup integer, intent(out) :: N_obsbias_max @@ -272,6 +272,8 @@ subroutine read_ens_upd_inputs( & character(len=400) :: err_msg character(len= 6) :: tmpstring6 logical :: file_exists + integer :: ios + character(len=400) :: iomsg ! ----------------------------------------------------------------- @@ -300,8 +302,12 @@ subroutine read_ens_upd_inputs( & if (logit) write (logunit,*) if (logit) write (logunit,'(400A)') 'reading *default* EnKF inputs from ' // trim(fname) if (logit) write (logunit,*) - - read (10, nml=ens_upd_inputs) + + read (10, nml=ens_upd_inputs, iostat=ios, iomsg=iomsg) + if (ios /= 0) then + err_msg = 'Error reading ens_upd_inputs. NOTE: "out_ObsFcstAna" must be integer. ' // trim(iomsg) + call ldas_abort(LDAS_GENERIC_ERROR, Iam, err_msg) + end if close(10,status='keep') @@ -325,7 +331,13 @@ subroutine read_ens_upd_inputs( & if (logit) write (logunit,'(400A)') 'reading *special* EnKF inputs from ' // trim(fname) if (logit) write (logunit,*) - read (10, nml=ens_upd_inputs) + read (10, nml=ens_upd_inputs, iostat=ios, iomsg=iomsg) + + ! if read error, exit gracefully; hint at possible cause (out_ObsFcstAna changed from logical to integer, Apr 2026) + if (ios /= 0) then + err_msg = 'Error reading ens_upd_inputs. NOTE: "out_ObsFcstAna" must be integer. ' // trim(iomsg) + call ldas_abort(LDAS_GENERIC_ERROR, Iam, err_msg) + end if close(10,status='keep') @@ -337,15 +349,63 @@ subroutine read_ens_upd_inputs( & ! ! none implemented so far (reichle, 19 Jul 2005) + + ! ----------------------------------------------------------------- + + ! fill obs_param%fcstvarname and obs_param%fcstunits + + do i=1,N_obs_species_nml + + ! expect "fcstvarname" and "fcstunits" to be 'NULL'; they are not meant to be filled by the user in the config nml file + + if (trim(obs_param_nml(i)%fcstvarname) /= 'NULL' .and. trim(obs_param_nml(i)%fcstvarname) /= 'null') & + call ldas_abort(LDAS_GENERIC_ERROR, Iam, 'obs_param_nml%fcstvarname must be NULL on input') + + if (trim(obs_param_nml(i)%fcstunits ) /= 'NULL' .and. trim(obs_param_nml(i)%fcstunits ) /= 'null') & + call ldas_abort(LDAS_GENERIC_ERROR, Iam, 'obs_param_nml%fcstunits must be NULL on input') + + ! IMPORTANT: Must maintain consistency in the mapping between obs and model variables + ! that is encoded here with that in get_obs_pred(). + + select case (trim(obs_param_nml(i)%varname)) + + case ('sfmc', 'sfds') + + ! NOTE: 'sfds' is scaled into volumetric units of sfmc + + obs_param_nml(i)%fcstvarname = 'sfmc' + obs_param_nml(i)%fcstunits = 'm3 m-3' + + case ('rzmc', 'tsurf', 'FT', 'Tb', 'asnow', 'NULL') + + obs_param_nml(i)%fcstvarname = obs_param_nml(i)%varname + obs_param_nml(i)%fcstunits = obs_param_nml(i)%units + + case default + + call ldas_abort(LDAS_GENERIC_ERROR, Iam, 'unknown obs_param_nml%varname') + + end select + + end do + ! ----------------------------------------------------------------- ! ! consistency checks etc + + if (out_ObsFcstAna<0 .or. out_ObsFcstAna>3) then + write(tmpstring6,'(I6)') out_ObsFcstAna + err_msg = 'Unknown integer value of "out_ObsFcstAna": "' // trim(tmpstring6) // '". Use 0=OFF, 1=BIN, 2=NC4, 3=BIN and NC4.' + call ldas_abort(LDAS_GENERIC_ERROR, Iam, err_msg) + end if if (update_type==0) then err_msg = 'executable was built for assimilation but update_type=0' call ldas_abort(LDAS_GENERIC_ERROR, Iam, err_msg) end if + ! obs bias init and checks + N_obsbias_max = 0 ! initialize do i=1,N_obs_species_nml @@ -1166,6 +1226,9 @@ subroutine get_obs_pred( & do i=1,N_obs_param + ! IMPORTANT: Must maintain consistency in the mapping between obs and model ("lH") variables + ! that is encoded here with that in read_ens_upd_inputs(). + select case (trim(obs_param(i)%varname)) case ('sfmc', 'sfds') diff --git a/GEOSldas_App/LDASsa_DEFAULT_inputs_ensupd.nml b/GEOSldas_App/LDASsa_DEFAULT_inputs_ensupd.nml index 0a31e143..13ac0538 100644 --- a/GEOSldas_App/LDASsa_DEFAULT_inputs_ensupd.nml +++ b/GEOSldas_App/LDASsa_DEFAULT_inputs_ensupd.nml @@ -43,7 +43,7 @@ update_type = 0 out_obslog = .true. -out_ObsFcstAna = .false. +out_ObsFcstAna = 0 ! 0=OFF, 1=BIN (legacy unformatted), 2=NC4, 3=BIN and NC4 out_smapL4SMaup = .false. ! --------------------------------------------------------------------- @@ -85,55 +85,59 @@ fcsterr_inflation_fac = -9999. ! Definition of obs_param_nml fields (see also enkf_types.F90): ! ! %descr = description -! %species = identifier for type of measurement +! %species = identifier for type of observation ! %orbit = type of (half-)orbit ! 0 = n/a [eg., in situ obs] ! 1 = ascending ! 2 = descending ! 3 = ascending or descending ! 4 = geostationary +! ! %pol = polarization ! 0 = n/a [eg., multi-pol. retrieval] ! 1 = horizontal ! 2 = vertical -! 3 = ... +! 3 = ... [could add 3rd/4th Stokes, HH, HV, VH, VV] +! ! %N_ang = # satellite viewing angles in species (radiance obs only) ! %ang = vector of satellite viewing angles ! %freq = frequency [Hz] ! %FOV = field-of-view *radius*, see NOTES below ! (if FOV==0. equate obs footprint w/ tile) -! %FOV_units = field-of-view units ('km' or 'deg'), see NOTES below -! %assim = Should this obs type be assimilated (state update)? (logical) -! %scale = Should this obs be scaled? (logical) -! %getinnov = Should innov be computed for this obs type (logical) -! (innovations are always computed if assim==.true.) -! %RTM_ID = ID of radiative transfer model to use for Tb forward modeling -! (subroutine get_obs_pred()) +! %FOV_units = FOV units ('km' or 'deg'), see NOTES below +! %assim = assimilate obs species: true/false (see also "obs_type") +! %scale = scale obs species: true/false +! %getinnov = compute innovs? (set to .T. if assim==.T.) +! %RTM_ID = ID of radiative transfer model to use for Tb forward modeling (see get_obs_pred()) ! 0 = none ! 1 = L-band tau-omega model as in De Lannoy et al. 2013 (doi:10.1175/JHM-D-12-092.1) (old SMOS preproc) ! 2 = same as 1 but without Pellarin atm corr (SMAP) ! 3 = same as 1 but with Mironov and SMAP L2_SM pol mixing -! 4 = same as 3 but without Pellarin atm corr (targeted for SMAP L4_SM Version 8) -! %bias_Npar = number of obs bias states tracked per day (integer) +! 4 = same as 3 but without Pellarin atm corr (SMAP L4_SM Version 8) +! +! %bias_Npar = number of obs bias states tracked per day ! %bias_trel = e-folding time scale of obs bias memory [s] ! %bias_tcut = cutoff time for confident obs bias estimate [s] ! %nodata = no-data-value -! %varname = equivalent model variable name (for "Obs_pred") -! %units = units (eg., 'K' or 'm3/m3') -! %path = path to measurement files -! %name = name identifier for file containing measurements +! %varname = observation native variable name (before scaling) +! %units = observation native units (before scaling; eg., 'K' or 'm3/m3') +! %fcstvarname = observation-equivalent model variable name (Obs_pred) [do not edit, filled by read_ens_upd_inputs()] +! %fcstunits = observation-equivalent model units (eg., 'K' or 'm3/m3') [do not edit, filled by read_ens_upd_inputs()] +! %path = path to observations file +! %name = name identifier for file containing observations ! %maskpath = path to obs mask file ! %maskname = filename for obs mask ! %scalepath = path to file(s) with scaling parameters ! %scalename = filename for scaling parameters ! %flistpath = path to file with list of obs file names ! %flistname = name of file with list of obs file names -! %errstd = default obs error std -! %std_normal_max = maximum allowed perturbation (relative to N(0,1)) -! %zeromean = enforce zero mean across ensemble -! %coarsen_pert = generate obs perturbations on coarser grid (see pert_param_type%coarsen) -! %xcorr = correlation length (deg) in longitude direction -! %ycorr = correlation length (deg) in latitude direction +! %errstd = default obs error std (before scaling) +! %std_normal_max = maximum allowed obs perturbation (relative to N(0,1)) (see pert_param_type) +! %zeromean = enforce zero mean across ensemble (see pert_param_type) +! %coarsen_pert = generate obs perturbations on coarser grid (see pert_param_type) +! %xcorr = obs error correlation length (deg) in longitude direction (see pert_param_type) +! %ycorr = obs error correlation length (deg) in latitude direction (see pert_param_type) +! %adapt = identifier for adaptive filtering ! ! For observation perturbations, always use: ! @@ -257,6 +261,8 @@ obs_param_nml( 1)%bias_tcut = 432000 obs_param_nml( 1)%nodata = -9999. obs_param_nml( 1)%varname = 'sfmc' obs_param_nml( 1)%units = 'm3/m3' +obs_param_nml( 1)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml( 1)%fcstunits = 'NULL' ! do not edit! obs_param_nml( 1)%path = '/land/l_data/AMSR/data/AMSR_E_L2_Land_V001/' obs_param_nml( 1)%name = 'AMSR_E_L2_Land_' obs_param_nml( 1)%maskpath = '' @@ -294,6 +300,8 @@ obs_param_nml( 2)%bias_tcut = 432000 obs_param_nml( 2)%nodata = -9999. obs_param_nml( 2)%varname = 'sfmc' obs_param_nml( 2)%units = 'm3/m3' +obs_param_nml( 2)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml( 2)%fcstunits = 'NULL' ! do not edit! obs_param_nml( 2)%path = '/land/l_data/AMSR/data/AMSR_E_L2_Land_V001/' obs_param_nml( 2)%name = 'AMSR_E_L2_Land_' obs_param_nml( 2)%maskpath = '' @@ -331,6 +339,8 @@ obs_param_nml( 3)%bias_tcut = 432000 obs_param_nml( 3)%nodata = -9999. obs_param_nml( 3)%varname = 'tsurf' obs_param_nml( 3)%units = 'K' +obs_param_nml( 3)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml( 3)%fcstunits = 'NULL' ! do not edit! obs_param_nml( 3)%path = '/land/l_data/ISCCP/GSWP2_1by1_V1/' obs_param_nml( 3)%name = 'isccpdx_tskin.' obs_param_nml( 3)%maskpath = '' @@ -368,6 +378,8 @@ obs_param_nml( 4)%bias_tcut = 432000 obs_param_nml( 4)%nodata = -999. obs_param_nml( 4)%varname = 'sfmc' obs_param_nml( 4)%units = 'm3/m3' +obs_param_nml( 4)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml( 4)%fcstunits = 'NULL' ! do not edit! obs_param_nml( 4)%path = '/land/l_data/RedArk/Retrievals_36km/retrievals_20060508/' obs_param_nml( 4)%name = 'SM_retrieval.' obs_param_nml( 4)%maskpath = '' @@ -405,6 +417,8 @@ obs_param_nml( 5)%bias_tcut = 432000 obs_param_nml( 5)%nodata = -999. obs_param_nml( 5)%varname = 'sfmc' obs_param_nml( 5)%units = 'm3/m3' +obs_param_nml( 5)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml( 5)%fcstunits = 'NULL' ! do not edit! obs_param_nml( 5)%path = '/land/l_data/RedArk/Truth/50mm_Soil_Moisture_Truth/' obs_param_nml( 5)%name = 'red_ark_50mm.sm.' obs_param_nml( 5)%maskpath = '' @@ -442,6 +456,8 @@ obs_param_nml( 6)%bias_tcut = 432000 obs_param_nml( 6)%nodata = -999. obs_param_nml( 6)%varname = 'rzmc' obs_param_nml( 6)%units = 'm3/m3' +obs_param_nml( 6)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml( 6)%fcstunits = 'NULL' ! do not edit! obs_param_nml( 6)%path = '/land/l_data/RedArk/Truth/400mm_Soil_Moisture_Truth/' obs_param_nml( 6)%name = 'red_ark_400mm.sm.' obs_param_nml( 6)%maskpath = '' @@ -479,6 +495,8 @@ obs_param_nml( 7)%bias_tcut = 432000 obs_param_nml( 7)%nodata = -9999. obs_param_nml( 7)%varname = 'sfmc' obs_param_nml( 7)%units = 'm3/m3' +obs_param_nml( 7)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml( 7)%fcstunits = 'NULL' ! do not edit! obs_param_nml( 7)%path = '/land/l_data/RedArk_OSSE/data/Retrievals_CLSM_synth/M0001_P0001_R0001_URI/std_synth_obs_0.020/' obs_param_nml( 7)%name = 'CLSM_synth_sm.' obs_param_nml( 7)%maskpath = '' @@ -516,6 +534,8 @@ obs_param_nml( 8)%bias_tcut = 432000 obs_param_nml( 8)%nodata = -9999. obs_param_nml( 8)%varname = 'sfmc' obs_param_nml( 8)%units = 'm3/m3' +obs_param_nml( 8)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml( 8)%fcstunits = 'NULL' ! do not edit! obs_param_nml( 8)%path = '/discover/nobackup/vmaggion/Synth_sfmc/radar_sim/std_synth_sfmc_0.040/' obs_param_nml( 8)%name = 'synth_sfmc_VivianaOK_' obs_param_nml( 8)%maskpath = '' @@ -553,6 +573,8 @@ obs_param_nml( 9)%bias_tcut = 432000 obs_param_nml( 9)%nodata = -9999. obs_param_nml( 9)%varname = 'sfmc' obs_param_nml( 9)%units = 'm3/m3' +obs_param_nml( 9)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml( 9)%fcstunits = 'NULL' ! do not edit! obs_param_nml( 9)%path = '/land/l_data/AMSR/data/AMSR_E_sm_LPRM/L2_EASE/bin/' obs_param_nml( 9)%name = 'AMSRsmUVA.EASE.v03.' obs_param_nml( 9)%maskpath = '' @@ -590,6 +612,8 @@ obs_param_nml(10)%bias_tcut = 432000 obs_param_nml(10)%nodata = -9999. obs_param_nml(10)%varname = 'sfmc' obs_param_nml(10)%units = 'm3/m3' +obs_param_nml(10)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(10)%fcstunits = 'NULL' ! do not edit! obs_param_nml(10)%path = '/land/l_data/AMSR/data/AMSR_E_sm_LPRM/L2_EASE/bin/' obs_param_nml(10)%name = 'AMSRsmUVA.EASE.v03.' obs_param_nml(10)%maskpath = '' @@ -627,6 +651,8 @@ obs_param_nml(11)%bias_tcut = 432000 obs_param_nml(11)%nodata = -9999. obs_param_nml(11)%varname = 'sfmc' obs_param_nml(11)%units = 'm3/m3' +obs_param_nml(11)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(11)%fcstunits = 'NULL' ! do not edit! obs_param_nml(11)%path = '/land/l_data/AMSR/data/AMSR_E_sm_LPRM/L2_EASE/bin/' obs_param_nml(11)%name = 'AMSRsmUVA.EASE.v03.' obs_param_nml(11)%maskpath = '' @@ -664,6 +690,8 @@ obs_param_nml(12)%bias_tcut = 432000 obs_param_nml(12)%nodata = -9999. obs_param_nml(12)%varname = 'sfmc' obs_param_nml(12)%units = 'm3/m3' +obs_param_nml(12)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(12)%fcstunits = 'NULL' ! do not edit! obs_param_nml(12)%path = '/land/l_data/AMSR/data/AMSR_E_sm_LPRM/L2_EASE/bin/' obs_param_nml(12)%name = 'AMSRsmUVA.EASE.v03.' obs_param_nml(12)%maskpath = '' @@ -705,6 +733,8 @@ obs_param_nml(13)%bias_tcut = 432000 obs_param_nml(13)%nodata = -9999. obs_param_nml(13)%varname = 'sfmc' obs_param_nml(13)%units = 'm3/m3' +obs_param_nml(13)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(13)%fcstunits = 'NULL' ! do not edit! obs_param_nml(13)%path = '/discover/nobackup/rreichle/l_data/ASCAT/TUW_W5.4/EASE/CONUS/bin/' obs_param_nml(13)%name = 'SDS_' obs_param_nml(13)%maskpath = '' @@ -746,6 +776,8 @@ obs_param_nml(14)%bias_tcut = 432000 obs_param_nml(14)%nodata = -9999. obs_param_nml(14)%varname = 'sfmc' obs_param_nml(14)%units = 'm3/m3' +obs_param_nml(14)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(14)%fcstunits = 'NULL' ! do not edit! obs_param_nml(14)%path = '/discover/nobackup/rreichle/l_data/ASCAT/TUW_W5.4/EASE/CONUS/bin/' obs_param_nml(14)%name = 'SDS_' obs_param_nml(14)%maskpath = '' @@ -783,6 +815,8 @@ obs_param_nml(15)%bias_tcut = 432000 obs_param_nml(15)%nodata = -9999. obs_param_nml(15)%varname = 'sfmc' obs_param_nml(15)%units = 'm3/m3' +obs_param_nml(15)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(15)%fcstunits = 'NULL' ! do not edit! obs_param_nml(15)%path = '/discover/nobackup/projects/gmao/ssd/land/l_data/SMOS/EASEv2/ESA_REPR/SMOS_M36_SMUDP2/' obs_param_nml(15)%name = '' obs_param_nml(15)%maskpath = '' @@ -820,6 +854,8 @@ obs_param_nml(16)%bias_tcut = 432000 obs_param_nml(16)%nodata = -9999. obs_param_nml(16)%varname = 'sfmc' obs_param_nml(16)%units = 'm3/m3' +obs_param_nml(16)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(16)%fcstunits = 'NULL' ! do not edit! obs_param_nml(16)%path = '/discover/nobackup/projects/gmao/ssd/land/l_data/SMOS/EASEv2/ESA_REPR/SMOS_M36_SMUDP2/' obs_param_nml(16)%name = '' obs_param_nml(16)%maskpath = '' @@ -905,6 +941,8 @@ obs_param_nml(17)%bias_tcut = 432000 obs_param_nml(17)%nodata = -9999. obs_param_nml(17)%varname = 'Tb' obs_param_nml(17)%units = 'K' +obs_param_nml(17)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(17)%fcstunits = 'NULL' ! do not edit! obs_param_nml(17)%path = '/discover/nobackup/projects/gmao/ssd/land/l_data/SMOS/EASEv2/ESA_REPR/SMOS_M36_SCLF1C_reg_nosky_noatm_v620_ESA_v102/' obs_param_nml(17)%name = '' obs_param_nml(17)%maskpath = '' @@ -949,6 +987,8 @@ obs_param_nml(18)%bias_tcut = 432000 obs_param_nml(18)%nodata = -9999. obs_param_nml(18)%varname = 'Tb' obs_param_nml(18)%units = 'K' +obs_param_nml(18)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(18)%fcstunits = 'NULL' ! do not edit! obs_param_nml(18)%path = '/discover/nobackup/projects/gmao/ssd/land/l_data/SMOS/EASEv2/ESA_REPR/SMOS_M36_SCLF1C_reg_nosky_noatm_v620_ESA_v102/' obs_param_nml(18)%name = '' obs_param_nml(18)%maskpath = '' @@ -993,6 +1033,8 @@ obs_param_nml(19)%bias_tcut = 432000 obs_param_nml(19)%nodata = -9999. obs_param_nml(19)%varname = 'Tb' obs_param_nml(19)%units = 'K' +obs_param_nml(19)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(19)%fcstunits = 'NULL' ! do not edit! obs_param_nml(19)%path = '/discover/nobackup/projects/gmao/ssd/land/l_data/SMOS/EASEv2/ESA_REPR/SMOS_M36_SCLF1C_reg_nosky_noatm_v620_ESA_v102/' obs_param_nml(19)%name = '' obs_param_nml(19)%maskpath = '' @@ -1037,6 +1079,8 @@ obs_param_nml(20)%bias_tcut = 432000 obs_param_nml(20)%nodata = -9999. obs_param_nml(20)%varname = 'Tb' obs_param_nml(20)%units = 'K' +obs_param_nml(20)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(20)%fcstunits = 'NULL' ! do not edit! obs_param_nml(20)%path = '/discover/nobackup/projects/gmao/ssd/land/l_data/SMOS/EASEv2/ESA_REPR/SMOS_M36_SCLF1C_reg_nosky_noatm_v620_ESA_v102/' obs_param_nml(20)%name = '' obs_param_nml(20)%maskpath = '' @@ -1075,6 +1119,8 @@ obs_param_nml(21)%bias_tcut = 432000 obs_param_nml(21)%nodata = -9999. obs_param_nml(21)%varname = 'Tb' obs_param_nml(21)%units = 'K' +obs_param_nml(21)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(21)%fcstunits = 'NULL' ! do not edit! obs_param_nml(21)%path = '/discover/nobackup/projects/gmao/smap/SMAP_Nature/SMOS/EASEv2/ESA/SMOS_M36_SCLF1C_fit_nosky_noatm_v724_ESA_v102/SMOS_fit_poly2/' obs_param_nml(21)%name = '' obs_param_nml(21)%maskpath = '' @@ -1113,6 +1159,8 @@ obs_param_nml(22)%bias_tcut = 432000 obs_param_nml(22)%nodata = -9999. obs_param_nml(22)%varname = 'Tb' obs_param_nml(22)%units = 'K' +obs_param_nml(22)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(22)%fcstunits = 'NULL' ! do not edit! obs_param_nml(22)%path = '/discover/nobackup/projects/gmao/smap/SMAP_Nature/SMOS/EASEv2/ESA/SMOS_M36_SCLF1C_fit_nosky_noatm_v724_ESA_v102/SMOS_fit_poly2/' obs_param_nml(22)%name = '' obs_param_nml(22)%maskpath = '' @@ -1151,6 +1199,8 @@ obs_param_nml(23)%bias_tcut = 432000 obs_param_nml(23)%nodata = -9999. obs_param_nml(23)%varname = 'Tb' obs_param_nml(23)%units = 'K' +obs_param_nml(23)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(23)%fcstunits = 'NULL' ! do not edit! obs_param_nml(23)%path = '/discover/nobackup/projects/gmao/smap/SMAP_Nature/SMOS/EASEv2/ESA/SMOS_M36_SCLF1C_fit_nosky_noatm_v724_ESA_v102/SMOS_fit_poly2/' obs_param_nml(23)%name = '' obs_param_nml(23)%maskpath = '' @@ -1189,6 +1239,8 @@ obs_param_nml(24)%bias_tcut = 432000 obs_param_nml(24)%nodata = -9999. obs_param_nml(24)%varname = 'Tb' obs_param_nml(24)%units = 'K' +obs_param_nml(24)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(24)%fcstunits = 'NULL' ! do not edit! obs_param_nml(24)%path = '/discover/nobackup/projects/gmao/smap/SMAP_Nature/SMOS/EASEv2/ESA/SMOS_M36_SCLF1C_fit_nosky_noatm_v724_ESA_v102/SMOS_fit_poly2/' obs_param_nml(24)%name = '' obs_param_nml(24)%maskpath = '' @@ -1227,6 +1279,8 @@ obs_param_nml(25)%bias_tcut = -9999 obs_param_nml(25)%nodata = -9999. obs_param_nml(25)%varname = 'NULL' obs_param_nml(25)%units = 'NULL' +obs_param_nml(25)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(25)%fcstunits = 'NULL' ! do not edit! obs_param_nml(25)%path = 'NULL' obs_param_nml(25)%name = 'NULL' obs_param_nml(25)%maskpath = 'NULL' @@ -1265,6 +1319,8 @@ obs_param_nml(26)%bias_tcut = 432000 obs_param_nml(26)%nodata = -9999. obs_param_nml(26)%varname = 'tsurf' obs_param_nml(26)%units = 'K' +obs_param_nml(26)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(26)%fcstunits = 'NULL' ! do not edit! obs_param_nml(26)%path = '/discover/nobackup/csdraper/LaRC_float/GOES-WEST/' obs_param_nml(26)%name = 'larc-v3.inst3_g15_Nch.' obs_param_nml(26)%maskpath = '' @@ -1303,6 +1359,8 @@ obs_param_nml(27)%bias_tcut = 432000 obs_param_nml(27)%nodata = -9999. obs_param_nml(27)%varname = 'tsurf' obs_param_nml(27)%units = 'K' +obs_param_nml(27)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(27)%fcstunits = 'NULL' ! do not edit! obs_param_nml(27)%path = '/discover/nobackup/csdraper/LaRC_float/v4/GOES-EAST/' obs_param_nml(27)%name = 'larc-v3.inst3_g13_Nch.' obs_param_nml(27)%maskpath = '' @@ -1341,6 +1399,8 @@ obs_param_nml(28)%bias_tcut = 432000 obs_param_nml(28)%nodata = -9999. obs_param_nml(28)%varname = 'tsurf' obs_param_nml(28)%units = 'K' +obs_param_nml(28)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(28)%fcstunits = 'NULL' ! do not edit! obs_param_nml(28)%path = '/discover/nobackup/csdraper/LaRC_float/MET09/' obs_param_nml(28)%name = 'larc-v3.inst3_mt9_Nch.' obs_param_nml(28)%maskpath = '' @@ -1379,6 +1439,8 @@ obs_param_nml(29)%bias_tcut = 432000 obs_param_nml(29)%nodata = -9999. obs_param_nml(29)%varname = 'tsurf' obs_param_nml(29)%units = 'K' +obs_param_nml(29)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(29)%fcstunits = 'NULL' ! do not edit! obs_param_nml(29)%path = '/discover/nobackup/csdraper/LaRC_float/FY2E/' obs_param_nml(29)%name = 'larc-v3.inst3_fye_Nch.' obs_param_nml(29)%maskpath = '' @@ -1417,6 +1479,8 @@ obs_param_nml(30)%bias_tcut = 432000 obs_param_nml(30)%nodata = -9999. obs_param_nml(30)%varname = 'tsurf' obs_param_nml(30)%units = 'K' +obs_param_nml(30)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(30)%fcstunits = 'NULL' ! do not edit! obs_param_nml(30)%path = '/discover/nobackup/csdraper/LaRC_float/MTSAT-2/' obs_param_nml(30)%name = 'larc-v3.inst3_mt2_Nch.' obs_param_nml(30)%maskpath = '' @@ -1465,6 +1529,8 @@ obs_param_nml(31)%bias_tcut = 432000 obs_param_nml(31)%nodata = -9999. obs_param_nml(31)%varname = 'Tb' obs_param_nml(31)%units = 'K' +obs_param_nml(31)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(31)%fcstunits = 'NULL' ! do not edit! obs_param_nml(31)%path = '/discover/nobackup/projects/gmao/smap/SMAP_L4/SMAP/OPS/L1C_TB/' obs_param_nml(31)%name = '' obs_param_nml(31)%maskpath = '' @@ -1503,6 +1569,8 @@ obs_param_nml(32)%bias_tcut = 432000 obs_param_nml(32)%nodata = -9999. obs_param_nml(32)%varname = 'Tb' obs_param_nml(32)%units = 'K' +obs_param_nml(32)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(32)%fcstunits = 'NULL' ! do not edit! obs_param_nml(32)%path = '/discover/nobackup/projects/gmao/smap/SMAP_L4/SMAP/OPS/L1C_TB/' obs_param_nml(32)%name = '' obs_param_nml(32)%maskpath = '' @@ -1541,6 +1609,8 @@ obs_param_nml(33)%bias_tcut = 432000 obs_param_nml(33)%nodata = -9999. obs_param_nml(33)%varname = 'Tb' obs_param_nml(33)%units = 'K' +obs_param_nml(33)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(33)%fcstunits = 'NULL' ! do not edit! obs_param_nml(33)%path = '/discover/nobackup/projects/gmao/smap/SMAP_L4/SMAP/OPS/L1C_TB/' obs_param_nml(33)%name = '' obs_param_nml(33)%maskpath = '' @@ -1579,6 +1649,8 @@ obs_param_nml(34)%bias_tcut = 432000 obs_param_nml(34)%nodata = -9999. obs_param_nml(34)%varname = 'Tb' obs_param_nml(34)%units = 'K' +obs_param_nml(34)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(34)%fcstunits = 'NULL' ! do not edit! obs_param_nml(34)%path = '/discover/nobackup/projects/gmao/smap/SMAP_L4/SMAP/OPS/L1C_TB/' obs_param_nml(34)%name = '' obs_param_nml(34)%maskpath = '' @@ -1628,6 +1700,8 @@ obs_param_nml(35)%bias_tcut = 432000 obs_param_nml(35)%nodata = -9999. obs_param_nml(35)%varname = 'Tb' obs_param_nml(35)%units = 'K' +obs_param_nml(35)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(35)%fcstunits = 'NULL' ! do not edit! obs_param_nml(35)%path = '/discover/nobackup/projects/gmao/smap/SMAP_L4/SMAP/OPS/L2_SM_AP/' obs_param_nml(35)%name = '' obs_param_nml(35)%maskpath = '' @@ -1666,6 +1740,8 @@ obs_param_nml(36)%bias_tcut = 432000 obs_param_nml(36)%nodata = -9999. obs_param_nml(36)%varname = 'Tb' obs_param_nml(36)%units = 'K' +obs_param_nml(36)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(36)%fcstunits = 'NULL' ! do not edit! obs_param_nml(36)%path = '/discover/nobackup/projects/gmao/smap/SMAP_L4/SMAP/OPS/L2_SM_AP/' obs_param_nml(36)%name = '' obs_param_nml(36)%maskpath = '' @@ -1704,6 +1780,8 @@ obs_param_nml(37)%bias_tcut = 432000 obs_param_nml(37)%nodata = -9999. obs_param_nml(37)%varname = 'Tb' obs_param_nml(37)%units = 'K' +obs_param_nml(37)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(37)%fcstunits = 'NULL' ! do not edit! obs_param_nml(37)%path = '/discover/nobackup/projects/gmao/smap/SMAP_L4/SMAP/OPS/L2_SM_AP/' obs_param_nml(37)%name = '' obs_param_nml(37)%maskpath = '' @@ -1742,6 +1820,8 @@ obs_param_nml(38)%bias_tcut = 432000 obs_param_nml(38)%nodata = -9999. obs_param_nml(38)%varname = 'Tb' obs_param_nml(38)%units = 'K' +obs_param_nml(38)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(38)%fcstunits = 'NULL' ! do not edit! obs_param_nml(38)%path = '/discover/nobackup/projects/gmao/smap/SMAP_L4/SMAP/OPS/L2_SM_AP/' obs_param_nml(38)%name = '' obs_param_nml(38)%maskpath = '' @@ -1779,6 +1859,8 @@ obs_param_nml(39)%bias_tcut = 432000 obs_param_nml(39)%nodata = -9999. obs_param_nml(39)%varname = 'FT' obs_param_nml(39)%units = '-' +obs_param_nml(39)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(39)%fcstunits = 'NULL' ! do not edit! obs_param_nml(39)%path = '/discover/nobackup/projects/gmao/smap/SMAP_L4/SMAP/OPS/L2_SM_AP/' obs_param_nml(39)%name = '' obs_param_nml(39)%maskpath = '' @@ -1816,6 +1898,8 @@ obs_param_nml(40)%bias_tcut = 432000 obs_param_nml(40)%nodata = -9999. obs_param_nml(40)%varname = 'FT' obs_param_nml(40)%units = '-' +obs_param_nml(40)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(40)%fcstunits = 'NULL' ! do not edit! obs_param_nml(40)%path = '/discover/nobackup/projects/gmao/smap/SMAP_L4/SMAP/OPS/L2_SM_AP/' obs_param_nml(40)%name = '' obs_param_nml(40)%maskpath = '' @@ -1876,6 +1960,8 @@ obs_param_nml(41)%bias_tcut = 432000 obs_param_nml(41)%nodata = -9999. obs_param_nml(41)%varname = 'Tb' obs_param_nml(41)%units = 'K' +obs_param_nml(41)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(41)%fcstunits = 'NULL' ! do not edit! obs_param_nml(41)%path = '/discover/nobackup/projects/gmao/smap/SMAP_L4/SMAP/OPS/L1C_TB_E/' obs_param_nml(41)%name = '' obs_param_nml(41)%maskpath = '' @@ -1914,6 +2000,8 @@ obs_param_nml(42)%bias_tcut = 432000 obs_param_nml(42)%nodata = -9999. obs_param_nml(42)%varname = 'Tb' obs_param_nml(42)%units = 'K' +obs_param_nml(42)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(42)%fcstunits = 'NULL' ! do not edit! obs_param_nml(42)%path = '/discover/nobackup/projects/gmao/smap/SMAP_L4/SMAP/OPS/L1C_TB_E/' obs_param_nml(42)%name = '' obs_param_nml(42)%maskpath = '' @@ -1952,6 +2040,8 @@ obs_param_nml(43)%bias_tcut = 432000 obs_param_nml(43)%nodata = -9999. obs_param_nml(43)%varname = 'Tb' obs_param_nml(43)%units = 'K' +obs_param_nml(43)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(43)%fcstunits = 'NULL' ! do not edit! obs_param_nml(43)%path = '/discover/nobackup/projects/gmao/smap/SMAP_L4/SMAP/OPS/L1C_TB_E/' obs_param_nml(43)%name = '' obs_param_nml(43)%maskpath = '' @@ -1990,6 +2080,8 @@ obs_param_nml(44)%bias_tcut = 432000 obs_param_nml(44)%nodata = -9999. obs_param_nml(44)%varname = 'Tb' obs_param_nml(44)%units = 'K' +obs_param_nml(44)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(44)%fcstunits = 'NULL' ! do not edit! obs_param_nml(44)%path = '/discover/nobackup/projects/gmao/smap/SMAP_L4/SMAP/OPS/L1C_TB_E/' obs_param_nml(44)%name = '' obs_param_nml(44)%maskpath = '' @@ -2028,6 +2120,8 @@ obs_param_nml(45)%bias_tcut = 432000 obs_param_nml(45)%nodata = -9999. obs_param_nml(45)%varname = 'Tb' obs_param_nml(45)%units = 'K' +obs_param_nml(45)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(45)%fcstunits = 'NULL' ! do not edit! obs_param_nml(45)%path = '/discover/nobackup/projects/gmao/smap/SMAP_L4/SMAP/OPS/L1C_TB_E/' obs_param_nml(45)%name = '' obs_param_nml(45)%maskpath = '' @@ -2066,6 +2160,8 @@ obs_param_nml(46)%bias_tcut = 432000 obs_param_nml(46)%nodata = -9999. obs_param_nml(46)%varname = 'Tb' obs_param_nml(46)%units = 'K' +obs_param_nml(46)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(46)%fcstunits = 'NULL' ! do not edit! obs_param_nml(46)%path = '/discover/nobackup/projects/gmao/smap/SMAP_L4/SMAP/OPS/L1C_TB_E/' obs_param_nml(46)%name = '' obs_param_nml(46)%maskpath = '' @@ -2104,6 +2200,8 @@ obs_param_nml(47)%bias_tcut = 432000 obs_param_nml(47)%nodata = -9999. obs_param_nml(47)%varname = 'Tb' obs_param_nml(47)%units = 'K' +obs_param_nml(47)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(47)%fcstunits = 'NULL' ! do not edit! obs_param_nml(47)%path = '/discover/nobackup/projects/gmao/smap/SMAP_L4/SMAP/OPS/L1C_TB_E/' obs_param_nml(47)%name = '' obs_param_nml(47)%maskpath = '' @@ -2142,6 +2240,8 @@ obs_param_nml(48)%bias_tcut = 432000 obs_param_nml(48)%nodata = -9999. obs_param_nml(48)%varname = 'Tb' obs_param_nml(48)%units = 'K' +obs_param_nml(48)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(48)%fcstunits = 'NULL' ! do not edit! obs_param_nml(48)%path = '/discover/nobackup/projects/gmao/smap/SMAP_L4/SMAP/OPS/L1C_TB_E/' obs_param_nml(48)%name = '' obs_param_nml(48)%maskpath = '' @@ -2187,6 +2287,8 @@ obs_param_nml(49)%bias_tcut = 432000 obs_param_nml(49)%nodata = -9999. obs_param_nml(49)%varname = 'sfds' obs_param_nml(49)%units = '%' +obs_param_nml(49)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(49)%fcstunits = 'NULL' ! do not edit! obs_param_nml(49)%path = '/discover/nobackup/projects/gmao/smap/SMAP_Nature/ASCAT_EUMETSAT/Metop_A/' obs_param_nml(49)%name = '' obs_param_nml(49)%maskpath = '' @@ -2226,6 +2328,8 @@ obs_param_nml(50)%bias_tcut = 432000 obs_param_nml(50)%nodata = -9999. obs_param_nml(50)%varname = 'sfds' obs_param_nml(50)%units = '%' +obs_param_nml(50)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(50)%fcstunits = 'NULL' ! do not edit! obs_param_nml(50)%path = '/discover/nobackup/projects/gmao/smap/SMAP_Nature/ASCAT_EUMETSAT/Metop_B/' obs_param_nml(50)%name = '' obs_param_nml(50)%maskpath = '' @@ -2265,6 +2369,8 @@ obs_param_nml(51)%bias_tcut = 432000 obs_param_nml(51)%nodata = -9999. obs_param_nml(51)%varname = 'sfds' obs_param_nml(51)%units = '%' +obs_param_nml(51)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(51)%fcstunits = 'NULL' ! do not edit! obs_param_nml(51)%path = '/discover/nobackup/projects/gmao/smap/SMAP_Nature/ASCAT_EUMETSAT/Metop_C/' obs_param_nml(51)%name = '' obs_param_nml(51)%maskpath = '' @@ -2306,6 +2412,8 @@ obs_param_nml(52)%bias_tcut = 432000 obs_param_nml(52)%nodata = -9999. obs_param_nml(52)%varname = 'asnow' obs_param_nml(52)%units = 'm2/m2' +obs_param_nml(52)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(52)%fcstunits = 'NULL' ! do not edit! obs_param_nml(52)%path = '/discover/nobackup/projects/S2SHMA/MODIS/MYD10C1_V61/' obs_param_nml(52)%name = 'MYD10C1.Ayyyyddd.061.hdf' obs_param_nml(52)%maskpath = '' @@ -2347,6 +2455,8 @@ obs_param_nml(53)%bias_tcut = 432000 obs_param_nml(53)%nodata = -9999. obs_param_nml(53)%varname = 'asnow' obs_param_nml(53)%units = 'm2/m2' +obs_param_nml(53)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(53)%fcstunits = 'NULL' ! do not edit! obs_param_nml(53)%path = '/discover/nobackup/projects/S2SHMA/MODIS/MOD10C1_V61/' obs_param_nml(53)%name = 'MOD10C1.Ayyyyddd.061.hdf' obs_param_nml(53)%maskpath = '' @@ -2386,6 +2496,8 @@ obs_param_nml(54)%bias_tcut = 432000 obs_param_nml(54)%nodata = -9999. obs_param_nml(54)%varname = 'sfmc' obs_param_nml(54)%units = 'm3 m-3' +obs_param_nml(54)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(54)%fcstunits = 'NULL' ! do not edit! obs_param_nml(54)%path = '/discover/nobackup/projects/gmao/smap/SMAP_Nature/CYGNSS/' obs_param_nml(54)%name = '' obs_param_nml(54)%maskpath = '' @@ -2425,6 +2537,8 @@ obs_param_nml(55)%bias_tcut = 432000 obs_param_nml(55)%nodata = -9999. obs_param_nml(55)%varname = 'sfmc' obs_param_nml(55)%units = 'm3 m-3' +obs_param_nml(55)%fcstvarname = 'NULL' ! do not edit! +obs_param_nml(55)%fcstunits = 'NULL' ! do not edit! obs_param_nml(55)%path = '/discover/nobackup/projects/gmao/smap/SMAP_Nature/CYGNSS/' obs_param_nml(55)%name = '' obs_param_nml(55)%maskpath = '' diff --git a/GEOSldas_App/LandAtmDAS_nml/LDASsa_SPECIAL_inputs_ensupd.nml b/GEOSldas_App/LandAtmDAS_nml/LDASsa_SPECIAL_inputs_ensupd.nml index e8ecb12e..e8bc5978 100644 --- a/GEOSldas_App/LandAtmDAS_nml/LDASsa_SPECIAL_inputs_ensupd.nml +++ b/GEOSldas_App/LandAtmDAS_nml/LDASsa_SPECIAL_inputs_ensupd.nml @@ -37,7 +37,7 @@ update_type = 10 out_obslog = .true. -out_ObsFcstAna = .true. +out_ObsFcstAna = 1 ! 0=OFF, 1=BIN (legacy unformatted), 2=NC4, 3=BIN and NC4 out_smapL4SMaup = .false. ! --------------------------------------------------------------------- diff --git a/GEOSldas_App/lenkf_j_template.py b/GEOSldas_App/lenkf_j_template.py index 60443ff0..53bfe5fb 100644 --- a/GEOSldas_App/lenkf_j_template.py +++ b/GEOSldas_App/lenkf_j_template.py @@ -471,7 +471,7 @@ # must be done before moving HISTORY files - set ObsFcses = `ls *.ldas_ObsFcstAna.*.bin` + set ObsFcses = `ls *.ldas_ObsFcstAna.*` foreach obsfcs ( $ObsFcses ) set ThisTime = `echo $obsfcs | rev | cut -d'.' -f2 | rev` set TY = `echo $ThisTime | cut -c1-4` diff --git a/GEOSldas_App/util/shared/matlab/read_obsparam.m b/GEOSldas_App/util/shared/matlab/read_obsparam.m index 19e66056..433a1f07 100644 --- a/GEOSldas_App/util/shared/matlab/read_obsparam.m +++ b/GEOSldas_App/util/shared/matlab/read_obsparam.m @@ -7,7 +7,9 @@ % % 1 Dec 2011 - reichle: minor modifications and check-in to CVS % -% 8 Jun 2017 - reichle: added "flistpath" and "flistname" +% 8 Jun 2017 - reichle: added "flistpath" and "flistname" +% +% 17 Apr 2026 - reichle: added "fcstvarname" and "fcstunits" % % ------------------------------------------------------------------ @@ -40,6 +42,8 @@ obs_param(i).nodata = fscanf(fid, '%f ', 1); obs_param(i).varname = fscanf(fid, '%s ', 1); obs_param(i).units = fscanf(fid, '%s ', 1); + obs_param(i).fcstvarname = fscanf(fid, '%s ', 1); + obs_param(i).fcstunits = fscanf(fid, '%s ', 1); obs_param(i).path = fscanf(fid, '%s ', 1); obs_param(i).name = fscanf(fid, '%s ', 1); obs_param(i).maskpath = fscanf(fid, '%s ', 1); @@ -58,18 +62,20 @@ % remove leading and trailing quotes from strings - obs_param(i).descr = obs_param(i).descr( 2:end-1); - obs_param(i).FOV_units = obs_param(i).FOV_units(2:end-1); - obs_param(i).varname = obs_param(i).varname( 2:end-1); - obs_param(i).units = obs_param(i).units( 2:end-1); - obs_param(i).path = obs_param(i).path( 2:end-1); - obs_param(i).name = obs_param(i).name( 2:end-1); - obs_param(i).maskpath = obs_param(i).maskpath( 2:end-1); - obs_param(i).maskname = obs_param(i).maskname( 2:end-1); - obs_param(i).scalepath = obs_param(i).scalepath(2:end-1); - obs_param(i).scalename = obs_param(i).scalename(2:end-1); - obs_param(i).flistpath = obs_param(i).flistpath(2:end-1); - obs_param(i).flistname = obs_param(i).flistname(2:end-1); + obs_param(i).descr = obs_param(i).descr( 2:end-1); + obs_param(i).FOV_units = obs_param(i).FOV_units( 2:end-1); + obs_param(i).varname = obs_param(i).varname( 2:end-1); + obs_param(i).units = obs_param(i).units( 2:end-1); + obs_param(i).fcstvarname = obs_param(i).fcstvarname(2:end-1); + obs_param(i).fcstunits = obs_param(i).fcstunits( 2:end-1); + obs_param(i).path = obs_param(i).path( 2:end-1); + obs_param(i).name = obs_param(i).name( 2:end-1); + obs_param(i).maskpath = obs_param(i).maskpath( 2:end-1); + obs_param(i).maskname = obs_param(i).maskname( 2:end-1); + obs_param(i).scalepath = obs_param(i).scalepath( 2:end-1); + obs_param(i).scalename = obs_param(i).scalename( 2:end-1); + obs_param(i).flistpath = obs_param(i).flistpath( 2:end-1); + obs_param(i).flistname = obs_param(i).flistname( 2:end-1); end diff --git a/LDAS_Shared/LDAS_ensdrv_mpi.F90 b/LDAS_Shared/LDAS_ensdrv_mpi.F90 index 1da40cdf..6b36a77a 100644 --- a/LDAS_Shared/LDAS_ensdrv_mpi.F90 +++ b/LDAS_Shared/LDAS_ensdrv_mpi.F90 @@ -739,6 +739,8 @@ subroutine init_MPI_types() ! real :: nodata ! block #7 (real) ! character(40) :: varname ! block #8 (character) ! character(40) :: units + ! character(40) :: fcstvarname + ! character(40) :: fcstunits ! character(200) :: path ! character(80) :: name ! character(200) :: maskpath @@ -781,7 +783,7 @@ subroutine init_MPI_types() iblock( 5) = 3 iblock( 6) = 4 iblock( 7) = 1 - iblock( 8) = 40+40+200+80+200+80+200+80+200+80 + iblock( 8) = 40+40+40+40+200+80+200+80+200+80+200+80 iblock( 9) = 2 iblock(10) = 2 iblock(11) = 2 diff --git a/LDAS_Shared/enkf_types.F90 b/LDAS_Shared/enkf_types.F90 index 566c75ef..cb1daf06 100644 --- a/LDAS_Shared/enkf_types.F90 +++ b/LDAS_Shared/enkf_types.F90 @@ -56,7 +56,8 @@ module enkf_types ! added "varname" field to "obs_param_type" - reichle 14 Jun 2011 ! major revisions to "obs_type" fields - reichle 16 Jun 2011 ! added "units" field to "obs_param_type" - reichle 22 Nov 2011 - + ! added "fcstvarname" and "fcstunits" to "obs_param_type" - reichle 17 Apr 2026 + type :: obs_type ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @@ -73,8 +74,8 @@ module enkf_types real*8 :: time ! time of obs (J2000 seconds w/ 'TT12' epoch; see date_time_util.F90) real :: lon ! longitude of obs real :: lat ! latitude of obs - real :: obs ! observed value - real :: obsvar ! obs error var + real :: obs ! observed value (after scaling) + real :: obsvar ! obs error var (after scaling) real :: fcst ! "forecast": value of obs pred before EnKF update (ens mean) real :: fcstvar ! forecast error var (in obs space), a.k.a. HPHt real :: ana ! "analysis": value of obs pred after EnKF update (ens mean) @@ -94,68 +95,74 @@ module enkf_types ! any subroutines or operators defined herein ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - character(40) :: descr ! description - integer :: species ! identifier for type of measurement + character(40) :: descr ! description + integer :: species ! identifier for type of observation - integer :: orbit ! type of (half-)orbit - ! 0 = n/a [eg., in situ obs] - ! 1 = ascending - ! 2 = descending - ! 3 = ascending or descending - ! 4 = geostationary + integer :: orbit ! type of (half-)orbit + ! 0 = n/a [eg., in situ obs] + ! 1 = ascending + ! 2 = descending + ! 3 = ascending or descending + ! 4 = geostationary - integer :: pol ! polarization - ! 0 = n/a [eg., multi-pol. retrieval] - ! 1 = horizontal - ! 2 = vertical - ! 3 = ... - ! [add 3rd/4th Stokes, HH, HV, VH, VV] + integer :: pol ! polarization + ! 0 = n/a [eg., multi-pol. retrieval] + ! 1 = horizontal + ! 2 = vertical + ! 3 = ... [cpuld add 3rd/4th Stokes, HH, HV, VH, VV] - integer :: N_ang ! # satellite viewing angles in species (radiance obs only) + integer :: N_ang ! # satellite viewing angles in species (radiance obs only) real, & - dimension(N_obs_ang_max) :: ang ! vector of satellite viewing angles + dimension(N_obs_ang_max) :: ang ! vector of satellite viewing angles - real :: freq ! frequency [Hz] - - real :: FOV ! field-of-view *radius* - ! if FOV==0. equate obs footprint w/ tile - ! for details see LDASsa_DEFAULT_inputs ensupd.nml - character(40) :: FOV_units ! FOV units ('km' or 'deg') - - logical :: assim ! assimilate yes/no? (see also "obs_type") - logical :: scale ! scale yes/no? - logical :: getinnov ! compute innovs? (.T. if assim==.T.) - - integer :: RTM_ID ! ID of radiative transfer model - - integer :: bias_Npar ! number of bias states tracked per day - integer :: bias_trel ! e-folding time scale of obs bias memory [s] - integer :: bias_tcut ! cutoff time for confident obs bias est [s] - - real :: nodata ! no-data-value - - character(40) :: varname ! equivalent model variable name (Obs_pred) - character(40) :: units ! units (eg., 'K' or 'm3/m3') - - character(200) :: path ! path to measurements file - character(80) :: name ! name identifier for measurements - character(200) :: maskpath ! path to obs mask file - character(80) :: maskname ! filename for obs mask - character(200) :: scalepath ! path to file with scaling parameters - character(80) :: scalename ! filename for scaling parameters - character(200) :: flistpath ! path to file with list of obs file names - character(80) :: flistname ! name of file with list of obs file names - - real :: errstd ! default obs error std + real :: freq ! frequency [Hz] + + real :: FOV ! field-of-view *radius* + ! if FOV==0. equate obs footprint w/ tile + ! for details see LDASsa_DEFAULT_inputs ensupd.nml + character(40) :: FOV_units ! FOV units ('km' or 'deg') + + logical :: assim ! assimilate obs species: yes/no? (see also "obs_type") + logical :: scale ! scale obs species: yes/no? + logical :: getinnov ! compute innovs? (set .T. if assim==.T.) + + integer :: RTM_ID ! ID of radiative transfer model to use for Tb forward modeling (see get_obs_pred()) + ! 0 = none + ! 1 = L-band tau-omega model as in De Lannoy et al. 2013 (doi:10.1175/JHM-D-12-092.1) (old SMOS preproc) + ! 2 = same as 1 but without Pellarin atm corr (SMAP) + ! 3 = same as 1 but with Mironov and SMAP L2_SM pol mixing + ! 4 = same as 3 but without Pellarin atm corr (SMAP L4_SM Version 8) + + integer :: bias_Npar ! number of obs bias states tracked per day + integer :: bias_trel ! e-folding time scale of obs bias memory [s] + integer :: bias_tcut ! cutoff time for confident obs bias estimate [s] + + real :: nodata ! no-data-value + + character(40) :: varname ! observation native variable name (before scaling) + character(40) :: units ! observation native units (before scaling; eg., 'K' or 'm3/m3') + character(40) :: fcstvarname ! observation-equivalent model variable name (Obs_pred) [do not edit, filled by read_ens_upd_inputs()] + character(40) :: fcstunits ! observation-equivalent model units (eg., 'K' or 'm3/m3') [do not edit, filled by read_ens_upd_inputs()] + + character(200) :: path ! path to observations file + character(80) :: name ! name identifier for file containing observations + character(200) :: maskpath ! path to obs mask file + character(80) :: maskname ! filename for obs mask + character(200) :: scalepath ! path to file(s) with scaling parameters + character(80) :: scalename ! filename for scaling parameters + character(200) :: flistpath ! path to file with list of obs file names + character(80) :: flistname ! name of file with list of obs file names + + real :: errstd ! default obs error std (before scaling) - real :: std_normal_max ! see pert_param_type - logical :: zeromean ! see pert_param_type - logical :: coarsen_pert ! see pert_param_type ("%coarsen") - real :: xcorr ! see pert_param_type - real :: ycorr ! see pert_param_type + real :: std_normal_max ! maximum allowed obs perturbation (relative to N(0,1)) (see pert_param_type) + logical :: zeromean ! enforce zero mean across ensemble (see pert_param_type) + logical :: coarsen_pert ! generate obs perturbations on coarser grid (see pert_param_type) + real :: xcorr ! obs error correlation length (deg) in longitude direction (see pert_param_type) + real :: ycorr ! obs error correlation length (deg) in latitude direction (see pert_param_type) - integer :: adapt ! identifier for adaptive filtering + integer :: adapt ! identifier for adaptive filtering end type obs_param_type @@ -201,16 +208,18 @@ subroutine write_obs_param(unitnumber, N_obs_param, obs_param) write (unitnumber, *) obs_param(i)%bias_trel write (unitnumber, *) obs_param(i)%bias_tcut write (unitnumber, *) obs_param(i)%nodata - write (unitnumber, '(42A)') "'" // trim(obs_param(i)%varname) // "'" - write (unitnumber, '(42A)') "'" // trim(obs_param(i)%units) // "'" - write (unitnumber,'(202A)') "'" // trim(obs_param(i)%path) // "'" - write (unitnumber, '(82A)') "'" // trim(obs_param(i)%name) // "'" - write (unitnumber,'(202A)') "'" // trim(obs_param(i)%maskpath) // "'" - write (unitnumber, '(82A)') "'" // trim(obs_param(i)%maskname) // "'" - write (unitnumber,'(202A)') "'" // trim(obs_param(i)%scalepath) // "'" - write (unitnumber, '(82A)') "'" // trim(obs_param(i)%scalename) // "'" - write (unitnumber,'(202A)') "'" // trim(obs_param(i)%flistpath) // "'" - write (unitnumber, '(82A)') "'" // trim(obs_param(i)%flistname) // "'" + write (unitnumber, '(42A)') "'" // trim(obs_param(i)%varname) // "'" + write (unitnumber, '(42A)') "'" // trim(obs_param(i)%units) // "'" + write (unitnumber, '(42A)') "'" // trim(obs_param(i)%fcstvarname) // "'" + write (unitnumber, '(42A)') "'" // trim(obs_param(i)%fcstunits) // "'" + write (unitnumber,'(202A)') "'" // trim(obs_param(i)%path) // "'" + write (unitnumber, '(82A)') "'" // trim(obs_param(i)%name) // "'" + write (unitnumber,'(202A)') "'" // trim(obs_param(i)%maskpath) // "'" + write (unitnumber, '(82A)') "'" // trim(obs_param(i)%maskname) // "'" + write (unitnumber,'(202A)') "'" // trim(obs_param(i)%scalepath) // "'" + write (unitnumber, '(82A)') "'" // trim(obs_param(i)%scalename) // "'" + write (unitnumber,'(202A)') "'" // trim(obs_param(i)%flistpath) // "'" + write (unitnumber, '(82A)') "'" // trim(obs_param(i)%flistname) // "'" write (unitnumber, *) obs_param(i)%errstd write (unitnumber, *) obs_param(i)%std_normal_max write (unitnumber, *) obs_param(i)%zeromean