diff --git a/src/soca/Fields/soca_fields_mod.F90 b/src/soca/Fields/soca_fields_mod.F90 index eb5cf225d..54fde11e0 100644 --- a/src/soca/Fields/soca_fields_mod.F90 +++ b/src/soca/Fields/soca_fields_mod.F90 @@ -27,9 +27,6 @@ module soca_fields_mod use fms_io_mod, only: fms_io_init, fms_io_exit, register_restart_field, & restart_file_type, restore_state, free_restart_type, save_restart use fms_mod, only: write_data, set_domain -use horiz_interp_mod, only : horiz_interp_type -use horiz_interp_spherical_mod, only : horiz_interp_spherical, horiz_interp_spherical_del, & - horiz_interp_spherical_new use MOM_remapping, only : remapping_CS, initialize_remapping, remapping_core_h, & end_remapping use mpp_domains_mod, only : mpp_update_domains, mpp_update_domains_ad @@ -38,6 +35,7 @@ module soca_fields_mod use soca_fields_metadata_mod, only : soca_field_metadata use soca_geom_mod, only : soca_geom use soca_utils, only: soca_mld +use soca_utils, only: soca_stencil_interp, soca_stencil_neighbors implicit none private @@ -199,8 +197,8 @@ module soca_fields_mod !> \copybrief soca_fields_update_halos \see soca_fields_update_halos procedure :: update_halos => soca_fields_update_halos - !> \copybrief soca_fields_colocate \see soca_fields_colocate - procedure :: colocate => soca_fields_colocate + !> \copybrief soca_fields_tohpoints \see soca_fields_tohpoints + procedure :: tohpoints => soca_fields_tohpoints !> \} !> \name serialization @@ -278,36 +276,6 @@ subroutine soca_field_update_halo(self, geom) call mpp_update_domains(self%val, geom%Domain%mpp_domain) end subroutine soca_field_update_halo - -! ------------------------------------------------------------------------------ -!> Perform spatial interpolation between two grids. -!! -!! Interpolation used is inverse distance weidghted, taking into -!! consideration the mask. -!! \param[in] geom: The geometry to interpolate to -!! \param[in] interp2d: interpolation object created by calling -!! \c horiz_interp_spherical_new() in FMS -!! -!! \relates soca_fields_mod::soca_field -subroutine soca_field_stencil_interp(self, geom, interp2d) - class(soca_field), intent(inout) :: self - type(soca_geom), pointer, intent(in) :: geom - type(horiz_interp_type), intent(in) :: interp2d - - integer :: k - real(kind=kind_real), allocatable :: val(:,:,:) - - allocate(val, mold=self%val) - val = self%val - do k = 1, self%nz - call horiz_interp_spherical(interp2d, & - & val(geom%isd:geom%ied, geom%jsd:geom%jed,k), & - & self%val(geom%isc:geom%iec, geom%jsc:geom%jec,k)) - end do - call self%update_halo(geom) -end subroutine soca_field_stencil_interp - - ! ------------------------------------------------------------------------------ !> Make sure the two fields are the same in terms of name, size, shape. !! @@ -328,6 +296,93 @@ subroutine soca_field_check_congruent(self, rhs) end do end subroutine soca_field_check_congruent +! ------------------------------------------------------------------------------ +!> Perform spatial interpolation between adjacent grid point in the same stencil +!! +!! Interpolation used is inverse distance weidghted, taking into +!! consideration the mask and using at most 6 neighbors. +subroutine soca_field_stencil_interp(self, geom, fromto) + class(soca_field), intent(inout) :: self + class(soca_geom), intent(in) :: geom !< geometry + character(len=4), intent(in) :: fromto !< "u2h", "v2h" + + integer :: i, j + real(kind=kind_real), allocatable :: val_tmp(:,:,:) + real(kind=kind_real) :: val_max = 9e8_kind_real + integer :: ij(2,6), sti, nn + real(kind_real) :: lon_src(6), lat_src(6) + real(kind=kind_real), allocatable :: val(:,:) + real(kind=kind_real), allocatable :: lonsrc_local(:,:), latsrc_local(:,:) + real(kind=kind_real), allocatable :: londst_local(:,:), latdst_local(:,:) + real(kind=kind_real), allocatable :: masksrc_local(:,:), maskdst_local(:,:) + + ! Initialize temporary arrays + allocate(val_tmp, mold=self%val) + val_tmp = 0_kind_real + + ! Identify source and destination grids + select case(fromto) + case("vtoh") + ! Horizontal interpolation: v-points to h-points + allocate(lonsrc_local, mold=geom%lonv); lonsrc_local = geom%lonv + allocate(latsrc_local, mold=geom%latv); latsrc_local = geom%latv + allocate(masksrc_local, mold=geom%mask2dv); masksrc_local = geom%mask2dv + allocate(londst_local, mold=geom%lon); londst_local = geom%lon + allocate(latdst_local, mold=geom%lat); latdst_local = geom%lat + allocate(maskdst_local, mold=geom%mask2d); maskdst_local = geom%mask2d + + case("utoh") + ! Horizontal interpolation: u-points to h-points + allocate(lonsrc_local, mold=geom%lonu); lonsrc_local = geom%lonu + allocate(latsrc_local, mold=geom%latu); latsrc_local = geom%latu + allocate(masksrc_local, mold=geom%mask2du); masksrc_local = geom%mask2du + allocate(londst_local, mold=geom%lon); londst_local = geom%lon + allocate(latdst_local, mold=geom%lat); latdst_local = geom%lat + allocate(maskdst_local, mold=geom%mask2d); maskdst_local = geom%mask2d + + case default + call abor1_ftn('soca_field::stencil_interp, option '//fromto//& + ' not implemented yet') + + end select + + ! Interpolate + allocate(val(6,self%nz)) + do j = geom%jsc, geom%jec + do i = geom%isc, geom%iec + ! destination on land, skip + if (maskdst_local(i,j) == 0_kind_real) cycle + + ! get the 6 or less src-point neighbors surrounding the (i,j) dst-point + call soca_stencil_neighbors(fromto, i, j, ij) + nn = 1 + val = 0_kind_real + do sti = 1, 6 + ! source point on land, skip + if (masksrc_local(ij(1,sti), ij(2,sti)) == 0_kind_real) cycle + + ! outcroping of layers, skip + if (abs(self%val(ij(1,sti), ij(2,sti),1)) > val_max) cycle + + ! store the valid neighbors + lon_src(nn) = lonsrc_local(ij(1,sti), ij(2,sti)) + lat_src(nn) = latsrc_local(ij(1,sti), ij(2,sti)) + val(nn,:) = self%val(ij(1,sti), ij(2,sti),:) + nn = nn + 1 + end do + nn = nn - 1 + + ! val_tmp: interpolated val at (i,j) dst-point along layers + if ( nn >=1 ) then + call soca_stencil_interp(lon_src, lat_src, & + londst_local(i,j), latdst_local(i,j), & + val, val_tmp(i,j,:), nn) + end if + end do + end do + self%val = val_tmp + +end subroutine soca_field_stencil_interp ! ------------------------------------------------------------------------------ !> Delete the soca_field object. @@ -1205,81 +1260,39 @@ subroutine soca_fields_write_rst(self, f_conf, vdate) end subroutine soca_fields_write_rst ! ------------------------------------------------------------------------------ -!> Colocate by interpolating from one c-grid location to another. +!> Interpolates from uv-points location to h-points. !! -!! \warning only works on the "h" grid currently (not the "u" or "v" grid) !! \relates soca_fields_mod::soca_fields -subroutine soca_fields_colocate(self, cgridlocout) - class(soca_fields), intent(inout) :: self !< self - character(len=1), intent(in) :: cgridlocout !< colocate to cgridloc (u, v or h) +subroutine soca_fields_tohpoints(self) + class(soca_fields), intent(inout) :: self !< self - integer :: i, k + integer :: i real(kind=kind_real), allocatable :: val(:,:,:) real(kind=kind_real), pointer :: lon_out(:,:) => null() real(kind=kind_real), pointer :: lat_out(:,:) => null() - type(soca_geom), pointer :: g => null() - type(horiz_interp_type) :: interp2d - - ! Associate lon_out and lat_out according to cgridlocout - select case(cgridlocout) - ! TODO: Test colocation to u and v grid - !case ('u') - ! lon_out => self%geom%lonu - ! lat_out => self%geom%latu - !case ('v') - ! lon_out => self%geom%lonv - ! lat_out => self%geom%latv - case ('h') - lon_out => self%geom%lon - lat_out => self%geom%lat - case default - call abor1_ftn('soca_fields::colocate(): unknown c-grid location '// cgridlocout) - end select + character(len=4) :: fromto + + ! Associate lon_out and lat_out with the h-grid + lon_out => self%geom%lon + lat_out => self%geom%lat ! Apply interpolation to all fields, when necessary do i=1,size(self%fields) + ! Check if already on h-points + if (self%fields(i)%metadata%grid == 'h') cycle - ! Check if already colocated - if (self%fields(i)%metadata%grid == cgridlocout) cycle - - ! Initialize fms spherical idw interpolation - g => self%geom - call horiz_interp_spherical_new(interp2d, & - & real(deg2rad*self%fields(i)%lon(g%isd:g%ied,g%jsd:g%jed), 8), & - & real(deg2rad*self%fields(i)%lat(g%isd:g%ied,g%jsd:g%jed), 8), & - & real(deg2rad*lon_out(g%isc:g%iec,g%jsc:g%jec), 8), & - & real(deg2rad*lat_out(g%isc:g%iec,g%jsc:g%jec), 8)) - - ! Make a temporary copy of field - if (allocated(val)) deallocate(val) - allocate(val, mold=self%fields(i)%val) - val = self%fields(i)%val - - ! Interpolate all levels - do k = 1, self%fields(i)%nz - call self%fields(i)%stencil_interp(self%geom, interp2d) - end do - - ! Update c-grid location - self%fields(i)%metadata%grid = cgridlocout - select case(cgridlocout) - ! TODO: Test colocation to u and v grid - !case ('u') - ! self%fields(i)%lon => self%geom%lonu - ! self%fields(i)%lat => self%geom%latu - !case ('v') - ! self%fields(i)%lon => self%geom%lonv - ! self%fields(i)%lat => self%geom%latv - case ('h') - self%fields(i)%lon => self%geom%lon - self%fields(i)%lat => self%geom%lat - end select + ! Interpolate to different location of the stencil + fromto = self%fields(i)%metadata%grid//'toh' + call self%fields(i)%stencil_interp(self%geom, fromto) + call self%fields(i)%update_halo(self%geom) + ! Update grid location to h-points + self%fields(i)%metadata%grid = 'h' + self%fields(i)%lon => self%geom%lon + self%fields(i)%lat => self%geom%lat end do - call horiz_interp_spherical_del(interp2d) - -end subroutine soca_fields_colocate +end subroutine soca_fields_tohpoints ! ------------------------------------------------------------------------------ !> Number of elements to return in the serialized array diff --git a/src/soca/State/State.cc b/src/soca/State/State.cc index 444bf43b7..92f7f3969 100644 --- a/src/soca/State/State.cc +++ b/src/soca/State/State.cc @@ -104,6 +104,21 @@ namespace soca { soca_state_rotate2grid_f90(toFortran(), u, v); } // ----------------------------------------------------------------------------- + /// Staggered grid interpolation + // ----------------------------------------------------------------------------- + void State::tohgrid(const oops::Variables & u, + const oops::Variables & v) const { + Log::trace() << "State::State interpolate vector to h-grid." + << std::endl; + soca_state_tohgrid_f90(toFortran()); + } + // ----------------------------------------------------------------------------- + void State::tocgrid(const oops::Variables & u, + const oops::Variables & v) const { + Log::trace() << "State::State interpolate vector to c-grid. NOT IMPLEMENTED" + << std::endl; + } + // ----------------------------------------------------------------------------- /// Interactions with Increments // ----------------------------------------------------------------------------- State & State::operator+=(const Increment & dx) { diff --git a/src/soca/State/State.h b/src/soca/State/State.h index b2fdc83e1..5e53f20e6 100644 --- a/src/soca/State/State.h +++ b/src/soca/State/State.h @@ -64,6 +64,10 @@ namespace soca { void rotate2north(const oops::Variables &, const oops::Variables &) const; void rotate2grid(const oops::Variables &, const oops::Variables &) const; + /// Staggered grid interpolation + void tohgrid(const oops::Variables &, const oops::Variables &) const; + void tocgrid(const oops::Variables &, const oops::Variables &) const; + /// Logarithmic and exponential transformations void logtrans(const oops::Variables &) const; void expontrans(const oops::Variables &) const; diff --git a/src/soca/State/StateFortran.h b/src/soca/State/StateFortran.h index 1221812ed..09b37bf11 100644 --- a/src/soca/State/StateFortran.h +++ b/src/soca/State/StateFortran.h @@ -43,6 +43,7 @@ namespace soca { void soca_state_rotate2north_f90(const F90flds &, const oops::Variables &, const oops::Variables &); + void soca_state_tohgrid_f90(const F90flds &); void soca_state_logtrans_f90(const F90flds &, const oops::Variables &); void soca_state_expontrans_f90(const F90flds &, const oops::Variables &); void soca_state_gpnorm_f90(const F90flds &, const int &, double &); diff --git a/src/soca/State/soca_state.interface.F90 b/src/soca/State/soca_state.interface.F90 index bbc4fdc36..3a5a48ca6 100644 --- a/src/soca/State/soca_state.interface.F90 +++ b/src/soca/State/soca_state.interface.F90 @@ -256,6 +256,17 @@ subroutine soca_state_rotate2north_c(c_key_self, c_uvars, c_vvars) bind(c,name=' end subroutine soca_state_rotate2north_c +! ------------------------------------------------------------------------------ +!> C++ interface for soca_state_mod::soca_state::tohgrid() +subroutine soca_state_tohgrid_c(c_key_self) bind(c,name='soca_state_tohgrid_f90') + integer(c_int), intent(in) :: c_key_self + + type(soca_state), pointer :: self + + call soca_state_registry%get(c_key_self,self) + call self%tohpoints() + +end subroutine soca_state_tohgrid_c ! ------------------------------------------------------------------------------ !> C++ interface to get soca_state_mod::soca_state dimensions sizes diff --git a/src/soca/State/soca_state_mod.F90 b/src/soca/State/soca_state_mod.F90 index 7158d8daf..5f3c8e3bf 100644 --- a/src/soca/State/soca_state_mod.F90 +++ b/src/soca/State/soca_state_mod.F90 @@ -44,7 +44,7 @@ module soca_state_mod !> \name misc !! \{ - !> \copybrief soca_state_rotate \see soca_state_rotate + !> \copybrief soca_state_rotate \see soca_state_rotate procedure :: rotate => soca_state_rotate !> \copybrief soca_state_convert \see soca_state_convert diff --git a/src/soca/Utils/soca_utils.F90 b/src/soca/Utils/soca_utils.F90 index 6ab26a635..e886da3c5 100644 --- a/src/soca/Utils/soca_utils.F90 +++ b/src/soca/Utils/soca_utils.F90 @@ -16,7 +16,8 @@ module soca_utils private public :: write2pe, soca_str2int, soca_adjust, & - soca_rho, soca_diff, soca_mld, nc_check, soca_remap_idw + soca_rho, soca_diff, soca_mld, nc_check, soca_remap_idw, & + soca_stencil_interp, soca_stencil_neighbors ! ------------------------------------------------------------------------------ contains @@ -232,6 +233,7 @@ subroutine soca_remap_idw(lon_src, lat_src, data_src, lon_dst, lat_dst, data_dst if (dist(n) /= dist(nn_max)) exit nn = n-1 end do + if (nn <= 0 ) call fckit_exception%abort( & "No valid points found in IDW remapping, uh oh.") @@ -257,5 +259,78 @@ subroutine soca_remap_idw(lon_src, lat_src, data_src, lon_dst, lat_dst, data_dst call kd%final() end subroutine soca_remap_idw +! ------------------------------------------------------------------------------ +!> find the 6 stencil neighbors +subroutine soca_stencil_neighbors(fromto, i, j, ij) + character(len=4), intent(in) :: fromto !< "u2h", "v2h" + integer, intent(in) :: i, j + integer, intent(out):: ij(2,6) + + integer :: sti + + select case(fromto) + case("vtoh") + ! return the 6 v-grid neighbors of h-grid(i,j) + ij(1,1) = i; ij(2,1) = j + ij(1,2) = i; ij(2,2) = j-1 + ij(1,3) = i+1; ij(2,3) = j + ij(1,4) = i-1; ij(2,4) = j + ij(1,5) = i-1; ij(2,5) = j-1 + ij(1,6) = i+1; ij(2,6) = j-1 + case("utoh") + ! return the 6 u-grid neighbors of h-grid(i,j) + ij(1,1) = i; ij(2,1) = j + ij(1,2) = i; ij(2,2) = j-1 + ij(1,3) = i; ij(2,3) = j+1 + ij(1,4) = i-1; ij(2,4) = j+1 + ij(1,5) = i-1; ij(2,5) = j + ij(1,6) = i-1; ij(2,6) = j-1 + case default + call fckit_exception%abort(fromto, "option is not implemented") + end select + +end subroutine soca_stencil_neighbors + +! ------------------------------------------------------------------------------ +!> idw interpolation given known neighbors +subroutine soca_stencil_interp(lon_src, lat_src, lon_dst, lat_dst, data, data_out, nn) + real(kind_real), intent(in) :: lon_src(:) + real(kind_real), intent(in) :: lat_src(:) + real(kind_real), intent(in) :: lon_dst + real(kind_real), intent(in) :: lat_dst + real(kind_real), intent(in) :: data(:,:) + real(kind_real), intent(inout) :: data_out(:) + integer, intent(in) :: nn + + type(atlas_geometry) :: ageometry + integer :: i, n, nz, k + real(kind_real) :: val + real(kind_real), allocatable :: w(:) + + ! Initialize atlas geometry on the sphere + ageometry = atlas_geometry("UnitSphere") + + ! nn cannot be larger than 6 + if (nn > 6 ) call fckit_exception%abort( "Using more than 6 neighbors is not allowed") + + ! compute idw weights + allocate(w(nn)) + do i = 1, nn + w(i) = 1_kind_real/ageometry%distance(lon_src(i), lat_src(i), lon_dst, lat_dst) + end do + + nz = size(data, dim=2) + data_out = 0_kind_real + + do k = 1, nz + val = 0_kind_real + do i = 1, nn + val = val + w(i)*data(i,k) + end do + data_out(k) = val/sum(w(1:nn)) + end do + +end subroutine soca_stencil_interp + ! ------------------------------------------------------------------------------ end module soca_utils diff --git a/src/soca/VariableChange/Ana2Model/Ana2Model.cc b/src/soca/VariableChange/Ana2Model/Ana2Model.cc deleted file mode 100644 index fb73fed3e..000000000 --- a/src/soca/VariableChange/Ana2Model/Ana2Model.cc +++ /dev/null @@ -1,73 +0,0 @@ -/* - * (C) Copyright 2017-2021 UCAR. - * - * This software is licensed under the terms of the Apache Licence Version 2.0 - * which can be obtained at http://www.apache.org/licenses/LICENSE-2.0. - */ - - -#include -#include - -#include "eckit/config/Configuration.h" -#include "eckit/exception/Exceptions.h" -#include "oops/util/Logger.h" -#include "soca/Geometry/Geometry.h" -#include "soca/State/State.h" -#include "soca/VariableChange/Ana2Model/Ana2Model.h" - -using oops::Log; - -namespace soca { - -// ----------------------------------------------------------------------------- - -static VariableChangeMaker makerVarChaA2M_("Ana2Model"); - -// ----------------------------------------------------------------------------- -Ana2Model::Ana2Model(const Geometry & resol, const eckit::Configuration & conf) -: uvars_(initRotate(conf, "u")), vvars_(initRotate(conf, "v")), - logvars_(initTrans(conf, "var")) -{ - Log::trace() << "Ana2Model::Ana2Model start" << std::endl; - ASSERT(uvars_.size() == vvars_.size()); - Log::trace() << "Ana2Model::Ana2Model Rotating:" - << " u = " << uvars_ << " v = " << vvars_ << std::endl; - Log::trace() << "Ana2Model::Log Transforming:" - << " var = " << logvars_ << std::endl; - Log::trace() << "Ana2Model::Ana2Model done" << std::endl; -} -// ----------------------------------------------------------------------------- -Ana2Model::~Ana2Model() { - oops::Log::trace() << "ChangeSOCA destructed" << std::endl; -} -// ----------------------------------------------------------------------------- -void Ana2Model::changeVar(const State & xa, - State & xm) const { - oops::Log::trace() << "Ana2Model::changeVar starting" << xa << - std::endl; - util::DateTime * vtime = &xm.validTime(); - xm = xa; - xm.rotate2grid(uvars_, vvars_); - xm.logtrans(logvars_); - xm.validTime() = xa.validTime(); - Log::trace() << "Ana2Model::changeVar done" << xm << std::endl; -} -// ----------------------------------------------------------------------------- -void Ana2Model::changeVarInverse(const State & xm, - State & xa) const { - Log::trace() << "Ana2Model::changeVarInverse starting" < +#include + +#include "eckit/config/Configuration.h" +#include "eckit/exception/Exceptions.h" +#include "oops/util/Logger.h" +#include "soca/Geometry/Geometry.h" +#include "soca/State/State.h" +#include "soca/VariableChange/Model2Ana/Model2Ana.h" + +using oops::Log; + +namespace soca { + +// ----------------------------------------------------------------------------- + +static VariableChangeMaker makerVarChaA2M_("Model2Ana"); + +// ----------------------------------------------------------------------------- +Model2Ana::Model2Ana(const Geometry & resol, const eckit::Configuration & conf) + : uvars_(initRotate(conf, "u")), vvars_(initRotate(conf, "v")), + interp_(initInterp(conf)), logvars_(initTrans(conf, "var")) +{ + Log::trace() << "Model2Ana::Model2Ana start" << std::endl; + ASSERT(uvars_.size() == vvars_.size()); + Log::trace() << "Model2Ana::Model2Ana Rotating:" + << " u = " << uvars_ << " v = " << vvars_ << std::endl; + Log::trace() << "Model2Ana::Model2Ana Interpolate:" + << interp_ << std::endl; + Log::trace() << "Model2Ana::Log Transforming:" + << " var = " << logvars_ << std::endl; + Log::trace() << "Model2Ana::Model2Ana done" << std::endl; +} +// ----------------------------------------------------------------------------- +Model2Ana::~Model2Ana() { + oops::Log::trace() << "ChangeSOCA destructed" << std::endl; +} +// ----------------------------------------------------------------------------- +void Model2Ana::changeVar(const State & xm, + State & xa) const { + Log::trace() << "Model2Ana::changeVar starting" < initTrans(const eckit::Configuration & conf, const std::string & trvar) const { @@ -56,6 +59,7 @@ class Ana2Model: public VariableChangeBase { void print(std::ostream &) const override; const oops::Variables uvars_; const oops::Variables vvars_; + const bool interp_; const oops::Variables logvars_; }; // ----------------------------------------------------------------------------- diff --git a/src/soca/VariableChange/VariableChange.cc b/src/soca/VariableChange/VariableChange.cc index 2d13dc288..48e3a956c 100644 --- a/src/soca/VariableChange/VariableChange.cc +++ b/src/soca/VariableChange/VariableChange.cc @@ -42,7 +42,7 @@ void VariableChange::changeVar(State & x, const oops::Variables & vars) const { << vars << std::endl; // TODO(travis) rename in/out variables so that skipping this - // works for Ana2Model (i.e. we need rotated/unrotate u/v renamed different) + // works for Model2Ana (i.e. we need rotated/unrotate u/v renamed different) // // If the variables are the same, don't bother doing anything! // if (!(x.variables() == vars)) { diff --git a/test/testinput/3dhybfgat.yml b/test/testinput/3dhybfgat.yml index 363efc186..c18739477 100644 --- a/test/testinput/3dhybfgat.yml +++ b/test/testinput/3dhybfgat.yml @@ -17,7 +17,7 @@ cost function: model variables: &model_vars [cicen, hicen, socn, uocn, vocn, tocn, ssh, hocn, sw, lhf, shf, lw, us, mld, layer_depth] variable change: - variable change name: Ana2Model + variable change name: Model2Ana rotate: u: [uocn] v: [vocn] diff --git a/test/testinput/3dvarfgat.yml b/test/testinput/3dvarfgat.yml index a44f8ca5e..5e56fa089 100644 --- a/test/testinput/3dvarfgat.yml +++ b/test/testinput/3dvarfgat.yml @@ -16,7 +16,7 @@ cost function: fields metadata: data_static/fields_metadata.yml variable change: - variable change name: Ana2Model + variable change name: Model2Ana rotate: u: [uocn] v: [vocn] diff --git a/test/testinput/3dvarfgat_pseudo.yml b/test/testinput/3dvarfgat_pseudo.yml index 963d46400..b457272df 100644 --- a/test/testinput/3dvarfgat_pseudo.yml +++ b/test/testinput/3dvarfgat_pseudo.yml @@ -16,7 +16,7 @@ cost function: fields metadata: data_static/fields_metadata.yml variable change: - variable change name: Ana2Model + variable change name: Model2Ana rotate: u: [uocn] v: [vocn] diff --git a/test/testinput/convertstate_changevar.yml b/test/testinput/convertstate_changevar.yml index 5cb8012b9..70b2c2a75 100644 --- a/test/testinput/convertstate_changevar.yml +++ b/test/testinput/convertstate_changevar.yml @@ -10,11 +10,12 @@ input geometry: &geom output geometry: *geom variable change: - variable change name: Ana2Model + variable change name: Model2Ana do inverse: false rotate: u: [uocn] v: [vocn] + interp: true log: var: [socn, chl] output variables: *soca_vars diff --git a/test/testinput/varchange_ana2model.yml b/test/testinput/varchange_ana2model.yml index 9a25b90ca..f7bf2151c 100644 --- a/test/testinput/varchange_ana2model.yml +++ b/test/testinput/varchange_ana2model.yml @@ -5,7 +5,7 @@ geometry: variable change tests: - variable change: - variable change name: Ana2Model + variable change name: Model2Ana state: read_from_file: 1 date: 2018-04-15T00:00:00Z diff --git a/test/testref/3dvarfgat.test b/test/testref/3dvarfgat.test index 6300630fa..0a8d73a61 100644 --- a/test/testref/3dvarfgat.test +++ b/test/testref/3dvarfgat.test @@ -1,23 +1,23 @@ CostJb : Nonlinear Jb = 0.0000000000000000e+00 -CostJo : Nonlinear Jo(CoolSkin) = 2.1763431423226604e+04, nobs = 48, Jo/n = 4.5340482131722092e+02, err = 3.8763170936361402e-01 -CostJo : Nonlinear Jo(SeaSurfaceTemp) = 3.3393856011145090e+02, nobs = 43, Jo/n = 7.7660130258476956e+00, err = 3.8581536591328441e-01 -CostJo : Nonlinear Jo(SeaSurfaceSalinity) = 1.7511224050621064e+00, nobs = 15, Jo/n = 1.1674149367080709e-01, err = 1.0000000000000000e+00 -CostJo : Nonlinear Jo(ADT) = 7.9943052246258645e+01, nobs = 29, Jo/n = 2.7566569740089188e+00, err = 1.0000000149011612e-01 +CostJo : Nonlinear Jo(CoolSkin) = 2.1763431423226597e+04, nobs = 48, Jo/n = 4.5340482131722075e+02, err = 3.8763170936361402e-01 +CostJo : Nonlinear Jo(SeaSurfaceTemp) = 3.3393856011141156e+02, nobs = 43, Jo/n = 7.7660130258467808e+00, err = 3.8581536591328441e-01 +CostJo : Nonlinear Jo(SeaSurfaceSalinity) = 1.7511224050619130e+00, nobs = 15, Jo/n = 1.1674149367079420e-01, err = 1.0000000000000000e+00 +CostJo : Nonlinear Jo(ADT) = 7.9943052302124300e+01, nobs = 29, Jo/n = 2.7566569759353206e+00, err = 1.0000000149011612e-01 CostJo : Nonlinear Jo(InsituTemperature) = 3.3968987169654760e+01, nobs = 30, Jo/n = 1.1322995723218254e+00, err = 9.1590378547213480e-01 CostJo : Nonlinear Jo(InsituSalinity) = 1.5845797862881717e+01, nobs = 33, Jo/n = 4.8017569281459749e-01, err = 6.1042996877640276e-01 CostJo : Nonlinear Jo(SeaIceFraction) = 2.3572786998945503e+02, nobs = 24, Jo/n = 9.8219945828939590e+00, err = 1.0000000149011612e-01 -CostJo : Nonlinear Jo(SurfaceU) = 2.0640750953699122e-01, nobs = 48, Jo/n = 4.3001564486873168e-03, err = 3.8763169345842458e-01 -CostFunction: Nonlinear J = 2.2464813220520904e+04 -RPCGMinimizer: reduction in residual norm = 1.0984350652758776e+00 +CostJo : Nonlinear Jo(SurfaceU) = 2.0640750940226621e-01, nobs = 48, Jo/n = 4.3001564458805463e-03, err = 3.8763169345842458e-01 +CostFunction: Nonlinear J = 2.2464813220576591e+04 +RPCGMinimizer: reduction in residual norm = 1.0984350652739356e+00 CostFunction::addIncrement: Analysis: Valid time: 2018-04-15T00:00:00Z - cicen min=-0.0286935633424766 max=1.0484531854277295 mean=0.1201656281496799 + cicen min=-0.0286935633424276 max=1.0484531854275954 mean=0.1201656281496798 hicen min=0.0000000000000000 max=4.0326673084246947 mean=0.4712515705773916 - socn min=10.7210460395083924 max=40.4416591897031168 mean=34.5458126924755220 - tocn min=-1.8883899372702533 max=65.1953082854073216 mean=6.0662982865278305 - uocn min=-0.8581688276833758 max=0.6999167056417410 mean=-0.0002591636654831 + socn min=10.7210460395083924 max=40.4416591897031168 mean=34.5458126924452742 + tocn min=-1.8883899372702533 max=65.1953082853449075 mean=6.0662982865389417 + uocn min=-0.8581688276833758 max=0.6999167056417414 mean=-0.0002591636654830 vocn min=-0.7661101215480253 max=1.4377766409421606 mean=0.0021972630399530 - ssh min=-1.9191227571636915 max=0.9230126424816942 mean=-0.2758628641082930 + ssh min=-1.9191227571581164 max=0.9230126425406677 mean=-0.2758628640973802 hocn min=0.0009999999999977 max=1345.6400000000003274 mean=128.6280642065023017 sw min=0.0000000000000000 max=0.0000000000000000 mean=0.0000000000000000 lhf min=0.0000000000000000 max=0.0000000000000000 mean=0.0000000000000000 @@ -36,13 +36,13 @@ layer_depth min=2.2854716757984130 max=5658.3057467114012979 mean=1200.522 -CostJb : Nonlinear Jb = 4900.5505544964244109 -CostJo : Nonlinear Jo(CoolSkin) = 25202.4168784580251668, nobs = 48, Jo/n = 525.0503516345421531, err = 0.3876317093636140 -CostJo : Nonlinear Jo(SeaSurfaceTemp) = 3569.3773977796836334, nobs = 43, Jo/n = 83.0087766925507822, err = 0.3858153659132844 -CostJo : Nonlinear Jo(SeaSurfaceSalinity) = 1.7426599076872593, nobs = 15, Jo/n = 0.1161773271791506, err = 1.0000000000000000 -CostJo : Nonlinear Jo(ADT) = 78.7688608809395703, nobs = 29, Jo/n = 2.7161676165841233, err = 0.1000000014901161 -CostJo : Nonlinear Jo(InsituTemperature) = 30.6442281578270510, nobs = 30, Jo/n = 1.0214742719275685, err = 0.9159037854721348 -CostJo : Nonlinear Jo(InsituSalinity) = 15.5656968088175152, nobs = 33, Jo/n = 0.4716877820853793, err = 0.6104299687764028 -CostJo : Nonlinear Jo(SeaIceFraction) = 200.2061460146289278, nobs = 24, Jo/n = 8.3419227506095392, err = 0.1000000014901161 -CostJo : Nonlinear Jo(SurfaceU) = 0.2049582352893839, nobs = 48, Jo/n = 0.0042699632351955, err = 0.3876316934584246 -CostFunction: Nonlinear J = 33999.4773807393212337 +CostJb : Nonlinear Jb = 4900.5505544617581108 +CostJo : Nonlinear Jo(CoolSkin) = 25202.4168780031068309, nobs = 48, Jo/n = 525.0503516250647635, err = 0.3876317093636140 +CostJo : Nonlinear Jo(SeaSurfaceTemp) = 3569.3773977347736945, nobs = 43, Jo/n = 83.0087766915063696, err = 0.3858153659132844 +CostJo : Nonlinear Jo(SeaSurfaceSalinity) = 1.7426599076836851, nobs = 15, Jo/n = 0.1161773271789123, err = 1.0000000000000000 +CostJo : Nonlinear Jo(ADT) = 78.7701054794583797, nobs = 29, Jo/n = 2.7162105337744267, err = 0.1000000014901161 +CostJo : Nonlinear Jo(InsituTemperature) = 30.6442281576848501, nobs = 30, Jo/n = 1.0214742719228282, err = 0.9159037854721348 +CostJo : Nonlinear Jo(InsituSalinity) = 15.5656968088179433, nobs = 33, Jo/n = 0.4716877820853922, err = 0.6104299687764028 +CostJo : Nonlinear Jo(SeaIceFraction) = 200.2061460146998684, nobs = 24, Jo/n = 8.3419227506124951, err = 0.1000000014901161 +CostJo : Nonlinear Jo(SurfaceU) = 0.2049585774054072, nobs = 48, Jo/n = 0.0042699703626126, err = 0.3876316934584246 +CostFunction: Nonlinear J = 33999.4786251453915611 diff --git a/test/testref/convertstate_changevar.test b/test/testref/convertstate_changevar.test index 7c384f827..9a7119ca8 100644 --- a/test/testref/convertstate_changevar.test +++ b/test/testref/convertstate_changevar.test @@ -9,7 +9,7 @@ Input state: cicen min=0.0000000000000000 max=1.0000000000000000 mean=0.1175125761012928 chl min=0.0000100000000000 max=4.6719899196170607 mean=0.1184829680403295 layer_depth min=2.2854716757984130 max=5658.3057467114012979 mean=1200.5229536158392420 -Variable transform: Ana2Model +Variable transform: Model2Ana Variable change from: 9 variables: ssh, tocn, socn, uocn, vocn, hocn, cicen, chl, layer_depth Variable change to: 9 variables: ssh, tocn, socn, uocn, vocn, hocn, cicen, chl, layer_depth State after variable transform: @@ -17,8 +17,8 @@ State after variable transform: ssh min=-1.9244847628277935 max=0.9272826517867588 mean=-0.2767903423591662 tocn min=-1.8883899372702533 max=31.7004645720658580 mean=6.0175649583533595 socn min=2.3722088224690281 max=3.6998604464674543 mean=3.5415386142981689 - uocn min=-0.8581693992488438 max=0.7000954286848975 mean=-0.0002180189769953 - vocn min=-0.7661101215480253 max=1.4375780595495007 mean=0.0020985163615815 + uocn min=-0.5496338818918315 max=0.3909715861579679 mean=0.0001603778035724 + vocn min=-0.6585294260094935 max=0.6907511741595780 mean=0.0020870869173532 hocn min=0.0009999999999977 max=1345.6400000000003274 mean=128.6280642065023017 cicen min=0.0000000000000000 max=1.0000000000000000 mean=0.1175125761012928 chl min=-11.4176152851659030 max=1.5415853019489218 mean=-4.1636979499256563 @@ -28,8 +28,8 @@ Output state: ssh min=-1.9244847628277935 max=0.9272826517867588 mean=-0.2767903423591662 tocn min=-1.8883899372702533 max=31.7004645720658580 mean=6.0175649583533595 socn min=2.3722088224690281 max=3.6998604464674543 mean=3.5415386142981689 - uocn min=-0.8581693992488438 max=0.7000954286848975 mean=-0.0002180189769953 - vocn min=-0.7661101215480253 max=1.4375780595495007 mean=0.0020985163615815 + uocn min=-0.5496338818918315 max=0.3909715861579679 mean=0.0001603778035724 + vocn min=-0.6585294260094935 max=0.6907511741595780 mean=0.0020870869173532 hocn min=0.0009999999999977 max=1345.6400000000003274 mean=128.6280642065023017 cicen min=0.0000000000000000 max=1.0000000000000000 mean=0.1175125761012928 chl min=-11.4176152851659030 max=1.5415853019489218 mean=-4.1636979499256563 @@ -45,7 +45,7 @@ Input state: cicen min=0.0000000000000000 max=1.0000000000000000 mean=0.1175125761012928 chl min=0.0000100000000000 max=4.6719899196170607 mean=0.1184829680403295 layer_depth min=2.2854716757984130 max=5658.3057467114012979 mean=1200.5229536158392420 -Variable transform: Ana2Model +Variable transform: Model2Ana Variable change from: 9 variables: ssh, tocn, socn, uocn, vocn, hocn, cicen, chl, layer_depth Variable change to: 9 variables: ssh, tocn, socn, uocn, vocn, hocn, cicen, chl, layer_depth State after variable transform: @@ -53,8 +53,8 @@ State after variable transform: ssh min=-1.9244847628277935 max=0.9272826517867588 mean=-0.2767903423591662 tocn min=-1.8883899372702533 max=31.7004645720658580 mean=6.0175649583533595 socn min=2.3722088224690281 max=3.6998604464674543 mean=3.5415386142981689 - uocn min=-0.8581693992488438 max=0.7000954286848975 mean=-0.0002180189769953 - vocn min=-0.7661101215480253 max=1.4375780595495007 mean=0.0020985163615815 + uocn min=-0.5496338818918315 max=0.3909715861579679 mean=0.0001603778035724 + vocn min=-0.6585294260094935 max=0.6907511741595780 mean=0.0020870869173532 hocn min=0.0009999999999977 max=1345.6400000000003274 mean=128.6280642065023017 cicen min=0.0000000000000000 max=1.0000000000000000 mean=0.1175125761012928 chl min=-11.4176152851659030 max=1.5415853019489218 mean=-4.1636979499256563 @@ -64,8 +64,8 @@ Output state: ssh min=-1.9244847628277935 max=0.9272826517867588 mean=-0.2767903423591662 tocn min=-1.8883899372702533 max=31.7004645720658580 mean=6.0175649583533595 socn min=2.3722088224690281 max=3.6998604464674543 mean=3.5415386142981689 - uocn min=-0.8581693992488438 max=0.7000954286848975 mean=-0.0002180189769953 - vocn min=-0.7661101215480253 max=1.4375780595495007 mean=0.0020985163615815 + uocn min=-0.5496338818918315 max=0.3909715861579679 mean=0.0001603778035724 + vocn min=-0.6585294260094935 max=0.6907511741595780 mean=0.0020870869173532 hocn min=0.0009999999999977 max=1345.6400000000003274 mean=128.6280642065023017 cicen min=0.0000000000000000 max=1.0000000000000000 mean=0.1175125761012928 chl min=-11.4176152851659030 max=1.5415853019489218 mean=-4.1636979499256563