diff --git a/io/module_write_internal_state.F90 b/io/module_write_internal_state.F90 index b370793e1..07fe9da3a 100644 --- a/io/module_write_internal_state.F90 +++ b/io/module_write_internal_state.F90 @@ -49,6 +49,7 @@ module write_internal_state integer :: lat_start, lon_start integer :: lat_end, lon_end real :: latstart, latlast, lonstart, lonlast + real :: latse, latnw, lonse, lonnw integer,dimension(:),allocatable :: lat_start_wrtgrp, lon_start_wrtgrp integer,dimension(:),allocatable :: lat_end_wrtgrp, lon_end_wrtgrp real,dimension(:,:),allocatable :: lonPtr, latPtr diff --git a/io/module_wrt_grid_comp.F90 b/io/module_wrt_grid_comp.F90 index 5bfb1486d..15ed9152e 100644 --- a/io/module_wrt_grid_comp.F90 +++ b/io/module_wrt_grid_comp.F90 @@ -756,6 +756,30 @@ subroutine wrt_initialize_p1(wrt_comp, imp_state_write, exp_state_write, clock, enddo enddo wrt_int_state%post_maptype = 207 + rot_lon = lon1(n) + rot_lat = lat1(n) + call rtll(rot_lon, rot_lat, geo_lon, geo_lat, dble(cen_lon(n)), dble(cen_lat(n))) + if (geo_lon < 0.0) geo_lon = geo_lon + 360.0 + wrt_int_state%lonstart = geo_lon + wrt_int_state%latstart = geo_lat + rot_lon = lon2(n) + rot_lat = lat1(n) + call rtll(rot_lon, rot_lat, geo_lon, geo_lat, dble(cen_lon(n)), dble(cen_lat(n))) + if (geo_lon < 0.0) geo_lon = geo_lon + 360.0 + wrt_int_state%lonse = geo_lon + wrt_int_state%latse = geo_lat + rot_lon = lon1(n) + rot_lat = lat2(n) + call rtll(rot_lon, rot_lat, geo_lon, geo_lat, dble(cen_lon(n)), dble(cen_lat(n))) + if (geo_lon < 0.0) geo_lon = geo_lon + 360.0 + wrt_int_state%lonnw = geo_lon + wrt_int_state%latnw = geo_lat + rot_lon = lon2(n) + rot_lat = lat2(n) + call rtll(rot_lon, rot_lat, geo_lon, geo_lat, dble(cen_lon(n)), dble(cen_lat(n))) + if (geo_lon < 0.0) geo_lon = geo_lon + 360.0 + wrt_int_state%lonlast = geo_lon + wrt_int_state%latlast = geo_lat else if ( trim(output_grid(n)) == 'rotated_latlon_moving' ) then ! Do not compute lonPtr, latPtr here. Will be done in the run phase wrt_int_state%post_maptype = 207 diff --git a/io/post_fv3.F90 b/io/post_fv3.F90 index 634e55910..c959c30bf 100644 --- a/io/post_fv3.F90 +++ b/io/post_fv3.F90 @@ -37,6 +37,9 @@ subroutine post_run_fv3(wrt_int_state,mypei,mpicomp,lead_write, & ! 6)read 3D cloud fraction from cld_amt for GFDL MP, ! and from cldfra for other MPs. ! Jun 2022 J. Meng 2D decomposition +! Jul 2022 W. Meng 1)output lat/lon of four corner point for rotated +! lat-lon grid. +! 2)read instant model top logwave ! !----------------------------------------------------------------------- !*** run post on write grid comp @@ -249,7 +252,8 @@ subroutine post_getattr_fv3(wrt_int_state,grid_id) lonstartv, lonlastv, cenlonv, & latstartv, latlastv, cenlatv, & latstart_r,latlast_r,lonstart_r, & - lonlast_r, STANDLON, maptype, gridtype + lonlast_r, STANDLON, maptype, gridtype, & + latse,lonse,latnw,lonnw ! implicit none ! @@ -356,6 +360,14 @@ subroutine post_getattr_fv3(wrt_int_state,grid_id) lonstart_r = lonstart latlast_r = latlast lonlast_r = lonlast + lonstart = nint(wrt_int_state%lonstart*gdsdegr) + latstart = nint(wrt_int_state%latstart*gdsdegr) + lonse = nint(wrt_int_state%lonse*gdsdegr) + latse = nint(wrt_int_state%latse*gdsdegr) + lonnw = nint(wrt_int_state%lonnw*gdsdegr) + latnw = nint(wrt_int_state%latnw*gdsdegr) + lonlast = nint(wrt_int_state%lonlast*gdsdegr) + latlast = nint(wrt_int_state%latlast*gdsdegr) if(dlon(grid_id) con_g, fv => con_fvirt, rgas => con_rd, & @@ -571,7 +583,7 @@ subroutine set_postvars_fv3(wrt_int_state,mpicomp,setvar_atmfile, & real,external::FPVSNEW real,dimension(:,:),allocatable :: dummy, p2d, t2d, q2d, qs2d, & cw2d, cfr2d - character(len=80) :: fieldname, wrtFBName + character(len=80) :: fieldname, wrtFBName, flatlon type(ESMF_Grid) :: wrtGrid type(ESMF_Field) :: theField type(ESMF_Field), allocatable :: fcstField(:) @@ -714,14 +726,13 @@ subroutine set_postvars_fv3(wrt_int_state,mpicomp,setvar_atmfile, & ! time averaged cloud fraction, set acfrst to spval, ncfrst to 1 ! UNDERGROUND RUNOFF, bgroff ! inst incoming sfc longwave -! inst model top outgoing longwave,rlwtoa ! inst incoming sfc shortwave, rswin ! inst incoming clear sky sfc shortwave, rswinc ! inst outgoing sfc shortwave, rswout ! snow phase change heat flux, snopcx ! GFS does not use total momentum flux,sfcuvx !$omp parallel do default(none),private(i,j),shared(jsta,jend,im,spval,ista,iend), & -!$omp& shared(acfrcv,ncfrcv,acfrst,ncfrst,bgroff,rlwtoa,rswin,rswinc,rswout,snopcx,sfcuvx) +!$omp& shared(acfrcv,ncfrcv,acfrst,ncfrst,bgroff,rswin,rswinc,rswout,snopcx,sfcuvx) do j=jsta,jend do i=ista,iend acfrcv(i,j) = spval @@ -729,7 +740,6 @@ subroutine set_postvars_fv3(wrt_int_state,mpicomp,setvar_atmfile, & acfrst(i,j) = spval ncfrst(i,j) = 1.0 bgroff(i,j) = spval - rlwtoa(i,j) = spval rswinc(i,j) = spval enddo enddo @@ -1678,6 +1688,17 @@ subroutine set_postvars_fv3(wrt_int_state,mpicomp,setvar_atmfile, & enddo endif + ! outgoing model top logwave + if(trim(fieldname)=='ulwrf_toa') then + !$omp parallel do default(none) private(i,j) shared(jsta,jend,ista,iend,rlwtoa,arrayr42d,fillValue,spval) + do j=jsta,jend + do i=ista, iend + rlwtoa(i,j) = arrayr42d(i,j) + if( abs(arrayr42d(i,j)-fillValue)< small) rlwtoa(i,j) = spval + enddo + enddo + endif + ! time averaged incoming sfc shortwave if(trim(fieldname)=='dswrf_ave') then !$omp parallel do default(none) private(i,j) shared(jsta,jend,ista,iend,aswin,arrayr42d,fillValue,spval) @@ -3265,6 +3286,18 @@ subroutine set_postvars_fv3(wrt_int_state,mpicomp,setvar_atmfile, & ! !more fields need to be computed ! + +! write lat/lon of the four corner point for rotated lat-lon grid + if(mype == 0 .and. maptype == 207)then + write(flatlon,1001)ifhr + open(112,file=trim(flatlon),form='formatted',status='unknown') + write(112,1002)latstart/1000,lonstart/1000,& + latse/1000,lonse/1000,latnw/1000,lonnw/1000, & + latlast/1000,lonlast/1000 + 1001 format('latlons_corners.txt.f',I3.3) + 1002 format(4(I6,I7,X)) + close(112) + endif end subroutine set_postvars_fv3