diff --git a/lis/configs/lis.config.adoc b/lis/configs/lis.config.adoc index 8a9eeb7c2..0f1c2d917 100644 --- a/lis/configs/lis.config.adoc +++ b/lis/configs/lis.config.adoc @@ -396,6 +396,7 @@ endif::devonly[] |"`COAMPSout`" | COAMPS output forcing |"`GALWEM forecast`" | GALWEM 17km or 0.25deg deterministic forecast reader |"`GALWEM-GE forecast`" | GALWEM-GE forecast reader +|"`Global EIS`" | Global Earth Information System (EIS) 5km forcing |==== .Example _lis.config_ entry @@ -6635,8 +6636,19 @@ of the GALWEM-GE forecast forcing data. GALWEM-GE forecast forcing directory: ./GALWEM_GE GALWEM-GE forecast run mode: forecast # forecast | analysis GALWEM-GE forecast number of ensemble members: 21 +.... + +[[sssec_geis,Global EIS]] +==== Global EIS + +NOTE: There is no spatial interpolation support for this dataset. +`Global EIS forcing directory:` specifies the location +of the Global EIS forcing files. +.Example _lis.config_ entry +.... +Global EIS forcing directory: ./geis .... [[ssec_lsm,Land surface models]] diff --git a/lis/make/default.cfg b/lis/make/default.cfg index 87a51f553..acee04406 100644 --- a/lis/make/default.cfg +++ b/lis/make/default.cfg @@ -267,6 +267,11 @@ enabled: True macro: MF_NLDAS2 path: metforcing/nldas2 +[GEIS] +enabled: True +macro: MF_GEIS +path: metforcing/geis + [GLDAS] enabled: True macro: MF_GLDAS diff --git a/lis/metforcing/geis/finalize_geis.F90 b/lis/metforcing/geis/finalize_geis.F90 new file mode 100644 index 000000000..6e0abd19b --- /dev/null +++ b/lis/metforcing/geis/finalize_geis.F90 @@ -0,0 +1,35 @@ +!-----------------------BEGIN NOTICE -- DO NOT EDIT----------------------- +! NASA Goddard Space Flight Center +! Land Information System Framework (LISF) +! Version 7.4 +! +! Copyright (c) 2022 United States Government as represented by the +! Administrator of the National Aeronautics and Space Administration. +! All Rights Reserved. +!-------------------------END NOTICE -- DO NOT EDIT----------------------- +#include "LIS_misc.h" +!BOP +! !MODULE: finalize_geis +! \label{finalize_geis} +! +! !REVISION HISTORY: +! 15 Nov 2023; Sujay Kumar, Initial Code +! +! !INTERFACE: +subroutine finalize_geis(findex) +! !USES: + use LIS_coreMod, only : LIS_rc + use geis_forcingMod, only : geis_struc +! +! !DESCRIPTION: +! Routine to cleanup geis forcing related memory allocations. +! +!EOP + implicit none + + integer :: n + integer :: findex + + deallocate(geis_struc) + +end subroutine finalize_geis diff --git a/lis/metforcing/geis/geis_forcingMod.F90 b/lis/metforcing/geis/geis_forcingMod.F90 new file mode 100644 index 000000000..5cb84a002 --- /dev/null +++ b/lis/metforcing/geis/geis_forcingMod.F90 @@ -0,0 +1,142 @@ +!-----------------------BEGIN NOTICE -- DO NOT EDIT----------------------- +! NASA Goddard Space Flight Center +! Land Information System Framework (LISF) +! Version 7.4 +! +! Copyright (c) 2022 United States Government as represented by the +! Administrator of the National Aeronautics and Space Administration. +! All Rights Reserved. +!-------------------------END NOTICE -- DO NOT EDIT----------------------- +module geis_forcingMod +!BOP +! !MODULE: geis_forcingMod +! +! !DESCRIPTION: +! This module contains variables and data structures that are used +! for the implementation of the global forcing data used in the +! Earth Information System (EIS). The variables are produced +! at 5km spatial resolution, and at hourly intervals. +! +! NOTE: no spatial interpolation of the code is supported currently +! +! !REVISION HISTORY: +! 15 Nov 2023: Sujay Kumar; Initial Specification +! +! !USES: + use LIS_constantsMod, only : LIS_CONST_PATH_LEN + implicit none + + PRIVATE +!----------------------------------------------------------------------------- +! !PUBLIC MEMBER FUNCTIONS: +!----------------------------------------------------------------------------- + public :: init_GEIS !defines the native resolution of + !the input data +!----------------------------------------------------------------------------- +! !PUBLIC TYPES: +!----------------------------------------------------------------------------- + public :: geis_struc +!EOP + + type, public :: geis_type_dec + real :: ts + integer :: ncold, nrold ! AWIPS 212 dimensions + character*50 :: geis_filesrc + character(len=LIS_CONST_PATH_LEN) :: geisdir ! NLDAS-2 Forcing Directory + real*8 :: geistime1,geistime2 + + integer :: findtime1, findtime2 + integer :: niter + integer :: c_off,r_off + real, allocatable :: metdata1(:,:,:) + real, allocatable :: metdata2(:,:,:) + + end type geis_type_dec + + type(geis_type_dec), allocatable :: geis_struc(:) +!EOP +contains + +!BOP +! +! !ROUTINE: init_GEIS +! \label{init_GEIS} +! +! !INTERFACE: + subroutine init_GEIS(findex) +! !USES: + use LIS_coreMod, only : LIS_rc, LIS_domain + use LIS_timeMgrMod, only : LIS_update_timestep + use LIS_logMod, only : LIS_logunit,LIS_endrun + use LIS_spatialDownscalingMod, only : LIS_init_pcpclimo_native + use map_utils, only : proj_latlon + use LIS_forecastMod + + implicit none + + integer, intent(in) :: findex +! +! !DESCRIPTION: +! Defines the native resolution of the input forcing for NLDAS-2 +! data. The grid description arrays are based on the decoding +! schemes used by NCEP and followed in the LIS interpolation +! schemes (see Section~\ref{interp}). +! +! The routines invoked are: +! \begin{description} +! \item[readcrd\_geis](\ref{readcrd_geis}) \newline +! reads the runtime options specified for NLDAS-2 data +! \item[bilinear\_interp\_input](\ref{bilinear_interp_input}) \newline +! computes the neighbor, weights for bilinear interpolation +! \item[conserv\_interp\_input](\ref{conserv_interp_input}) \newline +! computes the neighbor, weights for conservative interpolation +! \item[read\_geis\_elev](\ref{read_geis_elev}) \newline +! reads the native elevation of the NLDAS-2 data to be used +! for topographic adjustments to the forcing +! \end{description} +!EOP + + integer :: n + real :: lat1, lon1 + + allocate(geis_struc(LIS_rc%nnest)) + call readcrd_geis() + + do n=1, LIS_rc%nnest + geis_struc(n)%ts = 3600 + call LIS_update_timestep(LIS_rc, n, geis_struc(n)%ts) + enddo + + LIS_rc%met_nf(findex) = 7 + + ! Set NLDAS-2 grid dimensions and extent information: + geis_struc(:)%ncold = 7200 + geis_struc(:)%nrold = 3000 + + do n=1,LIS_rc%nnest + + allocate(geis_struc(n)%metdata1(1,& + LIS_rc%met_nf(findex),& + LIS_rc%ngrid(n))) + allocate(geis_struc(n)%metdata2(1,& + LIS_rc%met_nf(findex),& + LIS_rc%ngrid(n))) + + geis_struc(n)%metdata1 = 0 + geis_struc(n)%metdata2 = 0 + + geis_struc(n)%findtime1 = 0 + geis_struc(n)%findtime2 = 0 + + geis_struc(n)%nIter = 1 + + lat1 = LIS_domain(n)%lat(1) + lon1 = LIS_domain(n)%lon(1) + + geis_struc(n)%r_off = nint((lat1 + 59.975)/0.05) + 1 + geis_struc(n)%c_off = nint((lon1 + 179.975)/0.05) + 1 + + enddo + + end subroutine init_GEIS +end module geis_forcingMod diff --git a/lis/metforcing/geis/get_geis.F90 b/lis/metforcing/geis/get_geis.F90 new file mode 100644 index 000000000..f97449d58 --- /dev/null +++ b/lis/metforcing/geis/get_geis.F90 @@ -0,0 +1,268 @@ +!-----------------------BEGIN NOTICE -- DO NOT EDIT----------------------- +! NASA Goddard Space Flight Center +! Land Information System Framework (LISF) +! Version 7.4 +! +! Copyright (c) 2022 United States Government as represented by the +! Administrator of the National Aeronautics and Space Administration. +! All Rights Reserved. +!-------------------------END NOTICE -- DO NOT EDIT----------------------- +!BOP +! !ROUTINE: get_geis +! \label{get_geis} +! +! +! !REVISION HISTORY: +! 15 Nov 2023: Sujay Kumar; Initial Version in LIS +! +! !INTERFACE: +subroutine get_geis(n,findex) +! !USES: + use LIS_coreMod, only : LIS_rc, LIS_domain + use LIS_timeMgrMod, only : LIS_tick + use LIS_metforcingMod, only : LIS_forc + use LIS_logMod, only : LIS_logunit, LIS_endrun + use geis_forcingMod, only : geis_struc + use LIS_constantsMod, only : LIS_CONST_PATH_LEN + + implicit none +! !ARGUMENTS: + integer, intent(in) :: n + integer, intent(in) :: findex + +! +! !DESCRIPTION: +! Opens, reads, and interpolates hourly, GLBOAL EIS forcing. +! At the beginning of a simulation, the code reads the most recent +! past data (nearest hourly interval), and the nearest future data. +! These two datasets are used to temporally interpolate the data to +! the current model timestep. The strategy for missing data is to +! go backwards up to 10 days to get forcing at the same time of day. +! +!EOP + + integer :: c,f,ferror,try + integer :: order + real*8 :: time1,time2,timenow + real*8 :: dtime1, dtime2 + integer :: yr1,mo1,da1,hr1,mn1,ss1,doy1 + integer :: yr2,mo2,da2,hr2,mn2,ss2,doy2 + character(len=LIS_CONST_PATH_LEN) :: name + real :: gmt1,gmt2,ts1,ts2 + integer :: movetime ! 1=move time 2 data into time 1 + integer :: kk ! Forecast member index + +!=== End Variable Definition ============================================= + try=-999 + + kk = 1 + +!====Assumption will be not to find or move any data + geis_struc(n)%findtime1=0 + geis_struc(n)%findtime2=0 + movetime=0 + +!=== Determine Required NLDAS-2 Data Times (The previous hour and the future hour) + yr1=LIS_rc%yr + mo1=LIS_rc%mo + da1=LIS_rc%da + hr1=LIS_rc%hr + mn1=LIS_rc%mn + ss1=0 + ts1=0 + call LIS_tick(timenow,doy1,gmt1,yr1,mo1,da1,hr1,mn1,ss1,ts1) + + if(LIS_rc%ts.gt.3600) then + write(LIS_logunit,*) '[ERR] The model timestep is > forcing data timestep' + write(LIS_logunit,*) '[ERR] LIS does not support this mode currently' + write(LIS_logunit,*) '[ERR] Program stopping ...' + call LIS_endrun() + endif + + if(mod(nint(LIS_rc%ts),3600).eq.0) then + if(timenow.ge.geis_struc(n)%geistime2) then + yr1 = LIS_rc%yr + mo1=LIS_rc%mo + da1=LIS_rc%da + hr1=LIS_rc%hr + mn1=0 + ss1=0 + ts1=-60*60 + call LIS_tick(time1,doy1,gmt1,yr1,mo1,da1,hr1,mn1,ss1,ts1) + + yr2=LIS_rc%yr !next hour + mo2=LIS_rc%mo + da2=LIS_rc%da + hr2=LIS_rc%hr + mn2=0 + ss2=0 + ts2=0 + call LIS_tick(time2,doy2,gmt2,yr2,mo2,da2,hr2,mn2,ss2,ts2) + movetime = 1 + geis_struc(n)%findtime2 = 1 + endif + else + if(timenow.ge.geis_struc(n)%geistime2) then + yr1 = LIS_rc%yr + mo1=LIS_rc%mo + da1=LIS_rc%da + hr1=LIS_rc%hr + mn1=0 + ss1=0 + ts1=0 + call LIS_tick(time1,doy1,gmt1,yr1,mo1,da1,hr1,mn1,ss1,ts1) + + yr2=LIS_rc%yr !next hour + mo2=LIS_rc%mo + da2=LIS_rc%da + hr2=LIS_rc%hr + mn2=0 + ss2=0 + ts2=60*60 + call LIS_tick(time2,doy2,gmt2,yr2,mo2,da2,hr2,mn2,ss2,ts2) + + movetime = 1 + geis_struc(n)%findtime2 = 1 + endif + endif + + if(LIS_rc%tscount(n).eq.1 .or.LIS_rc%rstflag(n).eq.1 ) then !beginning of the run + geis_struc(n)%findtime1=1 + geis_struc(n)%findtime2=1 + movetime=0 + LIS_rc%rstflag(n) = 0 + endif + + if(movetime.eq.1) then + geis_struc(n)%geistime1=geis_struc(n)%geistime2 + do f=1,LIS_rc%met_nf(findex) + do c=1,LIS_rc%ngrid(n) + geis_struc(n)%metdata1(:,f,c)=geis_struc(n)%metdata2(:,f,c) + enddo + enddo + endif !end of movetime=1 + + if(geis_struc(n)%findtime1.eq.1) then +!=== the following looks back 10 days, at the same hour to fill data gaps. + ferror=0 + try=0 + ts1=-60*60*24 + + do + if ( ferror /= 0 ) exit + try=try+1 + + call get_geis_filename(n,kk,findex, & + name,geis_struc(n)%geisdir,& + yr1,mo1,da1,hr1) + + write(unit=LIS_logunit,fmt=*)'[INFO] getting file1.. ',trim(name) + order = 1 + call read_geis(n,kk,findex,order,mo1,name,ferror) + + + if(ferror.ge.1) geis_struc(n)%geistime1=time1 + call LIS_tick(dtime1,doy1,gmt1,yr1,mo1,da1,hr1,mn1,ss1,ts1) + if(try.gt.11)then + write(*,*)'error: GLOBAL EIS data gap exceeds 10 days on file 1' + stop + endif + enddo +!=== end of data search + endif !end of LIS_rc%findtime=1 + + + if(geis_struc(n)%findtime2.eq.1) then +!=== the following looks back 10 days, at the same hour to fill data gaps. + + ferror=0 + try=0 + ts2=-60*60*24 + do + if ( ferror /= 0 ) exit + try=try+1 + + call get_geis_filename(n,kk,findex,& + name,geis_struc(n)%geisdir,& + yr2,mo2,da2,hr2) + + write(unit=LIS_logunit,fmt=*)'[INFO] getting file2a.. ',trim(name) + order = 2 + + call read_geis(n,kk,findex,order,mo2,name,ferror) + + if(ferror.ge.1) then + geis_struc(n)%geistime2=time2 + endif + call LIS_tick(dtime2,doy2,gmt2,yr2,mo2,da2,hr2,mn2,ss2,ts2) + if(try.gt.11)then + write(*,*)'error: GLOBAL EIS data gap exceeds 10 days on file 2' + stop + endif + enddo + !=== end of data search + endif ! end of findtime2=1 + +end subroutine get_geis + + +!BOP +! !ROUTINE: get_geis_filename +! \label{get_geis_filename} +! +! !REVISION HISTORY: +! 15 Nov 2023: Sujay Kumar, initial specification +! +! !INTERFACE: + subroutine get_geis_filename(n,kk,findex,filename,odir,yr,mo,da,hr) + +! !USES: + use LIS_coreMod + use LIS_logMod + + implicit none +! !ARGUMENTS: + integer, intent(in) :: n + integer, intent(in) :: kk + integer, intent(in) :: findex + character(len=*), intent(out) :: filename + character(len=*), intent(in) :: odir + integer, intent(in) :: yr,mo,da,hr + +! !DESCRIPTION: +! This subroutine puts together file name for +! the Global EIS forcing data +! +! The arguments are: +! \begin{description} +! \item[odir] +! Name of the NLDAS-2 directory +! \item[yr] +! year +! \item[mo] +! month +! \item[da] +! day of month +! \item[hr] +! hour of day +! \item[filename] +! name of the timestamped Global EIS data +! \end{description} +! +!EOP + + character(len=6) :: fdir + character(len=10) :: ftime + + !=== end variable definition ============================================= + + !=== put together filename + + + write(UNIT=fdir,fmt='(i4,i2.2)') yr, mo + write(UNIT=ftime,fmt='(i4,i2.2,i2.2,i2.2)') yr, mo, da, hr + + filename = trim(odir) // '/' // fdir // '/MERRA2_CERES_IMERG_' // & + trim(ftime)//'00.d01.nc' + +end subroutine get_geis_filename diff --git a/lis/metforcing/geis/read_geis.F90 b/lis/metforcing/geis/read_geis.F90 new file mode 100644 index 000000000..086c57bc5 --- /dev/null +++ b/lis/metforcing/geis/read_geis.F90 @@ -0,0 +1,265 @@ +!-----------------------BEGIN NOTICE -- DO NOT EDIT----------------------- +! NASA Goddard Space Flight Center +! Land Information System Framework (LISF) +! Version 7.4 +! +! Copyright (c) 2022 United States Government as represented by the +! Administrator of the National Aeronautics and Space Administration. +! All Rights Reserved. +!-------------------------END NOTICE -- DO NOT EDIT----------------------- +#include "LIS_misc.h" +!BOP +! !ROUTINE: read_geis +! \label{read_geis} +! +! !REVISION HISTORY: +! 15 Nov 2023; Sujay Kumar, Initial Code +! +! !INTERFACE: +subroutine read_geis(n, kk, findex, order, month, name,ferror) +! !USES: + use LIS_coreMod + use LIS_logMod, only : LIS_logunit, LIS_verify, LIS_warning + use LIS_metforcingMod, only : LIS_forc + use geis_forcingMod, only : geis_struc +#if (defined USE_NETCDF3 || defined USE_NETCDF4) + use netcdf +#endif + + implicit none +! !ARGUMENTS: + integer, intent(in) :: n + integer, intent(in) :: kk ! Forecast member index + integer, intent(in) :: findex ! Forcing index + integer, intent(in) :: order + integer, intent(out) :: month + character(len=*), intent(in) :: name + integer, intent(out) :: ferror +! +! !DESCRIPTION: +! For the given time, reads parameters from +! the Global EIS data, transforms into 11 LIS forcing +! parameters +! +! The arguments are: +! \begin{description} +! \item[order] +! flag indicating which data to be read (order=1, read the previous +! hourly instance, order=2, read the next hourly instance) +! \item[n] +! index of the nest +! \item[name] +! name of the hourly GEIS forecast file +! \item[ferror] +! flag to indicate success of the call (=1 indicates success) +! \end{description} +! +! The routines invoked are: +! \begin{description} +! \item[interp\_geis](\ref{interp_geis}) \newline +! spatially interpolates a GEIS variable +! \end{description} +!EOP + + integer :: ftn + integer :: geis,varid + integer :: k,t,c,r,iret,rc + logical :: file_exists + real :: varfield(LIS_rc%lnc(n),LIS_rc%lnr(n)) + + ferror = 1 + +#if (defined USE_NETCDF4) + geis = (geis_struc(n)%ncold*geis_struc(n)%nrold) + + varfield = 0 + ferror = 1 + + inquire (file=name, exist=file_exists) + if (file_exists) then + + + call LIS_verify(nf90_open(path=trim(name), mode=NF90_NOWRITE, & + ncid=ftn), 'nf90_open failed in read_geis') + + call LIS_verify(nf90_inq_varid(ftn,'Tair',varId), & + 'nf90_inq_varid failed for Tair in read_geis') + call LIS_verify(nf90_get_var(ftn,varId, varfield,& + start=(/geis_struc(n)%c_off,& + geis_struc(n)%r_off/),& + count=(/LIS_rc%lnc(n),LIS_rc%lnr(n)/)),& + 'nf90_get_var failed for Tair in read_geis') + + do r=1,LIS_rc%lnr(n) + do c=1,LIS_rc%lnc(n) + if(LIS_domain(n)%gindex(c,r).ne.-1) then + if(order.eq.1) then + geis_struc(n)%metdata1(kk,1,& + LIS_domain(n)%gindex(c,r)) & + = varfield(c,r) + elseif(order.eq.2) then + geis_struc(n)%metdata2(kk,1,& + LIS_domain(n)%gindex(c,r))& + = varfield(c,r) + endif + endif + end do + enddo + + call LIS_verify(nf90_inq_varid(ftn,'Qair',varId), & + 'nf90_inq_varid failed for Qair in read_geis') + call LIS_verify(nf90_get_var(ftn,varId, varfield,& + start=(/geis_struc(n)%c_off,& + geis_struc(n)%r_off/),& + count=(/LIS_rc%lnc(n),LIS_rc%lnr(n)/)),& + 'nf90_get_var failed for Qair in read_geis') + + do r=1,LIS_rc%lnr(n) + do c=1,LIS_rc%lnc(n) + if(LIS_domain(n)%gindex(c,r).ne.-1) then + if(order.eq.1) then + geis_struc(n)%metdata1(kk,2,& + LIS_domain(n)%gindex(c,r)) & + = varfield(c,r) + elseif(order.eq.2) then + geis_struc(n)%metdata2(kk,2,& + LIS_domain(n)%gindex(c,r))& + = varfield(c,r) + endif + endif + end do + enddo + + + call LIS_verify(nf90_inq_varid(ftn,'SWdown',varId), & + 'nf90_inq_varid failed for SWdown in read_geis') + call LIS_verify(nf90_get_var(ftn,varId, varfield,& + start=(/geis_struc(n)%c_off,& + geis_struc(n)%r_off/),& + count=(/LIS_rc%lnc(n),LIS_rc%lnr(n)/)),& + 'nf90_get_var failed for SWdown in read_geis') + + do r=1,LIS_rc%lnr(n) + do c=1,LIS_rc%lnc(n) + if(LIS_domain(n)%gindex(c,r).ne.-1) then + if(order.eq.1) then + geis_struc(n)%metdata1(kk,3,& + LIS_domain(n)%gindex(c,r)) & + = varfield(c,r) + elseif(order.eq.2) then + geis_struc(n)%metdata2(kk,3,& + LIS_domain(n)%gindex(c,r))& + = varfield(c,r) + endif + endif + end do + enddo + + + call LIS_verify(nf90_inq_varid(ftn,'LWdown',varId), & + 'nf90_inq_varid failed for LWdown in read_geis') + call LIS_verify(nf90_get_var(ftn,varId, varfield,& + start=(/geis_struc(n)%c_off,& + geis_struc(n)%r_off/),& + count=(/LIS_rc%lnc(n),LIS_rc%lnr(n)/)),& + 'nf90_get_var failed for LWdown in read_geis') + + do r=1,LIS_rc%lnr(n) + do c=1,LIS_rc%lnc(n) + if(LIS_domain(n)%gindex(c,r).ne.-1) then + if(order.eq.1) then + geis_struc(n)%metdata1(kk,4,& + LIS_domain(n)%gindex(c,r)) & + = varfield(c,r) + elseif(order.eq.2) then + geis_struc(n)%metdata2(kk,4,& + LIS_domain(n)%gindex(c,r))& + = varfield(c,r) + endif + endif + end do + enddo + + + call LIS_verify(nf90_inq_varid(ftn,'Wind',varId), & + 'nf90_inq_varid failed for Wind in read_geis') + call LIS_verify(nf90_get_var(ftn,varId, varfield,& + start=(/geis_struc(n)%c_off,& + geis_struc(n)%r_off/),& + count=(/LIS_rc%lnc(n),LIS_rc%lnr(n)/)),& + 'nf90_get_var failed for Wind in read_geis') + + do r=1,LIS_rc%lnr(n) + do c=1,LIS_rc%lnc(n) + if(LIS_domain(n)%gindex(c,r).ne.-1) then + if(order.eq.1) then + geis_struc(n)%metdata1(kk,5,& + LIS_domain(n)%gindex(c,r)) & + = varfield(c,r) + elseif(order.eq.2) then + geis_struc(n)%metdata2(kk,5,& + LIS_domain(n)%gindex(c,r))& + = varfield(c,r) + endif + endif + end do + enddo + + call LIS_verify(nf90_inq_varid(ftn,'Psurf',varId), & + 'nf90_inq_varid failed for Psurf in read_geis') + call LIS_verify(nf90_get_var(ftn,varId, varfield,& + start=(/geis_struc(n)%c_off,& + geis_struc(n)%r_off/),& + count=(/LIS_rc%lnc(n),LIS_rc%lnr(n)/)),& + 'nf90_get_var failed for Psurf in read_geis') + + do r=1,LIS_rc%lnr(n) + do c=1,LIS_rc%lnc(n) + if(LIS_domain(n)%gindex(c,r).ne.-1) then + if(order.eq.1) then + geis_struc(n)%metdata1(kk,6,& + LIS_domain(n)%gindex(c,r)) & + = varfield(c,r) + elseif(order.eq.2) then + geis_struc(n)%metdata2(kk,6,& + LIS_domain(n)%gindex(c,r))& + = varfield(c,r) + endif + endif + end do + enddo + + call LIS_verify(nf90_inq_varid(ftn,'Rainf',varId), & + 'nf90_inq_varid failed for Rainf in read_geis') + call LIS_verify(nf90_get_var(ftn,varId, varfield,& + start=(/geis_struc(n)%c_off,& + geis_struc(n)%r_off/),& + count=(/LIS_rc%lnc(n),LIS_rc%lnr(n)/)),& + 'nf90_get_var failed for Rainf in read_geis') + + do r=1,LIS_rc%lnr(n) + do c=1,LIS_rc%lnc(n) + if(LIS_domain(n)%gindex(c,r).ne.-1) then + if(order.eq.1) then + geis_struc(n)%metdata1(kk,7,& + LIS_domain(n)%gindex(c,r)) & + = varfield(c,r) + elseif(order.eq.2) then + geis_struc(n)%metdata2(kk,7,& + LIS_domain(n)%gindex(c,r))& + = varfield(c,r) + endif + endif + end do + enddo + + + else + write(LIS_logunit,*) & + '[ERR] Could not find file: ',trim(name) + ferror = 0 + endif +#endif + +end subroutine read_geis + diff --git a/lis/metforcing/geis/readcrd_geis.F90 b/lis/metforcing/geis/readcrd_geis.F90 new file mode 100644 index 000000000..6fcb14c19 --- /dev/null +++ b/lis/metforcing/geis/readcrd_geis.F90 @@ -0,0 +1,51 @@ +!-----------------------BEGIN NOTICE -- DO NOT EDIT----------------------- +! NASA Goddard Space Flight Center +! Land Information System Framework (LISF) +! Version 7.4 +! +! Copyright (c) 2022 United States Government as represented by the +! Administrator of the National Aeronautics and Space Administration. +! All Rights Reserved. +!-------------------------END NOTICE -- DO NOT EDIT----------------------- +!BOP +! +! !ROUTINE: readcrd_geis +! \label{readcrd_geis} +! +! !REVISION HISTORY: +! 15 Nov 2023; Sujay Kumar, Initial Code +! +! !INTERFACE: +subroutine readcrd_geis() +! !USES: + use ESMF + use LIS_coreMod, only : LIS_rc, LIS_config + use LIS_logMod, only : LIS_logunit, LIS_verify, LIS_endrun + use geis_forcingMod, only : geis_struc + +! !DESCRIPTION: +! +! This routine reads the options specific to the global EIS forcing from +! the LIS configuration file. +! +!EOP + + implicit none + integer :: n,rc + + call ESMF_ConfigFindLabel(LIS_config,"Global EIS forcing directory:",rc=rc) + call LIS_verify(rc, 'Global EIS forcing directory: not defined') + do n=1,LIS_rc%nnest + call ESMF_ConfigGetAttribute(LIS_config,geis_struc(n)%geisdir,rc=rc) + enddo + + write(unit=LIS_logunit,fmt=*)'[INFO] Using Global EIS forcing' + + do n=1,LIS_rc%nnest + write(unit=LIS_logunit,fmt=*) '[INFO] Global EIS forcing directory : ',trim(geis_struc(n)%geisdir) + + geis_struc(n)%geistime1 = 3000.0 + geis_struc(n)%geistime2 = 0.0 + enddo + +end subroutine readcrd_geis diff --git a/lis/metforcing/geis/reset_geis.F90 b/lis/metforcing/geis/reset_geis.F90 new file mode 100644 index 000000000..ca45f9040 --- /dev/null +++ b/lis/metforcing/geis/reset_geis.F90 @@ -0,0 +1,38 @@ +!-----------------------BEGIN NOTICE -- DO NOT EDIT----------------------- +! NASA Goddard Space Flight Center +! Land Information System Framework (LISF) +! Version 7.4 +! +! Copyright (c) 2022 United States Government as represented by the +! Administrator of the National Aeronautics and Space Administration. +! All Rights Reserved. +!-------------------------END NOTICE -- DO NOT EDIT----------------------- +#include "LIS_misc.h" +!BOP +! !MODULE: reset_geis +! \label{reset_geis} +! +! !REVISION HISTORY: +! 15 Nov 2023; Sujay Kumar, Initial Code +! +! !INTERFACE: +subroutine reset_geis() +! !USES: + use LIS_coreMod, only : LIS_rc + use geis_forcingMod, only : geis_struc +! +! !DESCRIPTION: +! Routine to reset geis forcing related memory allocations. +! +!EOP + implicit none + + integer :: n + integer :: findex + + do n=1,LIS_rc%nnest + geis_struc(n)%geistime1 = 3000.0 + geis_struc(n)%geistime2 = 0.0 + enddo + +end subroutine reset_geis diff --git a/lis/metforcing/geis/timeinterp_geis.F90 b/lis/metforcing/geis/timeinterp_geis.F90 new file mode 100644 index 000000000..d0c05a928 --- /dev/null +++ b/lis/metforcing/geis/timeinterp_geis.F90 @@ -0,0 +1,306 @@ +!-----------------------BEGIN NOTICE -- DO NOT EDIT----------------------- +! NASA Goddard Space Flight Center +! Land Information System Framework (LISF) +! Version 7.4 +! +! Copyright (c) 2022 United States Government as represented by the +! Administrator of the National Aeronautics and Space Administration. +! All Rights Reserved. +!-------------------------END NOTICE -- DO NOT EDIT----------------------- +#include "LIS_misc.h" +!BOP +! !ROUTINE: timeinterp_geis +! \label{timeinterp_geis} +! +! !REVISION HISTORY: +! +! 15 Nov 2023: Sujay Kumar; Initial Version in LIS! +! !INTERFACE: +subroutine timeinterp_geis(n,findex) +! !USES: + use ESMF + use LIS_coreMod, only : LIS_rc,LIS_domain + use LIS_constantsMod, only : LIS_CONST_SOLAR + use LIS_metforcingMod, only : LIS_forc, LIS_FORC_Base_State + use LIS_FORC_AttributesMod + use LIS_timeMgrMod, only : LIS_tick, LIS_time2date + use LIS_logMod, only : LIS_logunit, LIS_verify + use geis_forcingMod, only : geis_struc + use LIS_forecastMod, only : LIS_get_iteration_index + + implicit none +! !ARGUMENTS: + integer, intent(in) :: n + integer, intent(in) :: findex +! +! !DESCRIPTION: +! Temporally interpolates the forcing data to the current model +! timestep. Downward shortwave radiation is interpolated using a +! zenith-angled based approach. Precipitation and longwave radiation +! are not temporally interpolated, and the previous 3 hourly value +! is used. All other variables are linearly interpolated between +! the 3 hourly blocks. +! +! The routines invoked are: +! \begin{description} +! \item[LIS\_time2date](\ref{LIS_time2date}) \newline +! converts the time to a date format +! \item[LIS\_tick](\ref{LIS_tick}) \newline +! advances or retracts time by the specified amount +! \item[zterp](\ref{zterp}) \newline +! zenith-angle based interpolation +! \end{description} +!EOP + integer :: zdoy + real :: zw1, zw2 + real :: czm, cze, czb + real :: wt1, wt2,swt1,swt2 + real :: gmt1, gmt2, tempbts + integer :: t,index1 + integer :: bdoy,byr,bmo,bda,bhr,bmn + real*8 :: btime,newtime1,newtime2 + real :: tempgmt1,tempgmt2 + integer :: tempbdoy,tempbyr,tempbmo,tempbda,tempbhr,tempbmn + integer :: tempbss + integer :: status + type(ESMF_Field) :: tmpField,q2Field,uField,vField,swdField,lwdField + type(ESMF_Field) :: psurfField,pcpField,cpcpField,fhgtField,acondField + type(ESMF_Field) :: PETField,CAPEField + real,pointer :: tmp(:),q2(:),uwind(:),vwind(:) + real,pointer :: swd(:),lwd(:),psurf(:),pcp(:),cpcp(:) + real,pointer :: fheight(:),acond(:),pet(:),cape(:) + logical :: forcing_z, forcing_ch, forcing_pet, forcing_cape + integer :: mfactor, m, k, kk +! ________________________________________ + + btime=geis_struc(n)%geistime1 + call LIS_time2date(btime,bdoy,gmt1,byr,bmo,bda,bhr,bmn) + + tempbdoy=bdoy + tempgmt1=gmt1 + tempbyr=byr + tempbmo=bmo + tempbda=bda + tempbhr=bhr + if (tempbhr.eq.24) tempbhr=0 + tempbmn=bmn + tempbss=0 + tempbts=0 + call LIS_tick(newtime1,tempbdoy,tempgmt1,& + tempbyr,tempbmo,tempbda,tempbhr,tempbmn, & + tempbss,tempbts) + + btime=geis_struc(n)%geistime2 + call LIS_time2date(btime,bdoy,gmt2,byr,bmo,bda,bhr,bmn) + tempbdoy=bdoy + tempgmt2=gmt2 + tempbyr=byr + tempbmo=bmo + tempbda=bda + tempbhr=bhr + if (tempbhr.eq.24) tempbhr=0 + tempbmn=bmn + tempbss=0 + tempbts=0 + call LIS_tick(newtime2,tempbdoy,tempgmt2,& + tempbyr,tempbmo,tempbda,tempbhr,tempbmn,& + tempbss,tempbts) + +!=== Interpolate Data in time + wt1=(geis_struc(n)%geistime2-LIS_rc%time)/ & + (geis_struc(n)%geistime2-geis_struc(n)%geistime1) + wt2=1.0-wt1 + swt1=(newtime2-LIS_rc%time)/(newtime2-newtime1) + swt2=1.0-swt1 + + + call ESMF_StateGet(LIS_FORC_Base_State(n,findex),LIS_FORC_Tair%varname(1),tmpField,& + rc=status) + call LIS_verify(status, 'Error: Enable Tair in the forcing variables list') + + call ESMF_StateGet(LIS_FORC_Base_State(n,findex),LIS_FORC_Qair%varname(1),q2Field,& + rc=status) + call LIS_verify(status, 'Error: Enable Qair in the forcing variables list') + + call ESMF_StateGet(LIS_FORC_Base_State(n,findex),LIS_FORC_SWdown%varname(1),swdField,& + rc=status) + call LIS_verify(status, 'Error: Enable SWdown in the forcing variables list') + + call ESMF_StateGet(LIS_FORC_Base_State(n,findex),LIS_FORC_LWdown%varname(1),lwdField,& + rc=status) + call LIS_verify(status, 'Error: Enable LWdown in the forcing variables list') + + call ESMF_StateGet(LIS_FORC_Base_State(n,findex),LIS_FORC_Wind_E%varname(1),uField,& + rc=status) + call LIS_verify(status, 'Error: Enable Wind_E in the forcing variables list') + + call ESMF_StateGet(LIS_FORC_Base_State(n,findex),LIS_FORC_Wind_N%varname(1),vField,& + rc=status) + call LIS_verify(status, 'Error: Enable Wind_N in the forcing variables list') + + call ESMF_StateGet(LIS_FORC_Base_State(n,findex),LIS_FORC_Psurf%varname(1),psurfField,& + rc=status) + call LIS_verify(status, 'Error: Enable Psurf in the forcing variables list') + + call ESMF_StateGet(LIS_FORC_Base_State(n,findex),LIS_FORC_Rainf%varname(1),pcpField,& + rc=status) + call LIS_verify(status, 'Error: Enable Rainf in the forcing variables list') + + + call ESMF_FieldGet(swdField,localDE=0,farrayPtr=swd,rc=status) + call LIS_verify(status) + +! Loop over number of forcing ensembles: + mfactor = LIS_rc%nensem(n)/geis_struc(n)%nIter + + do k=1,LIS_rc%ntiles(n)/mfactor + do m=1,mfactor + t = m + (k-1)*mfactor + index1 = LIS_domain(n)%tile(t)%index + zdoy=LIS_rc%doy + ! compute and apply zenith angle weights + call zterp(1,LIS_domain(n)%grid(index1)%lat,& + LIS_domain(n)%grid(index1)%lon,& + gmt1,gmt2,LIS_rc%gmt,zdoy,zw1,zw2,czb,cze,czm,LIS_rc) + + kk = LIS_get_iteration_index(n, k, index1, mfactor) + + if(geis_struc(n)%metdata1(kk,3,index1).ne.LIS_rc%udef.and.& + geis_struc(n)%metdata2(kk,3,index1).ne.LIS_rc%udef) then + swd(t) = geis_struc(n)%metdata1(kk,3,index1)*zw1+& + geis_struc(n)%metdata2(kk,3,index1)*zw2 + ! In cases of small cos(zenith) angles, use linear weighting + ! to avoid overly large weights + + if((swd(t).gt.geis_struc(n)%metdata1(kk,3,index1).and. & + swd(t).gt.geis_struc(n)%metdata2(kk,3,index1)).and. & + (czb.lt.0.1.or.cze.lt.0.1))then + swd(t) = geis_struc(n)%metdata1(kk,3,index1)*swt1+ & + geis_struc(n)%metdata2(kk,3,index1)*swt2 + endif + endif + + if (swd(t).gt.LIS_CONST_SOLAR) then + write(unit=LIS_logunit,fmt=*)'[WARN] sw radiation too high!!' + write(unit=LIS_logunit,fmt=*)'[WARN] it is', swd(t) + write(unit=LIS_logunit,fmt=*)'[WARN] data1=',geis_struc(n)%metdata1(kk,3,index1) + write(unit=LIS_logunit,fmt=*)'[WARN] data2=',geis_struc(n)%metdata2(kk,3,index1) + write(unit=LIS_logunit,fmt=*)'[WARN] zw1=',zw1,'zw2=',zw2 + write(unit=LIS_logunit,fmt=*)'[WARN] swt1=',swt1,'swt2=',swt2 + endif + enddo + enddo ! End for SWdown + + ! do block precipitation interpolation + call ESMF_FieldGet(pcpField,localDE=0,farrayPtr=pcp,rc=status) + call LIS_verify(status) + + do k=1,LIS_rc%ntiles(n)/mfactor + do m=1,mfactor + t = m + (k-1)*mfactor + index1 = LIS_domain(n)%tile(t)%index + kk = LIS_get_iteration_index(n, k, index1, mfactor) + if(geis_struc(n)%metdata2(kk,7,index1).ne.LIS_rc%udef) then + pcp(t)=geis_struc(n)%metdata2(kk,7,index1) + pcp(t) = pcp(t)/(60.0*60.0) + endif + end do + enddo + + + !linearly interpolate everything else + + call ESMF_FieldGet(tmpField,localDE=0,farrayPtr=tmp,rc=status) + call LIS_verify(status) + + do k=1,LIS_rc%ntiles(n)/mfactor + do m=1,mfactor + t = m + (k-1)*mfactor + index1 = LIS_domain(n)%tile(t)%index + kk = LIS_get_iteration_index(n, k, index1, mfactor) + if((geis_struc(n)%metdata1(kk,1,index1).ne.LIS_rc%udef).and.& + (geis_struc(n)%metdata2(kk,1,index1).ne.LIS_rc%udef)) then + tmp(t) = geis_struc(n)%metdata1(kk,1,index1)*wt1+ & + geis_struc(n)%metdata2(kk,1,index1)*wt2 + endif + enddo + enddo + + call ESMF_FieldGet(q2Field,localDE=0,farrayPtr=q2,rc=status) + call LIS_verify(status) + + do k=1,LIS_rc%ntiles(n)/mfactor + do m=1,mfactor + t = m + (k-1)*mfactor + index1 = LIS_domain(n)%tile(t)%index + kk = LIS_get_iteration_index(n, k, index1, mfactor) + if((geis_struc(n)%metdata1(kk,2,index1).ne.LIS_rc%udef).and.& + (geis_struc(n)%metdata2(kk,2,index1).ne.LIS_rc%udef)) then + q2(t) =geis_struc(n)%metdata1(kk,2,index1)*wt1+ & + geis_struc(n)%metdata2(kk,2,index1)*wt2 + endif + enddo + enddo + + call ESMF_FieldGet(lwdField,localDE=0,farrayPtr=lwd,rc=status) + call LIS_verify(status) + + do k=1,LIS_rc%ntiles(n)/mfactor + do m=1,mfactor + t = m + (k-1)*mfactor + index1 = LIS_domain(n)%tile(t)%index + kk = LIS_get_iteration_index(n, k, index1, mfactor) + if((geis_struc(n)%metdata1(kk,4,index1).ne.LIS_rc%udef).and.& + (geis_struc(n)%metdata2(kk,4,index1).ne.LIS_rc%udef)) then + lwd(t) =geis_struc(n)%metdata1(kk,4,index1)*wt1+ & + geis_struc(n)%metdata2(kk,4,index1)*wt2 + endif + enddo + enddo + + call ESMF_FieldGet(uField,localDE=0,farrayPtr=uwind,rc=status) + call LIS_verify(status) + + do k=1,LIS_rc%ntiles(n)/mfactor + do m=1,mfactor + t = m + (k-1)*mfactor + index1 = LIS_domain(n)%tile(t)%index + kk = LIS_get_iteration_index(n, k, index1, mfactor) + if((geis_struc(n)%metdata1(kk,5,index1).ne.LIS_rc%udef).and.& + (geis_struc(n)%metdata2(kk,5,index1).ne.LIS_rc%udef)) then + uwind(t) = geis_struc(n)%metdata1(kk,5,index1)*wt1+ & + geis_struc(n)%metdata2(kk,5,index1)*wt2 + endif + enddo + enddo + + call ESMF_FieldGet(vField,localDE=0,farrayPtr=vwind,rc=status) + call LIS_verify(status) + + do k=1,LIS_rc%ntiles(n)/mfactor + do m=1,mfactor + t = m + (k-1)*mfactor + index1 = LIS_domain(n)%tile(t)%index + + vwind(t) =0.0 + enddo + enddo + + call ESMF_FieldGet(psurfField,localDE=0,farrayPtr=psurf,rc=status) + call LIS_verify(status) + + do k=1,LIS_rc%ntiles(n)/mfactor + do m=1,mfactor + t = m + (k-1)*mfactor + index1 = LIS_domain(n)%tile(t)%index + kk = LIS_get_iteration_index(n, k, index1, mfactor) + if((geis_struc(n)%metdata1(kk,6,index1).ne.LIS_rc%udef).and.& + (geis_struc(n)%metdata2(kk,6,index1).ne.LIS_rc%udef)) then + psurf(t) =geis_struc(n)%metdata1(kk,6,index1)*wt1+ & + geis_struc(n)%metdata2(kk,6,index1)*wt2 + endif + enddo + enddo + +end subroutine timeinterp_geis + diff --git a/lis/plugins/LIS_metforcing_pluginMod.F90 b/lis/plugins/LIS_metforcing_pluginMod.F90 index a4f06a4cd..97c2c5ec7 100644 --- a/lis/plugins/LIS_metforcing_pluginMod.F90 +++ b/lis/plugins/LIS_metforcing_pluginMod.F90 @@ -208,6 +208,10 @@ subroutine LIS_metforcing_plugin use nldas2_forcingMod #endif +#if ( defined MF_GEIS ) + use geis_forcingMod +#endif + #if ( defined MF_NARR ) use narr_forcingMod #endif @@ -500,6 +504,13 @@ subroutine LIS_metforcing_plugin external reset_nldas2 #endif +#if ( defined MF_GEIS ) + external get_geis + external timeinterp_geis + external finalize_geis + external reset_geis +#endif + #if ( defined MF_NARR ) external get_narr external timeinterp_narr @@ -964,6 +975,16 @@ subroutine LIS_metforcing_plugin call registerresetmetforc(trim(LIS_nldas2Id)//char(0),reset_nldas2) #endif +#if ( defined MF_GEIS ) +! - GEIS Forcing: + call registerinitmetforc(trim(LIS_geisId)//char(0),init_geis) + call registerretrievemetforc(trim(LIS_geisId)//char(0),get_geis) + call registertimeinterpmetforc(trim(LIS_geisId)//char(0), & + timeinterp_geis) + call registerfinalmetforc(trim(LIS_geisId)//char(0),finalize_geis) + call registerresetmetforc(trim(LIS_geisId)//char(0),reset_geis) +#endif + #if ( defined MF_NARR ) ! - NARR profile data for CRTM call registerinitmetforc(trim(LIS_narrId)//char(0),init_NARR) diff --git a/lis/plugins/LIS_pluginIndices.F90 b/lis/plugins/LIS_pluginIndices.F90 index 76536bbc3..3c166e7ec 100644 --- a/lis/plugins/LIS_pluginIndices.F90 +++ b/lis/plugins/LIS_pluginIndices.F90 @@ -131,6 +131,7 @@ module LIS_pluginIndices character*50, public, parameter :: LIS_agrmetId = "AGRMET" character*50, public, parameter :: LIS_princetonId = "PRINCETON" character*50, public, parameter :: LIS_nldas2Id = "NLDAS2" + character*50, public, parameter :: LIS_geisId = "Global EIS" character*50, public, parameter :: LIS_gldasId = "GLDAS" character*50, public, parameter :: LIS_gfsId = "GFS"