From 82f8083c5e9561254a951c20bb9f65fef639a429 Mon Sep 17 00:00:00 2001 From: Eric Kemp Date: Thu, 1 Feb 2024 15:05:19 -0500 Subject: [PATCH] Some style updates to make the code easier to read. Identical results are given to older GNU and Intel compiler runs. --- lis/metforcing/mogreps_g/get_cdf_params.F90 | 144 ++++---- lis/metforcing/mogreps_g/get_mogrepsg.F90 | 202 +++++++----- .../mogreps_g/mogrepsg_forcingMod.F90 | 58 ++-- lis/metforcing/mogreps_g/readcrd_mogrepsg.F90 | 75 +++-- .../mogreps_g/timeinterp_mogrepsg.F90 | 308 ++++++++++-------- 5 files changed, 438 insertions(+), 349 deletions(-) diff --git a/lis/metforcing/mogreps_g/get_cdf_params.F90 b/lis/metforcing/mogreps_g/get_cdf_params.F90 index b63fb041e..5d9c0fc07 100644 --- a/lis/metforcing/mogreps_g/get_cdf_params.F90 +++ b/lis/metforcing/mogreps_g/get_cdf_params.F90 @@ -28,82 +28,85 @@ subroutine get_cdf_params (n, fname, month, param_a, param_b, mean, std) implicit none - !ARGUMENTS: + !ARGUMENTS: integer, intent(in) :: n character(len=*), intent(in) :: fname integer, intent(in) :: month ! month of year (1-12) - !EOP - logical :: file_exists - integer :: ftn - integer :: param_a_id, param_b_id, mean_id, std_id, ngrid_id, nlead_id - integer :: ngrid, nlead - integer :: r, c, nc, c1, r1, gid, l - real,allocatable :: param_hires(:,:) - real,allocatable :: param_hires_2d(:,:) - real :: param_a(LIS_rc%ngrid(n),8) !8: lead time - real :: param_b(LIS_rc%ngrid(n),8) - real :: mean(LIS_rc%ngrid(n),8) - real :: std(LIS_rc%ngrid(n),8) - - ! check file + !EOP + logical :: file_exists + integer :: ftn + integer :: param_a_id, param_b_id, mean_id, & + std_id, ngrid_id, nlead_id + integer :: ngrid, nlead + integer :: r, c, nc, c1, r1, gid, l + real, allocatable :: param_hires(:,:) + real, allocatable :: param_hires_2d(:,:) + real :: param_a(LIS_rc%ngrid(n),8) !8:lead time + real :: param_b(LIS_rc%ngrid(n),8) + real :: mean(LIS_rc%ngrid(n),8) + real :: std(LIS_rc%ngrid(n),8) + + ! check file inquire(file=fname, exist=file_exists) - if(.not. file_exists) then - write(LIS_logunit,*) '[ERR] ',trim(fname)//' does not exist' + if (.not. file_exists) then + write(LIS_logunit,*) '[ERR] ', trim(fname) // ' does not exist' call LIS_endrun() endif - write(LIS_logunit,*)'[INFO] Getting MOGREPS-G bias correction parameters ', trim(fname) + write(LIS_logunit,*) & + '[INFO] Getting MOGREPS-G bias correction parameters ', trim(fname) call LIS_verify(nf90_open(path=trim(fname), mode=NF90_NOWRITE, & ncid=ftn), 'nf90_open failed in get_cdf_params') call LIS_verify(nf90_inq_dimid(ftn, "ngrid", ngrid_id), & - 'nf90_inq_dimid failed for ngrid in get_cdf_params') + 'nf90_inq_dimid failed for ngrid in get_cdf_params') call LIS_verify(nf90_inquire_dimension(ftn, ngrid_id, len=ngrid),& - 'nf90_inquire_dimension failed for ngrid in get_cdf_params') + 'nf90_inquire_dimension failed for ngrid in get_cdf_params') - if(LIS_rc%gnc(n)*LIS_rc%gnr(n) .ne. ngrid) then - write(LIS_logunit,*) 'The input dimensions of the '//trim(fname) + if (LIS_rc%gnc(n)*LIS_rc%gnr(n) .ne. ngrid) then + write(LIS_logunit,*) '[ERR] The input dimensions of the '//trim(fname) write(LIS_logunit,*) '(',ngrid,')' - write(LIS_logunit,*) 'does not match the dimensions in the LIS parameter file' + write(LIS_logunit,*) & + 'does not match the dimensions in the LIS parameter file' write(LIS_logunit,*) '(',LIS_rc%gnc(n)*LIS_rc%gnr(n),')' call LIS_endrun() endif call LIS_verify(nf90_inq_dimid(ftn, "lead_time", nlead_id), & - 'nf90_inq_dimid failed for nlead in get_cdf_params') + 'nf90_inq_dimid failed for nlead in get_cdf_params') call LIS_verify(nf90_inquire_dimension(ftn, nlead_id, len=nlead),& - 'nf90_inquire_dimension failed for nlead in get_cdf_params') + 'nf90_inquire_dimension failed for nlead in get_cdf_params') - allocate(param_hires(LIS_rc%gnc(n)*LIS_rc%gnr(n),nlead)) + allocate(param_hires(LIS_rc%gnc(n)*LIS_rc%gnr(n),nlead)) allocate(param_hires_2d(LIS_rc%gnc(n),LIS_rc%gnr(n))) param_hires = -9999.0 param_hires_2d = -9999.0 - + ! read param_a - call LIS_verify(nf90_inq_varid(ftn,'cdf_param_a',param_a_id), & + call LIS_verify(nf90_inq_varid(ftn, 'cdf_param_a', param_a_id), & 'nf90_inq_varid failed for cdf_param_a in get_cdf_params') call LIS_verify(nf90_get_var(ftn, param_a_id, param_hires, & - start=(/1,1,month/), count=(/ngrid,nlead,1/)),& + start=(/1,1,month/), count=(/ngrid,nlead,1/)), & 'nf90_get_var failed for cdf_param_a in get_cdf_params') - - do l=1,nlead + + do l = 1, nlead ! 1D -> 2D - do r=1,LIS_rc%gnr(n) - do c=1,LIS_rc%gnc(n) - param_hires_2d(c,r)=param_hires(c+(r-1)*LIS_rc%gnc(n),l) + do r = 1, LIS_rc%gnr(n) + do c = 1, LIS_rc%gnc(n) + param_hires_2d(c,r) = param_hires(c+(r-1)*LIS_rc%gnc(n),l) enddo enddo !subsets the data for each processor's domain nc = (LIS_ewe_halo_ind(n,LIS_localPet+1)-LIS_ews_halo_ind(n,LIS_localPet+1))+1 - do r=LIS_nss_halo_ind(n,LIS_localPet+1),LIS_nse_halo_ind(n,LIS_localPet+1) - do c=LIS_ews_halo_ind(n,LIS_localPet+1),LIS_ewe_halo_ind(n,LIS_localPet+1) - c1 = c-LIS_ews_halo_ind(n,LIS_localPet+1)+1 - r1 = r-LIS_nss_halo_ind(n,LIS_localPet+1)+1 + do r = LIS_nss_halo_ind(n,LIS_localPet+1),LIS_nse_halo_ind(n,LIS_localPet+1) + do c = LIS_ews_halo_ind(n,LIS_localPet+1),LIS_ewe_halo_ind(n,LIS_localPet+1) + c1 = c - LIS_ews_halo_ind(n,LIS_localPet+1)+1 + r1 = r - LIS_nss_halo_ind(n,LIS_localPet+1)+1 gid = LIS_domain(n)%gindex(c1,r1) - if(gid.ne.-1) then + if (gid .ne. -1) then param_a(gid,l) = param_hires_2d(c,r) endif enddo @@ -111,28 +114,28 @@ subroutine get_cdf_params (n, fname, month, param_a, param_b, mean, std) enddo ! read param_b - call LIS_verify(nf90_inq_varid(ftn,'cdf_param_b',param_b_id), & + call LIS_verify(nf90_inq_varid(ftn, 'cdf_param_b', param_b_id), & 'nf90_inq_varid failed for cdf_param_b in get_cdf_params') call LIS_verify(nf90_get_var(ftn, param_b_id, param_hires, & start=(/1,1,month/), count=(/ngrid,nlead,1/)),& 'nf90_get_var failed for cdf_param_b in get_cdf_params') - do l=1,nlead + do l = 1, nlead ! 1D -> 2D - do r=1,LIS_rc%gnr(n) - do c=1,LIS_rc%gnc(n) - param_hires_2d(c,r)=param_hires(c+(r-1)*LIS_rc%gnc(n),l) + do r = 1, LIS_rc%gnr(n) + do c = 1,LIS_rc%gnc(n) + param_hires_2d(c,r) = param_hires(c+(r-1)*LIS_rc%gnc(n),l) enddo enddo !subsets the data for each processor's domain nc = (LIS_ewe_halo_ind(n,LIS_localPet+1)-LIS_ews_halo_ind(n,LIS_localPet+1))+1 - do r=LIS_nss_halo_ind(n,LIS_localPet+1),LIS_nse_halo_ind(n,LIS_localPet+1) - do c=LIS_ews_halo_ind(n,LIS_localPet+1),LIS_ewe_halo_ind(n,LIS_localPet+1) - c1 = c-LIS_ews_halo_ind(n,LIS_localPet+1)+1 - r1 = r-LIS_nss_halo_ind(n,LIS_localPet+1)+1 + do r = LIS_nss_halo_ind(n,LIS_localPet+1), LIS_nse_halo_ind(n,LIS_localPet+1) + do c = LIS_ews_halo_ind(n,LIS_localPet+1), LIS_ewe_halo_ind(n,LIS_localPet+1) + c1 = c - LIS_ews_halo_ind(n,LIS_localPet+1) + 1 + r1 = r - LIS_nss_halo_ind(n,LIS_localPet+1) + 1 gid = LIS_domain(n)%gindex(c1,r1) - if(gid.ne.-1) then + if (gid.ne.-1) then param_b(gid,l) = param_hires_2d(c,r) endif enddo @@ -140,28 +143,28 @@ subroutine get_cdf_params (n, fname, month, param_a, param_b, mean, std) enddo ! read mean - call LIS_verify(nf90_inq_varid(ftn,'mean',mean_id), & + call LIS_verify(nf90_inq_varid(ftn, 'mean', mean_id), & 'nf90_inq_varid failed for mean in get_cdf_params') call LIS_verify(nf90_get_var(ftn, mean_id, param_hires, & start=(/1,1,month/), count=(/ngrid,nlead,1/)),& 'nf90_get_var failed for mean in get_cdf_params') - do l=1,nlead + do l = 1, nlead ! 1D -> 2D - do r=1,LIS_rc%gnr(n) - do c=1,LIS_rc%gnc(n) - param_hires_2d(c,r)=param_hires(c+(r-1)*LIS_rc%gnc(n),l) + do r = 1, LIS_rc%gnr(n) + do c = 1, LIS_rc%gnc(n) + param_hires_2d(c,r) = param_hires(c+(r-1)*LIS_rc%gnc(n),l) enddo enddo !subsets the data for each processor's domain nc = (LIS_ewe_halo_ind(n,LIS_localPet+1)-LIS_ews_halo_ind(n,LIS_localPet+1))+1 - do r=LIS_nss_halo_ind(n,LIS_localPet+1),LIS_nse_halo_ind(n,LIS_localPet+1) - do c=LIS_ews_halo_ind(n,LIS_localPet+1),LIS_ewe_halo_ind(n,LIS_localPet+1) - c1 = c-LIS_ews_halo_ind(n,LIS_localPet+1)+1 - r1 = r-LIS_nss_halo_ind(n,LIS_localPet+1)+1 + do r = LIS_nss_halo_ind(n,LIS_localPet+1), LIS_nse_halo_ind(n,LIS_localPet+1) + do c = LIS_ews_halo_ind(n,LIS_localPet+1), LIS_ewe_halo_ind(n,LIS_localPet+1) + c1 = c - LIS_ews_halo_ind(n,LIS_localPet+1)+1 + r1 = r - LIS_nss_halo_ind(n,LIS_localPet+1)+1 gid = LIS_domain(n)%gindex(c1,r1) - if(gid.ne.-1) then + if (gid.ne.-1) then mean(gid,l) = param_hires_2d(c,r) endif enddo @@ -169,28 +172,28 @@ subroutine get_cdf_params (n, fname, month, param_a, param_b, mean, std) enddo ! read std - call LIS_verify(nf90_inq_varid(ftn,'std',std_id), & + call LIS_verify(nf90_inq_varid(ftn, 'std', std_id), & 'nf90_inq_varid failed for std in get_cdf_params') call LIS_verify(nf90_get_var(ftn, std_id, param_hires, & start=(/1,1,month/), count=(/ngrid,nlead,1/)),& 'nf90_get_var failed for std in get_cdf_params') - do l=1,nlead + do l = 1, nlead ! 1D -> 2D - do r=1,LIS_rc%gnr(n) - do c=1,LIS_rc%gnc(n) - param_hires_2d(c,r)=param_hires(c+(r-1)*LIS_rc%gnc(n),l) + do r = 1, LIS_rc%gnr(n) + do c = 1, LIS_rc%gnc(n) + param_hires_2d(c,r) =param_hires(c+(r-1)*LIS_rc%gnc(n),l) enddo enddo !subsets the data for each processor's domain nc = (LIS_ewe_halo_ind(n,LIS_localPet+1)-LIS_ews_halo_ind(n,LIS_localPet+1))+1 - do r=LIS_nss_halo_ind(n,LIS_localPet+1),LIS_nse_halo_ind(n,LIS_localPet+1) - do c=LIS_ews_halo_ind(n,LIS_localPet+1),LIS_ewe_halo_ind(n,LIS_localPet+1) - c1 = c-LIS_ews_halo_ind(n,LIS_localPet+1)+1 - r1 = r-LIS_nss_halo_ind(n,LIS_localPet+1)+1 + do r = LIS_nss_halo_ind(n,LIS_localPet+1), LIS_nse_halo_ind(n,LIS_localPet+1) + do c = LIS_ews_halo_ind(n,LIS_localPet+1), LIS_ewe_halo_ind(n,LIS_localPet+1) + c1 = c - LIS_ews_halo_ind(n,LIS_localPet+1) + 1 + r1 = r - LIS_nss_halo_ind(n,LIS_localPet+1) + 1 gid = LIS_domain(n)%gindex(c1,r1) - if(gid.ne.-1) then + if (gid.ne.-1) then std(gid,l) = param_hires_2d(c,r) endif enddo @@ -202,7 +205,8 @@ subroutine get_cdf_params (n, fname, month, param_a, param_b, mean, std) call LIS_verify(nf90_close(ftn),'failed to close in get_cdf_params') - write(LIS_logunit,*) 'Done reading MOGREPS-G bias correction parameters data ' + write(LIS_logunit,*) & + '[INFO] Done reading MOGREPS-G bias correction parameters data ' end subroutine get_cdf_params diff --git a/lis/metforcing/mogreps_g/get_mogrepsg.F90 b/lis/metforcing/mogreps_g/get_mogrepsg.F90 index 9877b0599..610802273 100644 --- a/lis/metforcing/mogreps_g/get_mogrepsg.F90 +++ b/lis/metforcing/mogreps_g/get_mogrepsg.F90 @@ -64,134 +64,156 @@ subroutine get_mogrepsg(n, findex) external :: read_mogrepsg ! MOGREPS-G cycles every 6 hours; ecch cycle provide up to 192 hours (8 days; 3-hour interval) forecast - if(LIS_rc%ts.gt.10800) then - write(LIS_logunit,*) '[WARN] The model timestep is > forcing data timestep ...' - write(LIS_logunit,*) '[WARN] LIS does not support this mode currently.' + if (LIS_rc%ts .gt. 10800) then + write(LIS_logunit,*) '[ERR] The model timestep is > forcing data timestep ...' + write(LIS_logunit,*) '[ERR] LIS does not support this mode currently.' call LIS_endrun() endif - openfile=0 + openfile = 0 - if(LIS_rc%tscount(n).eq.1 .or.LIS_rc%rstflag(n).eq.1) then !beginning of run + if (LIS_rc%tscount(n) .eq. 1 .or. LIS_rc%rstflag(n) .eq. 1) then !beginning of run LIS_rc%rstflag(n) = 0 endif ! First timestep of run - if(LIS_rc%tscount(n).eq.1 .or.LIS_rc%rstflag(n).eq.1) then - ! Bookend-time record 1 + if (LIS_rc%tscount(n) .eq. 1 .or. LIS_rc%rstflag(n) .eq. 1) then + ! Bookend-time record 1 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) + 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) ! Bookend-time record 2 - yr2=LIS_rc%yr !next hour - mo2=LIS_rc%mo - da2=LIS_rc%da - hr2=3 - mn2=0 - ss2=0 - ts2=0 - call LIS_tick(time2,doy2,gmt2,yr2,mo2,da2,hr2,mn2,ss2,ts2) + yr2 = LIS_rc%yr !next hour + mo2 = LIS_rc%mo + da2 = LIS_rc%da + hr2 = 3 + mn2 = 0 + ss2 = 0 + ts2 = 0 + call LIS_tick(time2, doy2, gmt2, yr2, mo2, da2, hr2, mn2, ss2, ts2) openfile=1 endif - ! 3 hourly interval + ! 3 hourly interval fcsthr_intv = 3 valid_hour = fcsthr_intv * (LIS_rc%hr/fcsthr_intv) - if((valid_hour==LIS_rc%hr .and. LIS_rc%mn==0) .or. & + if ((valid_hour == LIS_rc%hr .and. LIS_rc%mn == 0) .or. & openfile == 1) then ! Forecast hour condition within each file: - mogrepsg_struc(n)%fcst_hour = mogrepsg_struc(n)%fcst_hour + fcsthr_intv - + mogrepsg_struc(n)%fcst_hour = & + mogrepsg_struc(n)%fcst_hour + fcsthr_intv + ! Check if local forecast hour exceeds max grib file forecast hour: - if(mogrepsg_struc(n)%fcst_hour > 195 ) then + if (mogrepsg_struc(n)%fcst_hour > 195 ) then write(LIS_logunit,*) & - "[INFO] MOGREPS-G Forecast hour has exceeded the grib file's final" + "[ERR] MOGREPS-G Forecast hour has exceeded the grib file's final" write(LIS_logunit,*) & ' forecast hour (record). Run will end here for now ... ' call LIS_endrun endif - + ! Update bookend-time record 2: - if(LIS_rc%tscount(n).ne.1) then - mogrepsg_struc(n)%fcsttime1=mogrepsg_struc(n)%fcsttime2 - mogrepsg_struc(n)%metdata1(:,:,:)=mogrepsg_struc(n)%metdata2(:,:,:) - - yr2=LIS_rc%yr - mo2=LIS_rc%mo - da2=LIS_rc%da - hr2=valid_hour - mn2=fcsthr_intv*60 ! Backward looking - ss2=0 - ts2=0 - call LIS_tick(time2,doy2,gmt2,yr2,mo2,da2,hr2,mn2,ss2,ts2) + if (LIS_rc%tscount(n) .ne. 1) then + mogrepsg_struc(n)%fcsttime1 = mogrepsg_struc(n)%fcsttime2 + mogrepsg_struc(n)%metdata1(:,:,:) = & + mogrepsg_struc(n)%metdata2(:,:,:) + + yr2 = LIS_rc%yr + mo2 = LIS_rc%mo + da2 = LIS_rc%da + hr2 = valid_hour + mn2 = fcsthr_intv * 60 ! Backward looking + ss2 = 0 + ts2 = 0 + call LIS_tick(time2, doy2, gmt2, yr2, mo2, da2, hr2, mn2, ss2, & + ts2) endif - do m=1,mogrepsg_struc(n)%max_ens_members + do m = 1, mogrepsg_struc(n)%max_ens_members ! Read in file contents: - if(LIS_rc%tscount(n) == 1) then ! Read in first two book-ends - ferror=0 - order=1 - call get_mogrepsg_filename(mogrepsg_struc(n)%odir,mogrepsg_struc(n)%init_yr,& - mogrepsg_struc(n)%init_mo,mogrepsg_struc(n)%init_da,mogrepsg_struc(n)%init_hr,& - 0,m,fname) - - write(LIS_logunit,*)'[INFO] Getting MOGREPS-G forecast file1 ... ',trim(fname) + if (LIS_rc%tscount(n) == 1) then ! Read in first two book-ends + ferror = 0 + order = 1 + call get_mogrepsg_filename(mogrepsg_struc(n)%odir, & + mogrepsg_struc(n)%init_yr, & + mogrepsg_struc(n)%init_mo, & + mogrepsg_struc(n)%init_da, & + mogrepsg_struc(n)%init_hr, & + 0, m, fname) + + write(LIS_logunit,*)& + '[INFO] Getting MOGREPS-G forecast file1 ... ', & + trim(fname) call read_mogrepsg(n, m, findex, order, fname, ferror) - if(ferror.ge.1) mogrepsg_struc(n)%fcsttime1=time1 - - ferror=0 - order=2 - call get_mogrepsg_filename(mogrepsg_struc(n)%odir,mogrepsg_struc(n)%init_yr,& - mogrepsg_struc(n)%init_mo,mogrepsg_struc(n)%init_da,mogrepsg_struc(n)%init_hr,& - mogrepsg_struc(n)%fcst_hour,m,fname) - - write(LIS_logunit,*)'[INFO] Getting MOGREPS-G forecast file2 ... ',trim(fname) + if (ferror .ge. 1) mogrepsg_struc(n)%fcsttime1 = time1 + + ferror = 0 + order = 2 + call get_mogrepsg_filename(mogrepsg_struc(n)%odir, & + mogrepsg_struc(n)%init_yr, & + mogrepsg_struc(n)%init_mo, & + mogrepsg_struc(n)%init_da, & + mogrepsg_struc(n)%init_hr, & + mogrepsg_struc(n)%fcst_hour, m, fname) + + write(LIS_logunit,*) & + '[INFO] Getting MOGREPS-G forecast file2 ... ', & + trim(fname) call read_mogrepsg(n, m, findex, order, fname, ferror) - if(ferror.ge.1) mogrepsg_struc(n)%fcsttime2=time2 + if (ferror .ge. 1) mogrepsg_struc(n)%fcsttime2 = time2 - !only for T+0 due to mssing LW varaible - mogrepsg_struc(n)%metdata1(4,m,:) = mogrepsg_struc(n)%metdata2(4,m,:) + !only for T+0 due to mssing LW variable + mogrepsg_struc(n)%metdata1(4,m,:) = & + mogrepsg_struc(n)%metdata2(4,m,:) else - ferror=0 - order=2 + ferror = 0 + order = 2 ! met forcings except for pcp - call get_mogrepsg_filename(mogrepsg_struc(n)%odir,mogrepsg_struc(n)%init_yr,& - mogrepsg_struc(n)%init_mo,mogrepsg_struc(n)%init_da,mogrepsg_struc(n)%init_hr,& - mogrepsg_struc(n)%fcst_hour,m,fname) - - write(LIS_logunit,*)'[INFO] Getting MOGREPS-G forecast file2 ... ',trim(fname) + call get_mogrepsg_filename(mogrepsg_struc(n)%odir, & + mogrepsg_struc(n)%init_yr, & + mogrepsg_struc(n)%init_mo, & + mogrepsg_struc(n)%init_da, & + mogrepsg_struc(n)%init_hr, & + mogrepsg_struc(n)%fcst_hour, m, fname) + + write(LIS_logunit,*) & + '[INFO] Getting MOGREPS-G forecast file2 ... ', & + trim(fname) call read_mogrepsg(n, m, findex, order, fname, ferror) - if(ferror.ge.1) mogrepsg_struc(n)%fcsttime2=time2 + if (ferror .ge. 1) mogrepsg_struc(n)%fcsttime2 = time2 !only for T+141 due to mssing LW varaible - if(mogrepsg_struc(n)%fcst_hour == 141) then - mogrepsg_struc(n)%metdata2(4,m,:) = mogrepsg_struc(n)%metdata1(4,m,:) + if (mogrepsg_struc(n)%fcst_hour == 141) then + mogrepsg_struc(n)%metdata2(4,m,:) = & + mogrepsg_struc(n)%metdata1(4,m,:) endif endif ! apply precipitation bias correction (cdf from difference bewteen NAPFA and MOGREPS-G) if (mogrepsg_struc(n)%bc == 1) then - lead_time=floor((float(mogrepsg_struc(n)%fcst_hour))/24)+1 + lead_time=floor((float(mogrepsg_struc(n)%fcst_hour))/24) + 1 if (lead_time > 8) then lead_time = 8 endif - do t=1, LIS_rc%ngrid(n) - if(mogrepsg_struc(n)%metdata2(8,m,t).ne.LIS_rc%udef) then + do t = 1, LIS_rc%ngrid(n) + if (mogrepsg_struc(n)%metdata2(8,m,t) .ne. LIS_rc%udef) then ! only for land pixels if (mogrepsg_struc(n)%bc_param_a(t,lead_time) .ne. LIS_rc%udef) then ! perform centering and scaling - pcp1=mogrepsg_struc(n)%metdata2(8,m,t)-mogrepsg_struc(n)%metdata1(8,m,t) + pcp1= & + mogrepsg_struc(n)%metdata2(8,m,t) - & + mogrepsg_struc(n)%metdata1(8,m,t) if (mogrepsg_struc(n)%bc_std(t,lead_time) .ne. 0) then pcp2=(pcp1-mogrepsg_struc(n)%bc_mean(t,lead_time))/& mogrepsg_struc(n)%bc_std(t,lead_time) @@ -200,27 +222,29 @@ subroutine get_mogrepsg(n, findex) endif ! apply cdf params - pcp2=pcp2*mogrepsg_struc(n)%bc_param_a(t,lead_time)+mogrepsg_struc(n)%bc_param_b(t,lead_time) + pcp2 = pcp2 * & + mogrepsg_struc(n)%bc_param_a(t,lead_time)+mogrepsg_struc(n)%bc_param_b(t,lead_time) ! check for negative precipitation; if the corrected value has negative, keep the original value. - if(pcp2 >= 0) then - mogrepsg_struc(n)%pcp_bc(m,t)=pcp2 + if (pcp2 >= 0) then + mogrepsg_struc(n)%pcp_bc(m,t) = pcp2 else - mogrepsg_struc(n)%pcp_bc(m,t)=pcp1 + mogrepsg_struc(n)%pcp_bc(m,t) = pcp1 endif ! additionally, avoid bias correction for values that are too samll to reduce abnormal noise. if(pcp1 < 0.01) then - mogrepsg_struc(n)%pcp_bc(m,t)=pcp1 + mogrepsg_struc(n)%pcp_bc(m,t) = pcp1 endif else ! for water pixels - mogrepsg_struc(n)%pcp_bc(m,t)=mogrepsg_struc(n)%metdata2(8,m,t)-mogrepsg_struc(n)%metdata1(8,m,t) + mogrepsg_struc(n)%pcp_bc(m,t) = & + mogrepsg_struc(n)%metdata2(8,m,t) - & + mogrepsg_struc(n)%metdata1(8,m,t) endif endif enddo endif - enddo endif - openfile=0 + openfile = 0 end subroutine get_mogrepsg @@ -230,7 +254,8 @@ end subroutine get_mogrepsg ! \label{get_mogrepsg_filename} ! ! !INTERFACE: -subroutine get_mogrepsg_filename(rootdir,yr,mo,da,hr,fc_hr,ens_id,filename) +subroutine get_mogrepsg_filename(rootdir, yr, mo, da, hr, fc_hr, & + ens_id, filename) use LIS_logMod, only: LIS_endrun implicit none @@ -248,12 +273,11 @@ subroutine get_mogrepsg_filename(rootdir,yr,mo,da,hr,fc_hr,ens_id,filename) character(8) :: ftime character(2) :: chr character(3) :: fchr - character(2) :: ens - + character(2) :: ens character(len=36) :: fname write (UNIT=chr, FMT='(i2.2)') hr ! cycle 00/06/12/18 - write (UNIT=fchr, FMT='(i3.3)') fc_hr ! forecast time + write (UNIT=fchr, FMT='(i3.3)') fc_hr ! forecast time write (UNIT=ftime, FMT='(i4, i2.2, i2.2)') yr, mo, da fname = 'prods_op_mogreps-g_' @@ -265,7 +289,7 @@ subroutine get_mogrepsg_filename(rootdir,yr,mo,da,hr,fc_hr,ens_id,filename) else if (ens_id == 1) then write (UNIT=ens, FMT='(i2.2)') ens_id-1 ! start 00 - else + else write (UNIT=ens, FMT='(i2.2)') ens_id+16 ! start 18-34 endif endif diff --git a/lis/metforcing/mogreps_g/mogrepsg_forcingMod.F90 b/lis/metforcing/mogreps_g/mogrepsg_forcingMod.F90 index 98dc9d04a..dd27a7bab 100644 --- a/lis/metforcing/mogreps_g/mogrepsg_forcingMod.F90 +++ b/lis/metforcing/mogreps_g/mogrepsg_forcingMod.F90 @@ -47,14 +47,14 @@ module mogrepsg_forcingMod integer, allocatable :: gindex(:,:) integer :: mi - + integer, allocatable :: n111(:) integer, allocatable :: n121(:) integer, allocatable :: n211(:) integer, allocatable :: n221(:) real, allocatable :: w111(:),w121(:) real, allocatable :: w211(:),w221(:) - + integer, allocatable :: n112(:,:) integer, allocatable :: n122(:,:) integer, allocatable :: n212(:,:) @@ -63,14 +63,14 @@ module mogrepsg_forcingMod real, allocatable :: w212(:,:),w222(:,:) integer, allocatable :: n113(:) - + integer :: findtime1, findtime2 integer :: fcst_hour integer :: init_yr, init_mo, init_da, init_hr - real, allocatable :: metdata1(:,:,:) + real, allocatable :: metdata1(:,:,:) real, allocatable :: metdata2(:,:,:) - - integer :: nmodels + + integer :: nmodels ! only for v-wind due to difference resolution integer :: nrv @@ -78,15 +78,15 @@ module mogrepsg_forcingMod integer, allocatable :: nv121(:) integer, allocatable :: nv211(:) integer, allocatable :: nv221(:) - real, allocatable :: wv111(:),wv121(:) - real, allocatable :: wv211(:),wv221(:) + real, allocatable :: wv111(:), wv121(:) + real, allocatable :: wv211(:), wv221(:) integer, allocatable :: nv112(:,:) integer, allocatable :: nv122(:,:) integer, allocatable :: nv212(:,:) integer, allocatable :: nv222(:,:) - real, allocatable :: wv112(:,:),wv122(:,:) - real, allocatable :: wv212(:,:),wv222(:,:) + real, allocatable :: wv112(:,:), wv122(:,:) + real, allocatable :: wv212(:,:), wv222(:,:) integer, allocatable :: nv113(:) @@ -109,19 +109,19 @@ module mogrepsg_forcingMod ! ! !ROUTINE: init_mogrepsg ! \label{init_mogrepsg} -! +! ! !INTERFACE: subroutine init_mogrepsg(findex) -! !USES: +! !USES: use LIS_coreMod, only : LIS_rc use LIS_timeMgrMod, only : LIS_update_timestep use LIS_logMod, only : LIS_logunit, LIS_endrun implicit none -! !USES: +! !USES: integer, intent(in) :: findex -! -! !DESCRIPTION: +! +! !DESCRIPTION: ! Defines the native resolution of the input forcing for MOGREPS-G ! data. The grid description arrays are based on the decoding ! schemes used by NCEP and followed in the LIS interpolation @@ -139,12 +139,15 @@ subroutine init_mogrepsg(findex) external :: neighbor_interp_input external :: get_cdf_params - write(LIS_logunit,*) "[INFO] Initializing the MOGREPS-G forecast inputs " + write(LIS_logunit,*) & + "[INFO] Initializing the MOGREPS-G forecast inputs " ! Forecast mode -- NOT Available at this time for this forcing reader: if( LIS_rc%forecastMode.eq.1 ) then - write(LIS_logunit,*) '[ERR] Currently the MOGREPS-G forecast forcing reader' - write(LIS_logunit,*) '[ERR] is not set up to run in forecast mode.' + write(LIS_logunit,*) & + '[ERR] Currently the MOGREPS-G forecast forcing reader' + write(LIS_logunit,*) & + '[ERR] is not set up to run in forecast mode.' write(LIS_logunit,*) '[ERR] LIS forecast run-time ending.' call LIS_endrun() endif @@ -164,10 +167,10 @@ subroutine init_mogrepsg(findex) enddo ! 8 - key met field - LIS_rc%met_nf(findex) = 8 + LIS_rc%met_nf(findex) = 8 + + do n = 1, LIS_rc%nnest - do n=1,LIS_rc%nnest - ! Check if starting hour of LIS run matches 00/06/12/18 UTC: if((LIS_rc%shr .ne. 0) .and. (LIS_rc%shr .ne. 6) .and. & (LIS_rc%shr .ne. 12) .and. (LIS_rc%shr .ne. 18)) then @@ -177,10 +180,10 @@ subroutine init_mogrepsg(findex) write(LIS_logunit,*) "[ERR] your lis.config file.." call LIS_endrun() endif - + ! Allocate and initialize MOGREPS-G metforcing data structures: LIS_rc%met_nensem(findex) = mogrepsg_struc(n)%max_ens_members - + allocate(mogrepsg_struc(n)%metdata1(LIS_rc%met_nf(findex),& mogrepsg_struc(n)%max_ens_members,LIS_rc%ngrid(n))) allocate(mogrepsg_struc(n)%metdata2(LIS_rc%met_nf(findex),& @@ -196,7 +199,7 @@ subroutine init_mogrepsg(findex) mogrepsg_struc(n)%metdata1 = 0 mogrepsg_struc(n)%metdata2 = 0 gridDesci = 0 - + gridDesci(n,1) = 0 gridDesci(n,2) = real(mogrepsg_struc(n)%nc) !gnc gridDesci(n,3) = real(mogrepsg_struc(n)%nr) !gnr @@ -352,7 +355,7 @@ subroutine init_mogrepsg(findex) enddo ! precipitation bias correction - do n=1,LIS_rc%nnest + do n = 1, LIS_rc%nnest if (mogrepsg_struc(n)%bc == 1) then allocate(mogrepsg_struc(n)%pcp_bc(mogrepsg_struc(n)%max_ens_members,LIS_rc%ngrid(n))) allocate(mogrepsg_struc(n)%bc_param_a(LIS_rc%ngrid(n),8)) !8: lead time @@ -367,8 +370,9 @@ subroutine init_mogrepsg(findex) mogrepsg_struc(n)%bc_std = 0 ! read cdf parameters - call get_cdf_params(n,mogrepsg_struc(n)%cdf_fname,LIS_rc%mo,& - mogrepsg_struc(n)%bc_param_a, mogrepsg_struc(n)%bc_param_b,& + call get_cdf_params(n,mogrepsg_struc(n)%cdf_fname,LIS_rc%mo, & + mogrepsg_struc(n)%bc_param_a, & + mogrepsg_struc(n)%bc_param_b, & mogrepsg_struc(n)%bc_mean, mogrepsg_struc(n)%bc_std) endif enddo diff --git a/lis/metforcing/mogreps_g/readcrd_mogrepsg.F90 b/lis/metforcing/mogreps_g/readcrd_mogrepsg.F90 index 0ab191a4a..38cd984b4 100644 --- a/lis/metforcing/mogreps_g/readcrd_mogrepsg.F90 +++ b/lis/metforcing/mogreps_g/readcrd_mogrepsg.F90 @@ -16,7 +16,7 @@ ! 26 Jan 2023; Yeosang Yoon, Initial Code ! 01 Jan 2024; Yeosang Yoon; update codes for precpi. bias-correction ! -! !INTERFACE: +! !INTERFACE: subroutine readcrd_mogrepsg() ! !USES: use LIS_logMod @@ -26,55 +26,72 @@ subroutine readcrd_mogrepsg() ! ! !DESCRIPTION: ! -! This routine reads the options specific to MOGREPS-G forecast forcing from -! the LIS configuration file. -! +! This routine reads the options specific to MOGREPS-G forecast forcing from +! the LIS configuration file. +! !EOP implicit none - integer :: n,rc + integer :: n, rc - call ESMF_ConfigFindLabel(LIS_config,"MOGREPS-G forecast forcing directory:",rc=rc) - do n=1,LIS_rc%nnest - call ESMF_ConfigGetAttribute(LIS_config,mogrepsg_struc(n)%odir,rc=rc) - call LIS_verify(rc,'MOGREPS-G forecast forcing directory: not defined') + call ESMF_ConfigFindLabel(LIS_config, & + "MOGREPS-G forecast forcing directory:", rc=rc) + do n = 1, LIS_rc%nnest + call ESMF_ConfigGetAttribute(LIS_config, mogrepsg_struc(n)%odir, & + rc=rc) + call LIS_verify(rc, & + 'MOGREPS-G forecast forcing directory: not defined') enddo - call ESMF_ConfigFindLabel(LIS_config,"MOGREPS-G forecast run mode:",rc=rc) + call ESMF_ConfigFindLabel(LIS_config, "MOGREPS-G forecast run mode:", & + rc=rc) call LIS_verify(rc, 'MOGREPS-G forecast run mode: not defined ') - do n=1,LIS_rc%nnest - call ESMF_ConfigGetAttribute(LIS_config,mogrepsg_struc(n)%runmode,rc=rc) + do n = 1, LIS_rc%nnest + call ESMF_ConfigGetAttribute(LIS_config, & + mogrepsg_struc(n)%runmode, rc=rc) enddo - call ESMF_ConfigFindLabel(LIS_config,"MOGREPS-G forecast number of ensemble members:",rc=rc) - call LIS_verify(rc, 'MOGREPS-G forecast number of ensemble members: not defined') - do n=1,LIS_rc%nnest - call ESMF_ConfigGetAttribute(LIS_config,mogrepsg_struc(n)%max_ens_members,rc=rc) + call ESMF_ConfigFindLabel(LIS_config, & + "MOGREPS-G forecast number of ensemble members:", rc=rc) + call LIS_verify(rc, & + 'MOGREPS-G forecast number of ensemble members: not defined') + do n = 1, LIS_rc%nnest + call ESMF_ConfigGetAttribute(LIS_config, & + mogrepsg_struc(n)%max_ens_members, rc=rc) enddo - call ESMF_ConfigFindLabel(LIS_config,"Apply MOGREPS-G precipitation bias correction:",rc=rc) - call LIS_verify(rc, 'Apply MOGREPS-G precipitation bias correction: not defined') - do n=1,LIS_rc%nnest - call ESMF_ConfigGetAttribute(LIS_config,mogrepsg_struc(n)%bc,rc=rc) + call ESMF_ConfigFindLabel(LIS_config, & + "Apply MOGREPS-G precipitation bias correction:", rc=rc) + call LIS_verify(rc, & + 'Apply MOGREPS-G precipitation bias correction: not defined') + do n = 1, LIS_rc%nnest + call ESMF_ConfigGetAttribute(LIS_config, mogrepsg_struc(n)%bc, rc=rc) enddo - call ESMF_ConfigFindLabel(LIS_config,"MOGREPS-G model CDF file:",rc=rc) + call ESMF_ConfigFindLabel(LIS_config, "MOGREPS-G model CDF file:", rc=rc) call LIS_verify(rc, 'MOGREPS-G model CDF file: not defined') - do n=1,LIS_rc%nnest - call ESMF_ConfigGetAttribute(LIS_config,mogrepsg_struc(n)%cdf_fname,rc=rc) + do n = 1, LIS_rc%nnest + call ESMF_ConfigGetAttribute(LIS_config, & + mogrepsg_struc(n)%cdf_fname, rc=rc) enddo - do n=1,LIS_rc%nnest + do n = 1, LIS_rc%nnest write(LIS_logunit,*) '[INFO] Using MOGREPS-G forecast forcing' - write(LIS_logunit,*) '[INFO] MOGREPS-G forecast forcing directory: ', trim(mogrepsg_struc(n)%odir) - write(LIS_logunit,*) '[INFO] MOGREPS-G forecast run mode: ',mogrepsg_struc(n)%runmode - write(LIS_logunit,*) '[INFO] MOGREPS-G forecast number of ensemble members:',& + write(LIS_logunit,*) & + '[INFO] MOGREPS-G forecast forcing directory: ', & + trim(mogrepsg_struc(n)%odir) + write(LIS_logunit,*) '[INFO] MOGREPS-G forecast run mode: ', & + mogrepsg_struc(n)%runmode + write(LIS_logunit,*) & + '[INFO] MOGREPS-G forecast number of ensemble members:', & mogrepsg_struc(n)%max_ens_members - write(LIS_logunit,*) '[INFO] Using MOGREPS-G precipitation bias correction:',& + write(LIS_logunit,*) & + '[INFO] Using MOGREPS-G precipitation bias correction:',& mogrepsg_struc(n)%bc if (mogrepsg_struc(n)%bc == 1) then - write(LIS_logunit,*) '[INFO] MOGREPS-G model CDF file: ', trim(mogrepsg_struc(n)%cdf_fname) + write(LIS_logunit,*) '[INFO] MOGREPS-G model CDF file: ', & + trim(mogrepsg_struc(n)%cdf_fname) endif enddo end subroutine readcrd_mogrepsg diff --git a/lis/metforcing/mogreps_g/timeinterp_mogrepsg.F90 b/lis/metforcing/mogreps_g/timeinterp_mogrepsg.F90 index ea6537c25..c647d1c1a 100644 --- a/lis/metforcing/mogreps_g/timeinterp_mogrepsg.F90 +++ b/lis/metforcing/mogreps_g/timeinterp_mogrepsg.F90 @@ -59,150 +59,172 @@ subroutine timeinterp_mogrepsg(n,findex) real :: wt1, wt2, swt1, swt2 real :: gmt1, gmt2 integer :: t,index1 - integer :: bdoy,byr,bmo,bda,bhr,bmn - real*8 :: btime,newtime1,newtime2 - real :: tempgmt1,tempgmt2,tempbts - integer :: tempbdoy,tempbyr,tempbmo,tempbda,tempbhr,tempbmn,tempbss + integer :: bdoy, byr, bmo, bda, bhr, bmn + real*8 :: btime, newtime1, newtime2 + real :: tempgmt1, tempgmt2, tempbts + integer :: tempbdoy, tempbyr, tempbmo, tempbda, tempbhr, tempbmn, tempbss integer :: status - type(ESMF_Field) :: tmpField,q2Field,uField,vField,swdField,lwdField - type(ESMF_Field) :: psurfField,pcpField - real,pointer :: tmp(:),q2(:),uwind(:),vwind(:) - real,pointer :: swd(:),lwd(:),psurf(:),pcp(:) + type(ESMF_Field) :: tmpField, q2Field, uField, vField, swdField, & + lwdField + type(ESMF_Field) :: psurfField, pcpField + real, pointer :: tmp(:), q2(:), uwind(:), vwind(:) + real, pointer :: swd(:), lwd(:), psurf(:), pcp(:) integer :: mfactor, m, k, tid external :: zterp ! ________________________________________ btime=mogrepsg_struc(n)%fcsttime1 - 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=mogrepsg_struc(n)%fcsttime2 - 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) + 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 = mogrepsg_struc(n)%fcsttime2 + 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=real((mogrepsg_struc(n)%fcsttime2-LIS_rc%time)/ & + wt1 =real((mogrepsg_struc(n)%fcsttime2-LIS_rc%time)/ & (mogrepsg_struc(n)%fcsttime2-mogrepsg_struc(n)%fcsttime1)) - wt2=1.0-wt1 - swt1=real((newtime2-LIS_rc%time)/(newtime2-newtime1)) - swt2=1.0-swt1 + wt2 = 1.0 - wt1 + swt1 = real((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,& + 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 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,& + 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 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,& + 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 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,& + 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 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,& + 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 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,& + 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 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,& + 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 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,& + 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 LIS_verify(status, & + 'Error: Enable Rainf in the forcing variables list') - - call ESMF_FieldGet(tmpField,localDE=0,farrayPtr=tmp,rc=status) + call ESMF_FieldGet(tmpField, localDE=0, farrayPtr=tmp, rc=status) call LIS_verify(status) - call ESMF_FieldGet(q2Field,localDE=0,farrayPtr=q2,rc=status) + call ESMF_FieldGet(q2Field, localDE=0, farrayPtr=q2, rc=status) call LIS_verify(status) - call ESMF_FieldGet(swdField,localDE=0,farrayPtr=swd,rc=status) + call ESMF_FieldGet(swdField, localDE=0, farrayPtr=swd, rc=status) call LIS_verify(status) - call ESMF_FieldGet(lwdField,localDE=0,farrayPtr=lwd,rc=status) + call ESMF_FieldGet(lwdField, localDE=0, farrayPtr=lwd, rc=status) call LIS_verify(status) - call ESMF_FieldGet(uField,localDE=0,farrayPtr=uwind,rc=status) + call ESMF_FieldGet(uField, localDE=0, farrayPtr=uwind, rc=status) call LIS_verify(status) - call ESMF_FieldGet(vField,localDE=0,farrayPtr=vwind,rc=status) + call ESMF_FieldGet(vField, localDE=0, farrayPtr=vwind, rc=status) call LIS_verify(status) - call ESMF_FieldGet(psurfField,localDE=0,farrayPtr=psurf,rc=status) + call ESMF_FieldGet(psurfField, localDE=0, farrayPtr=psurf, rc=status) call LIS_verify(status) - call ESMF_FieldGet(pcpField,localDE=0,farrayPtr=pcp,rc=status) + call ESMF_FieldGet(pcpField, localDE=0, farrayPtr=pcp, rc=status) call LIS_verify(status) ! Metforcing ensemble member count factor - mfactor = LIS_rc%nensem(n)/mogrepsg_struc(n)%max_ens_members + mfactor = LIS_rc%nensem(n) / mogrepsg_struc(n)%max_ens_members ! Downward shortwave radiation (average): - do t=1,LIS_rc%ntiles(n)/LIS_rc%nensem(n) - do m=1,mogrepsg_struc(n)%max_ens_members - do k=1,mfactor - tid = (t-1)*LIS_rc%nensem(n)+(m-1)*mfactor+k + do t = 1, LIS_rc%ntiles(n) / LIS_rc%nensem(n) + do m = 1, mogrepsg_struc(n)%max_ens_members + do k = 1, mfactor + tid = (t - 1) * LIS_rc%nensem(n)+(m-1)*mfactor + k index1 = LIS_domain(n)%tile(tid)%index - zdoy=LIS_rc%doy + 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) - - if(mogrepsg_struc(n)%metdata1(3,m,index1).ne.LIS_rc%udef.and.& - mogrepsg_struc(n)%metdata2(3,m,index1).ne.LIS_rc%udef) then - swd(tid) = mogrepsg_struc(n)%metdata1(3,m,index1)*zw1+& - mogrepsg_struc(n)%metdata2(3,m,index1)*zw2 + 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) + + if (mogrepsg_struc(n)%metdata1(3,m,index1) .ne. LIS_rc%udef & + .and. & + mogrepsg_struc(n)%metdata2(3,m,index1).ne. LIS_rc%udef) & + then + swd(tid) = mogrepsg_struc(n)%metdata1(3,m,index1) * zw1 + & + mogrepsg_struc(n)%metdata2(3,m,index1) * zw2 !swd(tid) = zw1 * mogrepsg_struc(n)%metdata2(3,m,index1) ! In cases of small cos(zenith) angles, use linear weighting to avoid overly large weights - if((swd(tid).gt.mogrepsg_struc(n)%metdata1(3,m,index1).and. & - swd(tid).gt.mogrepsg_struc(n)%metdata2(3,m,index1)).and. & - (czb.lt.0.1.or.cze.lt.0.1))then + if ((swd(tid) .gt. mogrepsg_struc(n)%metdata1(3,m,index1) & + .and. & + swd(tid) .gt. mogrepsg_struc(n)%metdata2(3,m,index1)) & + .and. & + (czb .lt. 0.1 .or. cze .lt. 0.1)) then swd(tid) = mogrepsg_struc(n)%metdata1(3,m,index1)*swt1+ & mogrepsg_struc(n)%metdata2(3,m,index1)*swt2 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] sw radiation too high!!' write(unit=LIS_logunit,fmt=*)'[WARN] it is',swd(t) write(unit=LIS_logunit,fmt=*)'[WARN] mogrepsgdata1=',& mogrepsg_struc(n)%metdata1(3,m,index1) @@ -210,22 +232,24 @@ subroutine timeinterp_mogrepsg(n,findex) mogrepsg_struc(n)%metdata2(3,m,index1) write(unit=LIS_logunit,fmt=*)'[WARN] zw1=',zw1,'zw2=',zw2 swd(t) = LIS_CONST_SOLAR - write(unit=LIS_logunit,fmt=*)'[WARN] forcing set to ',swd(t) + write(unit=LIS_logunit,fmt=*) & + '[WARN] forcing set to ',swd(t) endif endif - if ((swd(t).ne.LIS_rc%udef).and.(swd(t).lt.0)) then - if (swd(t).gt.-0.00001) then + if ((swd(t) .ne. LIS_rc%udef) .and. (swd(t) .lt. 0)) then + if (swd(t) .gt. -0.00001) then swd(t) = 0.0 else - write(LIS_logunit,*) '[ERR] timeinterp_mogrepsg -- Stopping because ', & + write(LIS_logunit,*) & + '[ERR] timeinterp_mogrepsg -- Stopping because ', & 'forcing not udef but lt 0,' write(LIS_logunit,*)'[ERR] timeinterp_mogrepsg -- ', & - t,swd(t),mogrepsg_struc(n)%metdata2(3,m,index1), & - ' (',LIS_localPet,')' + t,swd(t),mogrepsg_struc(n)%metdata2(3,m,index1), & + ' (',LIS_localPet,')' call LIS_endrun endif - endif + endif enddo enddo enddo @@ -234,31 +258,35 @@ subroutine timeinterp_mogrepsg(n,findex) ! precip variable Block Interpolation !----------------------------------------------------------------------- ! Total precipitation field (accumulated): - do t=1,LIS_rc%ntiles(n)/LIS_rc%nensem(n) - do m=1,mogrepsg_struc(n)%max_ens_members - do k=1,mfactor - tid = (t-1)*LIS_rc%nensem(n)+(m-1)*mfactor+k + do t = 1, LIS_rc%ntiles(n) / LIS_rc%nensem(n) + do m = 1, mogrepsg_struc(n)%max_ens_members + do k = 1, mfactor + tid = (t - 1) * LIS_rc%nensem(n) + (m - 1) * mfactor + k index1 = LIS_domain(n)%tile(tid)%index ! apply precipitation bias correction if (mogrepsg_struc(n)%bc == 1) then !1 - use; or 0 - if(mogrepsg_struc(n)%pcp_bc(m,index1).ne.LIS_rc%udef) then + if (mogrepsg_struc(n)%pcp_bc(m,index1) .ne. & + LIS_rc%udef) then ! account for the accum fields - pcp(tid)=mogrepsg_struc(n)%pcp_bc(m,index1)/real(3600*3) - if(pcp(tid).lt.0) then + pcp(tid) = mogrepsg_struc(n)%pcp_bc(m,index1) / & + real(3600*3) + if (pcp(tid) .lt. 0) then pcp(tid) = 0.0 endif endif else !don't apply bias correction - if(mogrepsg_struc(n)%metdata2(8,m,index1).ne.LIS_rc%udef) then + if (mogrepsg_struc(n)%metdata2(8,m,index1) .ne. & + LIS_rc%udef) then ! account for the accum fields - pcp(tid)=(mogrepsg_struc(n)%metdata2(8,m,index1)-mogrepsg_struc(n)%metdata1(8,m,index1))/real(3600*3) - if(pcp(tid).lt.0) then + pcp(tid) = (mogrepsg_struc(n)%metdata2(8,m,index1) - & + mogrepsg_struc(n)%metdata1(8,m,index1)) / & + real(3600*3) + if (pcp(tid) .lt. 0) then pcp(tid) = 0.0 endif endif endif - enddo enddo enddo @@ -266,47 +294,59 @@ subroutine timeinterp_mogrepsg(n,findex) !----------------------------------------------------------------------- ! Linearly interpolate everything else !----------------------------------------------------------------------- - do t=1,LIS_rc%ntiles(n)/LIS_rc%nensem(n) - do m=1,mogrepsg_struc(n)%max_ens_members - do k=1,mfactor - tid = (t-1)*LIS_rc%nensem(n)+(m-1)*mfactor+k + do t = 1, LIS_rc%ntiles(n) / LIS_rc%nensem(n) + do m = 1, mogrepsg_struc(n)%max_ens_members + do k = 1, mfactor + tid = (t - 1)*LIS_rc%nensem(n) + (m - 1) * mfactor + k index1 = LIS_domain(n)%tile(tid)%index - + ! 2-meter air temp - if((mogrepsg_struc(n)%metdata1(1,m,index1).ne.LIS_rc%udef).and.& - (mogrepsg_struc(n)%metdata2(1,m,index1).ne.LIS_rc%udef)) then - tmp(tid) = mogrepsg_struc(n)%metdata1(1,m,index1)*wt1 + & - mogrepsg_struc(n)%metdata2(1,m,index1)*wt2 + if ((mogrepsg_struc(n)%metdata1(1,m,index1) .ne. LIS_rc%udef) & + .and. & + (mogrepsg_struc(n)%metdata2(1,m,index1) .ne. & + LIS_rc%udef)) then + tmp(tid) = mogrepsg_struc(n)%metdata1(1,m,index1) * wt1 + & + mogrepsg_struc(n)%metdata2(1,m,index1) * wt2 endif ! Specific humidity - if((mogrepsg_struc(n)%metdata1(2,m,index1).ne.LIS_rc%udef).and.& - (mogrepsg_struc(n)%metdata2(2,m,index1).ne.LIS_rc%udef)) then - q2(tid) = mogrepsg_struc(n)%metdata1(2,m,index1)*wt1 + & - mogrepsg_struc(n)%metdata2(2,m,index1)*wt2 + if ((mogrepsg_struc(n)%metdata1(2,m,index1) .ne. LIS_rc%udef) & + .and. & + (mogrepsg_struc(n)%metdata2(2,m,index1) .ne. & + LIS_rc%udef)) then + q2(tid) = mogrepsg_struc(n)%metdata1(2,m,index1) * wt1 + & + mogrepsg_struc(n)%metdata2(2,m,index1) * wt2 endif ! Downward longwave field - if((mogrepsg_struc(n)%metdata1(4,m,index1).ne.LIS_rc%udef).and.& - (mogrepsg_struc(n)%metdata2(4,m,index1).ne.LIS_rc%udef)) then - lwd(tid) = mogrepsg_struc(n)%metdata1(4,m,index1)*wt1 + & - mogrepsg_struc(n)%metdata2(4,m,index1)*wt2 + if ((mogrepsg_struc(n)%metdata1(4,m,index1) .ne. & + LIS_rc%udef) .and. & + (mogrepsg_struc(n)%metdata2(4,m,index1) .ne. & + LIS_rc%udef)) then + lwd(tid) = mogrepsg_struc(n)%metdata1(4,m,index1) * wt1 + & + mogrepsg_struc(n)%metdata2(4,m,index1) * wt2 endif ! U-wind component - if((mogrepsg_struc(n)%metdata1(5,m,index1).ne.LIS_rc%udef).and.& - (mogrepsg_struc(n)%metdata2(5,m,index1).ne.LIS_rc%udef)) then - uwind(tid) = mogrepsg_struc(n)%metdata1(5,m,index1)*wt1+& - mogrepsg_struc(n)%metdata2(5,m,index1)*wt2 + if ((mogrepsg_struc(n)%metdata1(5,m,index1) .ne. LIS_rc%udef) & + .and. & + (mogrepsg_struc(n)%metdata2(5,m,index1) .ne. & + LIS_rc%udef)) then + uwind(tid) = mogrepsg_struc(n)%metdata1(5,m,index1) * wt1 & + + mogrepsg_struc(n)%metdata2(5,m,index1) * wt2 endif ! V-wind component - if((mogrepsg_struc(n)%metdata1(6,m,index1).ne.LIS_rc%udef).and.& - (mogrepsg_struc(n)%metdata2(6,m,index1).ne.LIS_rc%udef)) then - vwind(tid) = mogrepsg_struc(n)%metdata1(6,m,index1)*wt1 + & - mogrepsg_struc(n)%metdata2(6,m,index1)*wt2 + if ((mogrepsg_struc(n)%metdata1(6,m,index1) .ne. LIS_rc%udef) & + .and. & + (mogrepsg_struc(n)%metdata2(6,m,index1) .ne. & + LIS_rc%udef)) then + vwind(tid) = mogrepsg_struc(n)%metdata1(6,m,index1) * wt1 & + + mogrepsg_struc(n)%metdata2(6,m,index1) * wt2 endif ! Surface pressure field - if((mogrepsg_struc(n)%metdata1(7,m,index1).ne.LIS_rc%udef).and.& - (mogrepsg_struc(n)%metdata2(7,m,index1).ne.LIS_rc%udef)) then - psurf(tid) = mogrepsg_struc(n)%metdata1(7,m,index1)*wt1 + & - mogrepsg_struc(n)%metdata2(7,m,index1)*wt2 + if ((mogrepsg_struc(n)%metdata1(7,m,index1) .ne. LIS_rc%udef) & + .and. & + (mogrepsg_struc(n)%metdata2(7,m,index1) .ne. & + LIS_rc%udef)) then + psurf(tid) = mogrepsg_struc(n)%metdata1(7,m,index1) * wt1 & + + mogrepsg_struc(n)%metdata2(7,m,index1) * wt2 endif enddo enddo