diff --git a/tools/cubed_sphere_inc_mod.F90 b/tools/cubed_sphere_inc_mod.F90 index 2b169b6e..19f80742 100644 --- a/tools/cubed_sphere_inc_mod.F90 +++ b/tools/cubed_sphere_inc_mod.F90 @@ -4,8 +4,9 @@ module cubed_sphere_inc_mod use field_manager_mod, only: MODEL_ATMOS use fv_arrays_mod, only: fv_atmos_type use fms2_io_mod, only: open_file, close_file, read_data, variable_exists, & - FmsNetcdfDomainFile_t, register_axis - + FmsNetcdfDomainFile_t, register_axis, FmsNetcdfFile_t + use mpp_mod, only: mpp_error, NOTE + implicit none type increment_data_type real, allocatable :: ua_inc(:,:,:) @@ -22,24 +23,37 @@ module cubed_sphere_inc_mod !---------------------------------------------------------------------------------------- - subroutine read_cubed_sphere_inc(fname, increment_data, Atm) + subroutine read_cubed_sphere_inc(fname, increment_data, Atm, IAU_regional) character(*), intent(in) :: fname type(increment_data_type), intent(inout) :: increment_data type(fv_atmos_type), intent(in) :: Atm + logical, intent(in) :: IAU_regional type(FmsNetcdfDomainFile_t) :: fileobj + type(FmsNetcdfFile_t) :: fileobj_regional integer :: itracer, ntracers character(len=64) :: tracer_name - + character(len=:), allocatable :: fname_base + integer :: ipos ! Get various dimensions call get_number_tracers(MODEL_ATMOS, num_tracers=ntracers) ! Open file - if ( open_file(fileobj, trim(fname), 'read', Atm%domain) ) then + if (IAU_regional) then + ipos = index(fname, '.') + if (ipos > 0) then + fname_base = fname(1:ipos-1) + else + fname_base = trim(fname) + endif + fname_base=trim(fname_base)//'.nc' + else + fname_base=trim(fname) + end if + if ( open_file(fileobj, trim(fname_base), 'read', Atm%domain) ) then ! Register axes call register_axis(fileobj, 'xaxis_1', 'x') call register_axis(fileobj, 'yaxis_1', 'y') - ! Read increments call read_data(fileobj, 'u_inc', increment_data%ua_inc) call read_data(fileobj, 'v_inc', increment_data%va_inc) @@ -54,12 +68,13 @@ subroutine read_cubed_sphere_inc(fname, increment_data, Atm) call get_tracer_names(MODEL_ATMOS, itracer, tracer_name) if ( variable_exists(fileobj, trim(tracer_name)//'_inc') ) then call read_data(fileobj, trim(tracer_name)//'_inc', increment_data%tracer_inc(:,:,:,itracer)) + else + call mpp_error(NOTE, 'No Increment for '//trim(tracer_name)//' found, assuming zero') + increment_data%tracer_inc(:,:,:,itracer) = 0.0 end if end do - call close_file(fileobj) end if - end subroutine read_cubed_sphere_inc !---------------------------------------------------------------------------------------- diff --git a/tools/fv_iau_mod.F90 b/tools/fv_iau_mod.F90 index 19375a37..38aa3e14 100644 --- a/tools/fv_iau_mod.F90 +++ b/tools/fv_iau_mod.F90 @@ -182,9 +182,15 @@ subroutine IAU_initialize (IPD_Control, IAU_Data, Init_parm, Atm) if( file_exists(fname) ) then call open_ncfile( fname, ncid ) ! open the file - call get_ncdim1( ncid, 'lon', im) - call get_ncdim1( ncid, 'lat', jm) - call get_ncdim1( ncid, 'lev', km) + if (IPD_Control%iau_regional) then + call get_ncdim1( ncid, 'xaxis_1', im) + call get_ncdim1( ncid, 'yaxis_1', jm) + call get_ncdim1( ncid, 'zaxis_1', km) + else + call get_ncdim1( ncid, 'lon', im) + call get_ncdim1( ncid, 'lat', jm) + call get_ncdim1( ncid, 'lev', km) + end if if (km.ne.npz) then if (is_master()) print *, 'km = ', km @@ -194,20 +200,22 @@ subroutine IAU_initialize (IPD_Control, IAU_Data, Init_parm, Atm) if(is_master()) write(*,*) fname, ' DA increment dimensions:', im,jm,km - allocate ( lon(im) ) - allocate ( lat(jm) ) + if (.not.IPD_Control%iau_regional) then + allocate ( lon(im) ) + allocate ( lat(jm) ) - call _GET_VAR1 (ncid, 'lon', im, lon ) - call _GET_VAR1 (ncid, 'lat', jm, lat ) - call close_ncfile(ncid) + call _GET_VAR1 (ncid, 'lon', im, lon ) + call _GET_VAR1 (ncid, 'lat', jm, lat ) + call close_ncfile(ncid) - ! Convert to radians - do i=1,im - lon(i) = lon(i) * deg2rad - enddo - do j=1,jm - lat(j) = lat(j) * deg2rad - enddo + ! Convert to radians + do i=1,im + lon(i) = lon(i) * deg2rad + enddo + do j=1,jm + lat(j) = lat(j) * deg2rad + enddo + end if else call mpp_error(FATAL,'==> Error in IAU_initialize: Expected file '& @@ -226,10 +234,12 @@ subroutine IAU_initialize (IPD_Control, IAU_Data, Init_parm, Atm) agrid(is-1+i,js-1+j,2)=Init_parm%xlat(i,j) enddo enddo - call remap_coef( is, ie, js, je, is, ie, js, je, & + if (.not.IPD_Control%iau_regional) then + call remap_coef( is, ie, js, je, is, ie, js, je, & im, jm, lon, lat, id1, id2, jdc, s2c, & agrid) - deallocate ( lon, lat,agrid ) + deallocate ( lon, lat,agrid ) + end if allocate(IAU_Data%ua_inc(is:ie, js:je, km)) @@ -247,11 +257,12 @@ subroutine IAU_initialize (IPD_Control, IAU_Data, Init_parm, Atm) allocate (iau_state%inc1%tracer_inc(is:ie, js:je, km,ntracers)) iau_state%hr1=IPD_Control%iaufhrs(1) iau_state%wt = 1.0 ! IAU increment filter weights (default 1.0) + iau_state%wt=iau_state%wt*IPD_Control%iau_inc_scale !Increase the weight iau_state%wt_normfact = 1.0 if (IPD_Control%iau_filter_increments) then ! compute increment filter weights, sum to obtain normalization factor dtp=IPD_control%dtp - nstep = 0.5*IPD_Control%iau_delthrs*3600/dtp + nstep = nint(0.5*IPD_Control%iau_delthrs*3600/dtp) ! compute normalization factor for filter weights normfact = 0. do k=1,2*nstep+1 @@ -268,8 +279,9 @@ subroutine IAU_initialize (IPD_Control, IAU_Data, Init_parm, Atm) enddo iau_state%wt_normfact = (2*nstep+1)/normfact endif + if (IPD_Control%iau_regional) iau_state%wt_normfact=iau_state%wt_normfact*IPD_Control%iau_inc_scale if ( Atm%flagstruct%increment_file_on_native_grid ) then - call read_cubed_sphere_inc('INPUT/'//trim(IPD_Control%iau_inc_files(1)), iau_state%inc1, Atm) + call read_cubed_sphere_inc('INPUT/'//trim(IPD_Control%iau_inc_files(1)), iau_state%inc1, Atm, IPD_Control%iau_regional) else call read_iau_forcing(IPD_Control,iau_state%inc1,'INPUT/'//trim(IPD_Control%iau_inc_files(1))) endif @@ -285,7 +297,7 @@ subroutine IAU_initialize (IPD_Control, IAU_Data, Init_parm, Atm) allocate (iau_state%inc2%tracer_inc(is:ie, js:je, km,ntracers)) iau_state%hr2=IPD_Control%iaufhrs(2) if ( Atm%flagstruct%increment_file_on_native_grid ) then - call read_cubed_sphere_inc('INPUT/'//trim(IPD_Control%iau_inc_files(2)), iau_state%inc2, Atm) + call read_cubed_sphere_inc('INPUT/'//trim(IPD_Control%iau_inc_files(2)), iau_state%inc2, Atm, IPD_Control%iau_regional) else call read_iau_forcing(IPD_Control,iau_state%inc2,'INPUT/'//trim(IPD_Control%iau_inc_files(2))) endif @@ -323,9 +335,9 @@ subroutine getiauforcing(IPD_Control,IAU_Data,Atm) ! in window kstep=-nstep,nstep (2*nstep+1 total) ! time step IPD_control%dtp dtp=IPD_control%dtp - nstep = 0.5*IPD_Control%iau_delthrs*3600/dtp + nstep = nint(0.5*IPD_Control%iau_delthrs*3600/dtp) ! compute normalized filter weight - kstep = ((IPD_Control%fhour-t1) - 0.5*IPD_Control%iau_delthrs)*3600./dtp + kstep = nint(((IPD_Control%fhour-t1) - 0.5*IPD_Control%iau_delthrs)*3600./dtp) if (IPD_Control%fhour >= t1 .and. IPD_Control%fhour < t2) then sx = acos(-1.)*kstep/nstep wx = acos(-1.)*kstep/(nstep+1) @@ -375,7 +387,7 @@ subroutine getiauforcing(IPD_Control,IAU_Data,Atm) iau_state%inc1=iau_state%inc2 if (is_master()) print *,'reading next increment file',trim(IPD_Control%iau_inc_files(itnext)) if ( Atm%flagstruct%increment_file_on_native_grid ) then - call read_cubed_sphere_inc('INPUT/'//trim(IPD_Control%iau_inc_files(itnext)), iau_state%inc2, Atm) + call read_cubed_sphere_inc('INPUT/'//trim(IPD_Control%iau_inc_files(itnext)), iau_state%inc2, Atm, IPD_Control%iau_regional) else call read_iau_forcing(IPD_Control,iau_state%inc2,'INPUT/'//trim(IPD_Control%iau_inc_files(itnext))) endif @@ -481,7 +493,12 @@ subroutine read_iau_forcing(IPD_Control,increments,fname) enddo enddo - allocate ( wk3(1:im,jbeg:jend, 1:km) ) + if (IPD_Control%iau_regional) then + allocate ( wk3(1:im,js:je, 1:km) ) + else + allocate ( wk3(1:im,jbeg:jend, 1:km) ) + endif + ! read in 1 time level call interp_inc('T_inc',increments%temp_inc(:,:,:),jbeg,jend) call interp_inc('delp_inc',increments%delp_inc(:,:,:),jbeg,jend) diff --git a/tools/fv_treat_da_inc.F90 b/tools/fv_treat_da_inc.F90 index 7e985fb1..6d749bf1 100644 --- a/tools/fv_treat_da_inc.F90 +++ b/tools/fv_treat_da_inc.F90 @@ -522,7 +522,7 @@ subroutine read_da_inc_cubed_sphere(Atm, fv_domain, bd, npz_in, nq, & ! Read increments fname = trim(fname_prefix) // '.nc' - call read_cubed_sphere_inc(fname, increment_data, Atm) + call read_cubed_sphere_inc(fname, increment_data, Atm, .False.) ! Wind increments ! ---------------