SUBROUTINE INITPOST_NMM !$$$ SUBPROGRAM DOCUMENTATION BLOCK ! . . . ! SUBPROGRAM: INITPOST INITIALIZE POST FOR RUN ! PRGRMMR: RUSS TREADON ORG: W/NP2 DATE: 93-11-10 ! ! ABSTRACT: THIS ROUTINE INITIALIZES CONSTANTS AND ! VARIABLES AT THE START OF AN ETA MODEL OR POST ! PROCESSOR RUN. ! ! THIS ROUTINE ASSUMES THAT INTEGERS AND REALS ARE THE SAME SIZE ! . ! ! PROGRAM HISTORY LOG: ! 93-11-10 RUSS TREADON - ADDED DOCBLOC ! 98-05-29 BLACK - CONVERSION OF POST CODE FROM 1-D TO 2-D ! 99-01 20 TUCCILLO - MPI VERSION ! 01-10-25 H CHUANG - MODIFIED TO PROCESS HYBRID MODEL OUTPUT ! 02-06-19 MIKE BALDWIN - WRF VERSION ! 02-08-15 H CHUANG - UNIT CORRECTION AND GENERALIZE PROJECTION OPTIONS ! 03-07-25 H CHUANG - MODIFIED TO PROCESS NMM WRF ! 05-12-05 H CHUANG - ADD CAPABILITY TO OUTPUT OFF-HOUR FORECAST WHICH HAS ! NO INPACTS ON ON-HOUR FORECAST ! ! USAGE: CALL INIT ! INPUT ARGUMENT LIST: ! NONE ! ! OUTPUT ARGUMENT LIST: ! NONE ! ! OUTPUT FILES: ! NONE ! ! SUBPROGRAMS CALLED: ! UTILITIES: ! NONE ! LIBRARY: ! COMMON - CTLBLK ! LOOKUP ! SOILDEPTH ! ! ! ATTRIBUTES: ! LANGUAGE: FORTRAN ! MACHINE : CRAY C-90 !$$$ use vrbls3d, only: t, u, uh, v, vh, q, cwm, f_ice, f_rain, f_rimef, q,& qqw, qqr, qqs, qqi, qqg, qqw, cwm , q2, wh, pint, alpint, pmid,& omga, pmidv, zmid, rlwtt, rswtt, ttnd, tcucn, train, exch_h,& el_pbl, cfr, zint use vrbls2d, only: fis, cfrach, cfracl, cfracm, u10h, u10, v10h, v10,th10,& q10, tshltr, qshltr, pshltr, smstav, smstot, acfrcv, acfrst, ncfrcv,& ncfrst, ssroff, bgroff, sfcevp, sfcexc, vegfrc, acsnow, acsnom,& cmc, sst, mdltaux, mdltauy, thz0, qz0, uz0, vz0, qs, z0, pblh, mixht,& ustar, akhs, akms, ths, prec, cuprec, acprec, ancprc, cprate, cuppt,& lspa, cldefi, htop, hbot, htopd, czmean, rswout, rlwin, rlwtoa, sigt4,& radot, aswin, aswout, alwin, alwout, alwtoa, aswtoa, hbotd, htops,& hbots, sr, rswin, rswinc, czen, tg, soiltb, twbs, sfcshx, qwbs,& sfclhx, grnflx, subshx, potevp, sno, si, pctsno, ivgtyp, isltyp,& islope, albedo, albase, mxsnal, epsr, f use soil, only: smc, sh2o, stc, sldpth, sllevel use masks, only: lmv, lmh, htm, vtm, hbm2, sm, sice, gdlat, gdlon, dx, dy use params_mod, only: tfrz, g, rd, d608, rtd, dtr, erad use lookup_mod, only: thl, plq, ptbl, ttbl, rdq, rdth, rdp, rdthe, pl,& qs0, sqs, sthe, the0, ttblq, rdpq, rdtheq, stheq, the0q use ctlblk_mod, only: jsta, jend, nprec, jsta_2l, jend_2u, filename,& datahandle, datestr, ihrst, imin, sdat, spval, imp_physics, pt,& icu_physics, pdtop, nsoil, isf_surface_physics, jsta_m, jend_m,& avrain, avcnvc, ardsw, ardlw, asrfc, me, mpi_comm_comp, nphs, spl,& lsm, dt, dtq2,tsrfc, trdlw, trdsw, idat, ifhr, ifmin, restrt,& theat, tclod, tprec, alsl, lm, im, jm , submodelname use gridspec_mod, only: latstart, latlast, cenlat, lonstart, lonlast,& cenlon, dxval, dyval, maptype, gridtype, truelat1, truelat2,& psmapf ! use wrf_io_flags_mod ! implicit none ! ! INCLUDE/SET PARAMETERS. ! INCLUDE "mpif.h" ! ! This version of INITPOST shows how to initialize, open, read from, and ! close a NetCDF dataset. In order to change it to read an internal (binary) ! dataset, do a global replacement of _ncd_ with _int_. real :: dcenlat, dcenlon character(len=31) :: VarName integer :: Status, cen1, cen2 character startdate*19,SysDepInfo*80 ! ! NOTE: SOME INTEGER VARIABLES ARE READ INTO DUMMY ( A REAL ). THIS IS OK ! AS LONG AS REALS AND INTEGERS ARE THE SAME SIZE. ! ! ALSO, EXTRACT IS CALLED WITH DUMMY ( A REAL ) EVEN WHEN THE NUMBERS ARE ! INTEGERS - THIS IS OK AS LONG AS INTEGERS AND REALS ARE THE SAME SIZE. CHARACTER*4 RESTHR INTEGER IDATE(8),JDATE(8) INTEGER :: i_parent_start, j_parent_start ! ! DECLARE VARIABLES. ! REAL RINC(5) REAL ETA1(LM), ETA2(LM) REAL DUMMY ( IM, JM ) REAL DUMMY2 ( IM, JM ) REAL FI(IM,JM,2) REAL DUM3D ( IM+1, JM+1, LM+1 ) REAL DUM3D2 ( IM+1, JM+1, LM+1 ) !mp INTEGER IDUMMY ( IM, JM ) ! !jw integer ii,jj,js,je,jev,iyear,imn,iday,itmp,ioutcount,istatus, & nsrfc,nrdlw,nrdsw,nheat,nclod, & I,J,L,LL,N,LONEND,LATEND,IMM,INAV,IRTN, & IFDX,IFDY,IGDOUT,ICEN,JCEN ! integer iw, ie real TSPH,fact,dumcst,tstart,tmp real LAT ! ! Declarations for : ! putting 10 m wind on V points because copygb assume such INTEGER IE, IW !code from R.Rozumalski INTEGER latnm, latsm, lonem, lonwm, idxave, dlat, dlon, nlat, nlon !*********************************************************************** ! START INIT HERE. ! WRITE(6,*)'INITPOST: ENTER INITPOST' print*,'im,jm,lm= ',im,jm,lm ii=im/2 ! diagnostic print indices jj=(jsta+jend)/2 ll=lm/2 ! ! STEP 1. READ MODEL OUTPUT FILE ! ! !*** ! ! set default to not empty buket NSRFC=0 NRDLW=0 NRDSW=0 NHEAT=0 NCLOD=0 NPREC=0 ! LMH always = LM for sigma-type vert coord ! LMV always = LM for sigma-type vert coord do j = jsta_2l, jend_2u do i = 1, im LMV ( i, j ) = lm LMH ( i, j ) = lm end do end do ! HTM VTM all 1 for sigma-type vert coord do l = 1, lm do j = jsta_2l, jend_2u do i = 1, im HTM ( i, j, l ) = 1.0 VTM ( i, j, l ) = 1.0 end do end do end do ! ! how do I get the filename? ! fileName = '/ptmp/wx20mb/wrfout_01_030500' ! DateStr = '2002-03-05_18:00:00' ! how do I get the filename? call ext_ncd_ioinit(SysDepInfo,Status) print*,'called ioinit', Status call ext_ncd_open_for_read( trim(fileName), 0, 0, " ", & DataHandle, Status) print*,'called open for read', Status if ( Status /= 0 ) then print*,'error opening ',fileName, ' Status = ', Status ; stop endif ! get date/time info ! this routine will get the next time from the file, not using it print *,'DateStr before calling ext_ncd_get_next_time=',DateStr ! call ext_ncd_get_next_time(DataHandle, DateStr, Status) print *,'DateStri,Status,DataHandle = ',DateStr,Status,DataHandle ! The end j row is going to be jend_2u for all variables except for V. JS=JSTA_2L JE=JEND_2U IF (JEND_2U.EQ.JM) THEN JEV=JEND_2U+1 ELSE JEV=JEND_2U ENDIF ! ! Getting start time call ext_ncd_get_dom_ti_char(DataHandle,'START_DATE',startdate, & status ) ! patch for NMM WRF because it does not have start-date in output yet ! startdate='2003-04-17T00:00:00' print*,'startdate= ',startdate ! jdate=0 idate=0 read(startdate,15)iyear,imn,iday,ihrst,imin 15 format(i4,1x,i2,1x,i2,1x,i2,1x,i2) print*,'start yr mo day hr min =',iyear,imn,iday,ihrst,imin print*,'processing yr mo day hr min=',idat(3),idat(1),idat(2), & idat(4),idat(5) idate(1)=iyear idate(2)=imn idate(3)=iday idate(5)=ihrst idate(6)=imin SDAT(1)=imn SDAT(2)=iday SDAT(3)=iyear ! jdate(1)=idat(3) jdate(2)=idat(1) jdate(3)=idat(2) jdate(5)=idat(4) jdate(6)=idat(5) ! CALL W3DIFDAT(JDATE,IDATE,2,RINC) ! ifhr=nint(rinc(2)) CALL W3DIFDAT(JDATE,IDATE,0,RINC) ifhr=nint(rinc(2)+rinc(1)*24.) ifmin=nint(rinc(3)) print*,' in INITPOST ifhr ifmin fileName=',ifhr,ifmin,fileName ! Getting tstart call ext_ncd_get_dom_ti_real(DataHandle,'TSTART',tmp,1,ioutcount, & istatus) if(istatus==0)then tstart=tmp else tstart=0. end if print*,'status for getting TSTART= ',istatus print*,'TSTART= ',TSTART ! Getting restart RESTRT=.TRUE. ! set RESTRT default call ext_ncd_get_dom_ti_integer(DataHandle,'RESTARTBIN',itmp,1, & ioutcount,istatus) IF(itmp .LT. 1)THEN RESTRT=.FALSE. ELSE RESTRT=.TRUE. END IF print*,'status for getting RESTARTBIN= ',istatus print*,'Is this a restrt run? ',RESTRT ! IF(RESTRT)THEN ! ifhr=ifhr+NINT(tstart) ! print*,'new forecast hours for restrt run= ',ifhr ! END IF IF(tstart .GT. 1.0E-2)THEN ifhr=ifhr+NINT(tstart) rinc=0 idate=0 rinc(2)=-1.0*ifhr call w3movdat(rinc,jdate,idate) SDAT(1)=idate(2) SDAT(2)=idate(3) SDAT(3)=idate(1) IHRST=idate(5) print*,'new forecast hours for restrt run= ',ifhr print*,'new start yr mo day hr min =',sdat(3),sdat(1), & sdat(2),ihrst,imin END IF VarName='HBM2' HBM2=SPVAL call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im HBM2 ( i, j ) = dummy ( i, j ) end do end do ! OK, since all of the variables are dimensioned/allocated to be ! the same size, this means we have to be careful int getVariable ! to not try to get too much data. For example, ! DUM3D is dimensioned IM+1,JM+1,LM+1 but there might actually ! only be im,jm,lm points of data available for a particular variable. ! get 3-D variables VarName='T' call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D, & IM+1,1,JM+1,LM+1,IM,JS,JE,LM) do l = 1, lm do j = jsta_2l, jend_2u do i = 1, im t ( i, j, l ) = dum3d ( i, j, l ) ! if(l.eq.1)print*,'Debug: I,J,T= ',i,j,t ( i, j, l ) ! t ( i, j, l ) = dum3d ( i, j, l ) + 300. ! th ( i, j, l ) = dum3d ( i, j, l ) + 300. end do end do end do do l=1,lm if(jj.ge. jsta .and. jj.le.jend)print*,'sample L,T= ',L,T(ii,jj,l) end do ! VarName='T_ADJ' ! call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D, ! & IM+1,1,JM+1,LM+1,IM,JS,JE,LM) ! do l = 1, lm ! do j = jsta_2l, jend_2u ! do i = 1, im ! t_ADJ ( i, j, l ) = dum3d ( i, j, l ) ! end do ! end do ! end do ! do l=1,lm ! if(jj.ge. jsta .and. jj.le.jend)print*,'sample L,T_ADJ= ',L ! &,T_ADJ(ii,jj,l) ! end do VarName='U' call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D, & IM+1,1,JM+1,LM+1,IM,JS,JE,LM) do l = 1, lm do j = jsta_2l, jend_2u do i = 1, im u ( i, j, l ) = dum3d ( i, j, l ) UH( i, j, l ) = dum3d ( i, j, l ) ! if(l.eq.1)print*,'Debug: I,J,U= ',i,j,u( i, j, l ) end do end do ! fill up UH which is U at P-points including 2 row halo ! do j = jsta_2l, jend_2u ! do i = 1, im ! UH (I,J,L) = (dum3d(I,J,L)+dum3d(I+1,J,L))*0.5 ! end do ! end do end do if(jj.ge. jsta .and. jj.le.jend)print*,'sample U= ',U(ii,jj,ll) VarName='V' call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D, & IM+1,1,JM+1,LM+1,IM,JS,JE,LM) do l = 1, lm do j = jsta_2l, jend_2u do i = 1, im v ( i, j, l ) = dum3d ( i, j, l ) VH( i, j, l ) = dum3d ( i, j, l ) end do end do ! fill up VH which is V at P-points including 2 row halo ! do j = jsta_2l, jend_2u ! do i = 1, im ! VH(I,J,L) = (dum3d(I,J,L)+dum3d(I,J+1,L))*0.5 ! end do ! end do end do if(jj.ge. jsta .and. jj.le.jend)print*,'sample V= ',V(ii,jj,ll) call ext_ncd_get_dom_ti_integer(DataHandle,'MP_PHYSICS' & ,itmp,1,ioutcount,istatus) imp_physics=itmp ! Chuang: will initialize microphysics constants differently for 85 now ! if(imp_physics .eq. 85) imp_physics=5 !HWRF print*,'MP_PHYSICS= ',imp_physics ! Initializes constants for Ferrier microphysics if(imp_physics==5 .or. imp_physics==85 .or. imp_physics==95)then CALL MICROINIT(imp_physics) end if call ext_ncd_get_dom_ti_integer(DataHandle,'CU_PHYSICS' & ,itmp,1,ioutcount,istatus) icu_physics=itmp if (icu_physics .eq. 84 .or. icu_physics .eq. 85) icu_physics = 4 ! HWRF print*,'CU_PHYSICS= ',icu_physics if(imp_physics==5 .or. imp_physics==85 .or. imp_physics==95)then VarName='Q' call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D, & IM+1,1,JM+1,LM+1,IM,JS,JE,LM) do l = 1, lm do j = jsta_2l, jend_2u do i = 1, im if (dum3d(i,j,l) .lt. 10E-12) dum3d(i,j,l) = 10E-12 q ( i, j, l ) = dum3d ( i, j, l ) end do end do end do print*,'finish reading specific humidity' if(jj.ge. jsta .and. jj.le.jend)print*,'sample Q= ',Q(ii,jj,ll) VarName='CWM' !????? call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D, & IM+1,1,JM+1,LM+1,IM,JS,JE,LM) do l = 1, lm do j = jsta_2l, jend_2u do i = 1, im cwm ( i, j, l ) = dum3d ( i, j, l ) end do end do end do print*,'finish reading cloud mixing ratio' VarName='F_ICE' call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D, & IM+1,1,JM+1,LM+1,IM,JS,JE,LM) do l = 1, lm do j = jsta_2l, jend_2u do i = 1, im F_ICE ( i, j, l ) = dum3d ( i, j, l ) end do end do end do VarName='F_RAIN' call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D, & IM+1,1,JM+1,LM+1,IM,JS,JE,LM) do l = 1, lm do j = jsta_2l, jend_2u do i = 1, im F_RAIN ( i, j, l ) = dum3d ( i, j, l ) end do end do end do VarName='F_RIMEF' call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D, & IM+1,1,JM+1,LM+1,IM,JS,JE,LM) do l = 1, lm do j = jsta_2l, jend_2u do i = 1, im F_RIMEF ( i, j, l ) = dum3d ( i, j, l ) end do end do end do else ! retrieve hydrometeo fields directly for non-Ferrier VarName='QVAPOR' call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D, & IM+1,1,JM+1,LM+1,IM,JS,JE,LM) do l = 1, lm do j = jsta_2l, jend_2u do i = 1, im ! q ( i, j, l ) = dum3d ( i, j, l ) ! if(l.eq.1)print*,'Debug: I,J,Q= ',i,j,q( i, j, l ) !CHC CONVERT MIXING RATIO TO SPECIFIC HUMIDITY if (dum3d(i,j,l) .lt. 10E-12) dum3d(i,j,l) = 10E-12 q ( i, j, l ) = dum3d ( i, j, l )/(1.0+dum3d ( i, j, l )) end do end do end do print*,'finish reading specific humidity' if(jj.ge. jsta .and. jj.le.jend)print*,'sample Q= ',Q(ii,jj,ll) qqw=spval qqr=spval qqs=spval qqi=spval qqg=spval cwm=spval f_rimef=spval if(imp_physics/=0)then VarName='QCLOUD' call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D, & IM+1,1,JM+1,LM+1,IM,JS,JE,LM) do l = 1, lm do j = jsta_2l, jend_2u do i = 1, im ! partition cloud water and ice for WSM3 if(imp_physics.eq.3)then if(t(i,j,l) .ge. TFRZ)then qqw ( i, j, l ) = dum3d ( i, j, l ) else qqi ( i, j, l ) = dum3d ( i, j, l ) end if else ! bug fix provided by J CASE qqw ( i, j, l ) = dum3d ( i, j, l ) end if cwm(i,j,l)=dum3d(i,j,l) end do end do end do end if if(jj.ge. jsta .and. jj.le.jend)print*,'sample qqw= ' & ,Qqw(ii,jj,ll) if(imp_physics.ne.1 .and. imp_physics.ne.3 & .and. imp_physics.ne.0)then VarName='QICE' call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D, & IM+1,1,JM+1,LM+1,IM,JS,JE,LM) do l = 1, lm do j = jsta_2l, jend_2u do i = 1, im qqi ( i, j, l ) = dum3d ( i, j, l ) cwm(i,j,l)=cwm(i,j,l)+dum3d(i,j,l) end do end do end do end if if(jj.ge. jsta .and. jj.le.jend)print*,'sample qqi= ' & ,Qqi(ii,jj,ll) if(imp_physics.ne.0)then VarName='QRAIN' call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D, & IM+1,1,JM+1,LM+1,IM,JS,JE,LM) do l = 1, lm do j = jsta_2l, jend_2u do i = 1, im ! partition rain and snow for WSM3 if(imp_physics .eq. 3)then if(t(i,j,l) .ge. TFRZ)then qqr ( i, j, l ) = dum3d ( i, j, l ) else qqs ( i, j, l ) = dum3d ( i, j, l ) end if else ! bug fix provided by J CASE qqr ( i, j, l ) = dum3d ( i, j, l ) end if cwm(i,j,l)=cwm(i,j,l)+dum3d(i,j,l) end do end do end do end if if(jj.ge. jsta .and. jj.le.jend)print*,'sample qqr= ' & ,Qqr(ii,jj,ll) if(imp_physics.ne.1 .and. imp_physics.ne.3 & .and. imp_physics.ne.0)then VarName='QSNOW' call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D, & IM+1,1,JM+1,LM+1,IM,JS,JE,LM) do l = 1, lm do j = jsta_2l, jend_2u do i = 1, im qqs ( i, j, l ) = dum3d ( i, j, l ) cwm(i,j,l)=cwm(i,j,l)+dum3d(i,j,l) end do end do end do end if if(jj.ge. jsta .and. jj.le.jend)print*,'sample qqs= ' & ,Qqs(ii,jj,ll) if(imp_physics.eq.2 .or. imp_physics.eq.6 & .or. imp_physics.eq.8 .or. imp_physics.eq.28)then VarName='QGRAUP' call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D, & IM+1,1,JM+1,LM+1,IM,JS,JE,LM) do l = 1, lm do j = jsta_2l, jend_2u do i = 1, im qqg ( i, j, l ) = dum3d ( i, j, l ) cwm(i,j,l)=cwm(i,j,l)+dum3d(i,j,l) end do end do end do end if if(jj.ge. jsta .and. jj.le.jend)print*,'sample qqg= ' & ,Qqg(ii,jj,ll) end if ! end of retrieving hydrometeo for different MP options ! call getVariable(fileName,DateStr,DataHandle,'TKE_PBL',DUM3D, call getVariable(fileName,DateStr,DataHandle,'Q2',DUM3D, & IM+1,1,JM+1,LM+1,IM,JS,JE,LM) do l = 1, lm do j = jsta_2l, jend_2u do i = 1, im q2 ( i, j, l ) = dum3d ( i, j, l ) end do end do end do VarName='W' call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D, & IM+1,1,JM+1,LM+1,IM,JS,JE,LM+1) ! & IM+1,1,JM+1,LM+1,IM,JS,JE,LM) ! do l = 1, lm+1 ! do j = jsta_2l, jend_2u ! do i = 1, im ! w ( i, j, l ) = dum3d ( i, j, l ) ! end do ! end do ! end do ! fill up WH which is W at P-points including 2 row halo DO L=1,LM DO I=1,IM DO J=JSTA_2L,JEND_2U ! WH(I,J,L) = (W(I,J,L)+W(I,J,L+1))*0.5 wh ( i, j, l ) = dum3d ( i, j, l+1 ) ENDDO ENDDO ENDDO print*,'finish reading W' !MEB call getVariable(fileName,DateStr,DataHandle,'QRAIN',new) VarName='PINT' call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D, & IM+1,1,JM+1,LM+1,IM,JS,JE,LM+1) ! VarName='P' ! call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D2, ! & IM+1,1,JM+1,LM+1,IM,JS,JE,LM) do l = 1, lm+1 do j = jsta_2l, jend_2u do i = 1, im ! PMID(I,J,L)=DUM3D(I,J,L)+DUM3D2(I,J,L) PINT(I,J,L)=DUM3D(I,J,L) ALPINT(I,J,L)=ALOG(PINT(I,J,L)) end do end do end do ! do l = 1, lm+1 ! if(jj.ge. jsta .and. jj.le.jend)print*,'sample PINT= ' ! & ,PINT(ii,jj,l) ! end do ! DO L=1,LM DO I=1,IM DO J=JSTA_2L,JEND_2U PMID(I,J,L)=(PINT(I,J,L)+PINT(I,J,L+1))*0.5 ! TH(I,J,L)=T(I,J,L)*(1.E5/PMID(I,J,L))**CAPA IF(ABS(T(I,J,L)).GT.1.0E-3) & OMGA(I,J,L) = -WH(I,J,L)*PMID(I,J,L)*G/ & (RD*T(I,J,L)*(1.+D608*Q(I,J,L))) ! ! PINT(I,J,L)=EXP((ALOG(PMID(I,J,L-1))+ ! & ALOG(PMID(I,J,L)))*0.5) ! ave of ln p ! ALPINT(I,J,L)=ALOG(PINT(I,J,L)) ENDDO ENDDO ENDDO ! do l = 1, lm do j = jsta, jend do i = 1, im-MOD(J,2) IF(J .EQ. 1 .AND. I .LT. IM)THEN !SOUTHERN BC PMIDV(I,J,L)=0.5*(PMID(I,J,L)+PMID(I+1,J,L)) ELSE IF(J.EQ.JM .AND. I.LT.IM)THEN !NORTHERN BC PMIDV(I,J,L)=0.5*(PMID(I,J,L)+PMID(I+1,J,L)) ELSE IF(I .EQ. 1 .AND. MOD(J,2) .EQ. 0) THEN !WESTERN EVEN BC PMIDV(I,J,L)=0.5*(PMID(I,J-1,L)+PMID(I,J+1,L)) ELSE IF(I .EQ. IM .AND. MOD(J,2) .EQ. 0 & .AND. J .LT. JM) THEN !EASTERN EVEN BC PMIDV(I,J,L)=0.5*(PMID(I,J-1,L)+PMID(I,J+1,L)) ELSE IF (MOD(J,2) .LT. 1) THEN PMIDV(I,J,L)=0.25*(PMID(I,J,L)+PMID(I-1,J,L) & +PMID(I,J+1,L)+PMID(I,J-1,L)) ELSE PMIDV(I,J,L)=0.25*(PMID(I,J,L)+PMID(I+1,J,L) & +PMID(I,J+1,L)+PMID(I,J-1,L)) END IF end do end do end do ! get sfc pressure ! VarName='MU' ! call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, ! & IM,1,JM,1,IM,JS,JE,1) ! VarName='MUB' ! call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY2, ! & IM,1,JM,1,IM,JS,JE,1) VarName='PT' call getVariable(fileName,DateStr,DataHandle,VarName,PT, & 1,1,1,1,1,1,1,1) VarName='PDTOP' call getVariable(fileName,DateStr,DataHandle,VarName,PDTOP, & 1,1,1,1,1,1,1,1) ! patch for no pt in wrf output ! pt=2500. ! print*,'PT= ',pt ! do j = jsta_2l, jend_2u ! do i = 1, im ! PINT(I,J,1)=pt ! ALPINT(I,J,1)=ALOG(PINT(I,J,1)) ! end do ! end do ! ! if(jj.ge. jsta .and. jj.le.jend)then do l = 1, lm+1 print*,'sample PINT= ',ii,jj,l,PINT(ii,jj,l) if(l.le.lm)print*,'sample PMID=',l,PMID(II,JJ,L) end do end if ! DO I=1,IM ! DO J=JS,JE ! PINT (I,J,LM+1) = DUMMY(I,J)+DUMMY2(I,J)+PT ! PINT (I,J,1) = PT ! ALPINT(I,J,LM+1)=ALOG(PINT(I,J,LM+1)) ! ALPINT(I,J,1)=ALOG(PINT(I,J,1)) ! ENDDO ! ENDDO ! NO HEIGHT OUTPUT IN NMM -> DERIVE IT FROM HYDROSTATIC RELATION ! VarName='PHB' ! call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D, ! & IM+1,1,JM+1,LM+1,IM,JS,JE,LM+1) ! VarName='PH' ! call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D2, ! & IM+1,1,JM+1,LM+1,IM,JS,JE,LM+1) ! FIRST, OBTAIN TERRAIN HEIGHT VarName='FIS' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im FIS ( i, j ) = dummy ( i, j ) ZINT(I,J,LM+1)=FIS(I,J)/G FI(I,J,1)=FIS(I,J) end do end do print*,'FIS at ',ii,jj,' = ',FIS(ii,jj) ! SECOND, INTEGRATE HEIGHT HYDROSTATICLY DO L=LM,1,-1 do j = jsta_2l, jend_2u do i = 1, im FI(I,J,2)=HTM(I,J,L)*T(I,J,L)*(Q(I,J,L)*D608+1.0)*RD* & (ALPINT(I,J,L+1)-ALPINT(I,J,L))+FI(I,J,1) ZINT(I,J,L)=FI(I,J,2)/G if(i.eq.ii.and.j.eq.jj) & print*,'L,sample HTM,T,Q,ALPINT(L+1),ALPINT(l),ZINT= ', & l,HTM(I,J,L),T(I,J,L),Q(I,J,L),ALPINT(I,J,L+1), & ALPINT(I,J,L),ZINT(I,J,L) FI(I,J,1)=FI(I,J,2) ENDDO ENDDO END DO print*,'finish deriving geopotential in nmm' ! DO L=1,LM DO I=1,IM DO J=JS,JE ! ZMID(I,J,L)=(ZINT(I,J,L+1)+ZINT(I,J,L))*0.5 ! ave of z FACT=(ALOG(PMID(I,J,L))-ALOG(PINT(I,J,L)))/ & (ALOG(PINT(I,J,L+1))-ALOG(PINT(I,J,L))) ZMID(I,J,L)=ZINT(I,J,L)+(ZINT(I,J,L+1)-ZINT(I,J,L))*FACT ENDDO ENDDO ENDDO ! get 3-d soil variables VarName='SMC' call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D, & IM+1,1,JM+1,LM+1,IM,JS,JE,NSOIL) do l = 1, nsoil do j = jsta_2l, jend_2u do i = 1, im ! smc ( i, j, l ) = dum3d ( i, j, l ) ! flip soil layer again because wrf soil variable vertical indexing ! is the same with eta and vertical indexing was flipped for both ! atmospheric and soil layers within getVariable smc ( i, j, l ) = dum3d ( i, j, nsoil-l+1) end do end do end do VarName='SH2O' call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D, & IM+1,1,JM+1,LM+1,IM,JS,JE,NSOIL) do l = 1, nsoil do j = jsta_2l, jend_2u do i = 1, im sh2o( i, j, l ) = dum3d ( i, j, nsoil-l+1) end do end do end do VarName='STC' call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D, & IM+1,1,JM+1,LM+1,IM,JS,JE,NSOIL) do l = 1, nsoil do j = jsta_2l, jend_2u do i = 1, im ! stc ( i, j, l ) = dum3d ( i, j, l ) stc ( i, j, l ) = dum3d ( i, j, nsoil-l+1) end do end do end do VarName='CFRACH' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im CFRACH ( i, j ) = dummy ( i, j ) ! print*,'Debug: I,J,TSHLTR= ',i,j,TSHLTR(i,j) end do end do print*,'CFRACH at ',ii,jj,' = ',CFRACH(ii,jj) VarName='CFRACL' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im CFRACL ( i, j ) = dummy ( i, j ) ! print*,'Debug: I,J,TSHLTR= ',i,j,TSHLTR(i,j) end do end do print*,'CFRACL at ',ii,jj,' = ',CFRACL(ii,jj) VarName='CFRACM' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im CFRACM ( i, j ) = dummy ( i, j ) ! print*,'Debug: I,J,TSHLTR= ',i,j,TSHLTR(i,j) end do end do print*,'CFRACM at ',ii,jj,' = ',CFRACM(ii,jj) ! Soil Layer/depth SLDPTH = 0.0 ! RUC LSM - use depths of center of soil layer if (iSF_SURFACE_PHYSICS==3)then !RUC LSM VarName='SLDPTH' call getVariable(fileName,DateStr,DataHandle,VarName,SLLEVEL, & NSOIL,1,1,1,NSOIL,1,1,1) print*,'SLLEVEL= ', (SLLEVEL(N),N=1,NSOIL) ! other LSM - use thickness of soil layer else VarName='DZSOIL' call getVariable(fileName,DateStr,DataHandle,VarName,SLDPTH, & NSOIL,1,1,1,NSOIL,1,1,1) print*,'SLDPTH= ',(SLDPTH(N),N=1,NSOIL) end if ! get 10m variables ! Chuang Aug 2012: 10 m winds are computed on mass points in the model ! post interpolates them onto V points because copygb interpolates ! wind points differently and 10 m winds are identified as 33/34 VarName='U10' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) U10 = SPVAL ! Wind on V point for output to copygb DO J=JSTA_M,JEND_M DO I=2,IM-1 u10h(i,j)=dummy(i,j) IE=I+MOD(J,2) IW=IE-1 u10(i,j)=(dummy(IW,J)+dummy(IE,J) & ! assuming e grid +dummy(I,J+1)+dummy(I,J-1))/4.0 END DO u10(1,j)=0.5*(dummy(1,j)+dummy(1,j+1)) u10h(1,j)=dummy(1,j) u10(im,j)=0.5*(dummy(im,j)+dummy(im,j+1)) u10h(im,j)=dummy(im,j) END DO ! Complete first row IF (JSTA_M.EQ.2) THEN DO I=1, IM-1 u10(I,1)=0.5*(dummy(I,1)+dummy(I+1,1)) u10h(I,1)=dummy(I,1) END DO u10(im,1) = dummy(im,1) u10h(im,1) = dummy(im,1) END IF ! Complete last row IF (JEND_M.EQ.(JM-1)) THEN DO I=1, IM-1 u10(I,jm)=0.5*(dummy(I,jm)+dummy(I+1,jm)) u10h(I,jm)=dummy(I,jm) END DO u10(im,jm) = dummy(im,jm) u10h(im,jm) = dummy(im,jm) END IF VarName='V10' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) V10 = SPVAL ! Wind on V point for output to copygb DO J=JSTA_M,JEND_M DO I=2,IM-1 v10h(i,j)=dummy(i,j) IE=I+MOD(J,2) IW=IE-1 v10(i,j)=(dummy(IW,J)+dummy(IE,J) & ! assuming e grid +dummy(I,J+1)+dummy(I,J-1))/4.0 END DO v10(1,j)=0.5*(dummy(1,j-1)+dummy(1,j+1)) v10h(1,j)=dummy(1,j) v10(im,j)=0.5*(dummy(im,j-1)+dummy(im,j+1)) v10h(im,j)=dummy(im,j) END DO ! Complete first row IF (JSTA_M.EQ.2) THEN DO I=1, IM-1 v10(I,1)=0.5*(dummy(I,1)+dummy(I+1,1)) v10h(I,1)=dummy(I,1) END DO v10(im,1) = dummy(im,1) v10h(im,1) = dummy(im,1) END IF ! Complete last row IF (JEND_M.EQ.(JM-1)) THEN DO I=1, IM-1 v10(I,jm)=0.5*(dummy(I,jm)+dummy(I+1,jm)) v10h(I,jm)=dummy(I,jm) END DO v10(im,jm) = dummy(im,jm) v10h(im,jm) = dummy(im,jm) END IF VarName='TH10' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im TH10 ( i, j ) = dummy ( i, j ) end do end do VarName='Q10' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im Q10 ( i, j ) = dummy ( i, j ) end do end do print*,'Q10 at ',ii,jj,' = ',Q10(ii,jj) ! get 2-m theta ! VarName='TH2' VarName='TSHLTR' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im TSHLTR ( i, j ) = dummy ( i, j ) ! print*,'Debug: I,J,TSHLTR= ',i,j,TSHLTR(i,j) end do end do print*,'TSHLTR at ',ii,jj,' = ',TSHLTR(ii,jj) ! get 2-m specific humidity ! VarName='Q2' VarName='QSHLTR' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im QSHLTR ( i, j ) = dummy ( i, j ) end do end do print*,'QSHLTR at ',ii,jj,' = ',QSHLTR(ii,jj) VarName='PSHLTR' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im PSHLTR ( i, j ) = dummy ( i, j ) end do end do print*,'PSHLTR at ',ii,jj,' = ',QSHLTR(ii,jj) VarName='SMSTAV' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im SMSTAV ( i, j ) = dummy ( i, j ) end do end do VarName='SMSTOT' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im SMSTOT ( i, j ) = dummy ( i, j ) end do end do print*,'SMSTOT at ',ii,jj,' = ',SMSTOT(ii,jj) VarName='ACFRCV' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im ACFRCV ( i, j ) = dummy ( i, j ) end do end do write(6,*) 'MAX ACFRCV: ', maxval(DUMMY) VarName='ACFRST' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im ACFRST ( i, j ) = dummy ( i, j ) end do end do write(6,*) 'max ACFRST ', maxval(DUMMY) varname='RLWTT' write(6,*) 'call getVariable for : ', VarName call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D, & IM+1,1,JM+1,LM+1,IM,JS,JE,LM) do l = 1, lm do j = jsta_2l, jend_2u do i = 1, im rlwtt( i, j, l ) = dum3d ( i, j, l ) end do end do end do ! write(6,*) 'RLWTT(II,jj,ll): ', DUM3D(ii,jj,ll) varname='RSWTT' write(6,*) 'call getVariable for : ', VarName call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D, & IM+1,1,JM+1,LM+1,IM,JS,JE,LM) do l = 1, lm do j = jsta_2l, jend_2u do i = 1, im rswtt ( i, j, l ) = dum3d ( i, j, l ) end do end do end do ! write(6,*) 'RSWTT(II,jj,ll): ', DUM3D(ii,jj,ll) do l = 1, lm do j = jsta_2l, jend_2u do i = 1, im ttnd ( i, j, l ) = rswtt(i,j,l) + rlwtt(i,j,l) end do end do end do VarName='AVRAIN' call getVariable(fileName,DateStr,DataHandle,VarName,AVRAIN, & 1,1,1,1,1,1,1,1) VarName='AVCNVC' call getVariable(fileName,DateStr,DataHandle,VarName,AVCNVC, & 1,1,1,1,1,1,1,1) print*,'AVRAIN,AVCNVC= ',AVRAIN,AVCNVC varname='TCUCN' write(6,*) 'call getVariable for : ', VarName call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D, & IM+1,1,JM+1,LM+1,IM,JS,JE,LM) do l=1,lm do j = jsta_2l, jend_2u do i = 1, im tcucn ( i, j,l ) = dum3d ( i, j,l ) end do end do end do print*,'max tcucn= ',maxval(tcucn) varname='TRAIN' write(6,*) 'call getVariable for : ', VarName call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D, & IM+1,1,JM+1,LM+1,IM,JS,JE,LM) do l=1,lm do j = jsta_2l, jend_2u do i = 1, im train ( i, j, l ) = dum3d ( i, j,l ) end do end do end do print*,'max train= ',maxval(train) VarName='NCFRCV' write(6,*) 'call getIVariable for : ', VarName call getIVariableN(fileName,DateStr,DataHandle,VarName,IDUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im ncfrcv ( i, j ) = float(idummy ( i, j )) ! if(ncfrcv(i,j).gt.1.0e-5)print*,'nonzero ncfrcv',ncfrcv(i,j) end do end do VarName='NCFRST' write(6,*) 'call getIVariable for : ', VarName call getIVariableN(fileName,DateStr,DataHandle,VarName,IDUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im ncfrst ( i, j ) = float(idummy ( i, j )) ! if(ncfrst(i,j).gt.1.0e-5)print*,'nonzero ncfrst',ncfrst(i,j) end do end do VarName='SSROFF' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im SSROFF ( i, j ) = dummy ( i, j ) end do end do VarName='UDROFF' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im BGROFF ( i, j ) = dummy ( i, j ) end do end do VarName='SFCEVP' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im SFCEVP( i, j ) = dummy ( i, j ) end do end do ! print*,'SFCEVP at ',ii,jj,' = ',SFCEVP(ii,jj) VarName='SFCEXC' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im SFCEXC( i, j ) = dummy ( i, j ) end do end do ! print*,'SFCEXC at ',ii,jj,' = ',SFCEXC(ii,jj) VarName='VEGFRC' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im VEGFRC ( i, j ) = dummy ( i, j ) end do end do print*,'VEGFRC at ',ii,jj,' = ',VEGFRC(ii,jj) VarName='ACSNOW' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im ACSNOW ( i, j ) = dummy ( i, j ) end do end do VarName='ACSNOM' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im ACSNOM ( i, j ) = dummy ( i, j ) end do end do ! VarName='CANWAT' VarName='CMC' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im CMC ( i, j ) = dummy ( i, j ) end do end do VarName='SST' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im SST ( i, j ) = dummy ( i, j ) end do end do print*,'SST at ',ii,jj,' = ',sst(ii,jj) VarName='TAUX' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im MDLTAUX ( i, j ) = dummy ( i, j ) end do end do print*,'MDLTAUX at ',ii,jj,' = ',mdltaux(ii,jj) VarName='TAUY' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im MDLTAUY ( i, j ) = dummy ( i, j ) end do end do print*,'MDLTAUY at ',ii,jj,' = ',mdltauy(ii,jj) VarName='EXCH_H' call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D, & IM+1,1,JM+1,LM+1,IM,JS,JE,LM) do l = 1, lm do j = jsta_2l, jend_2u do i = 1, im EXCH_H ( i, j, l ) = dum3d ( i, j, l ) dummy(i,j)=dum3d ( i, j, l ) end do end do print*,'l, max exch = ',l,maxval(dummy) end do do l=1,lm print*,'sample EXCH_H= ',EXCH_H(ii,jj,l) end do VarName='EL_PBL' call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D, & IM+1,1,JM+1,LM+1,IM,JS,JE,LM) do l = 1, lm do j = jsta_2l, jend_2u do i = 1, im EL_PBL ( i, j, l ) = dum3d ( i, j, l ) dummy(i,j)=dum3d ( i, j, l ) end do end do print*,'l, max EL_PBL = ',l,maxval(dummy) end do VarName='THZ0' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im THZ0 ( i, j ) = dummy ( i, j ) end do end do print*,'THZ0 at ',ii,jj,' = ',THZ0(ii,jj) VarName='QZ0' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im QZ0 ( i, j ) = dummy ( i, j ) end do end do print*,'QZ0 at ',ii,jj,' = ',QZ0(ii,jj) VarName='UZ0' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im UZ0 ( i, j ) = dummy ( i, j ) end do end do print*,'UZ0 at ',ii,jj,' = ',UZ0(ii,jj) VarName='VZ0' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im VZ0 ( i, j ) = dummy ( i, j ) end do end do print*,'VZ0 at ',ii,jj,' = ',VZ0(ii,jj) ! VarName='QSFC' VarName='QS' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im QS ( i, j ) = dummy ( i, j ) ! if(qs(i,j).gt.1.0e-7)print*,'nonzero qsfc' end do end do VarName='Z0' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im Z0 ( i, j ) = dummy ( i, j ) end do end do VarName='PBLH' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im PBLH( i, j ) = dummy ( i, j ) end do end do ! write(6,*) 'PBLH(ii,jj): ', DUMMY(ii,jj) VarName='MIXHT' !PLee (3/07) MIXHT=SPVAL !Init value to detect read failure call getVariable(filename,DateStr,DataHandle,Varname,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im MIXHT( i, j ) = dummy ( i, j ) end do end do VarName='USTAR' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im USTAR( i, j ) = dummy ( i, j ) end do end do print*,'USTAR at ',ii,jj,' = ',USTAR(ii,jj) VarName='AKHS_OUT' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im AKHS ( i, j ) = dummy ( i, j ) end do end do print*,'max akhs= ',maxval(akhs) VarName='AKMS_OUT' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im AKMS ( i, j ) = dummy ( i, j ) end do end do print*,'max akms= ',maxval(akms) ! ! In my version, variable is TSK (skin temp, not skin pot temp) ! !mp call getVariable(fileName,DateStr,DataHandle,'THSK',DUMMY, ! VarName='TSK' VarName='THS' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im THS ( i, j ) = dummy ( i, j ) end do end do print*,'THS at ',ii,jj,' = ',THS(ii,jj) !C !CMP !C !C RAINC is "ACCUMULATED TOTAL CUMULUS PRECIPITATION" !C RAINNC is "ACCUMULATED TOTAL GRID SCALE PRECIPITATION" write(6,*) 'getting RAINC' VarName='PREC' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im ! CUPREC ( i, j ) = dummy ( i, j ) * 0.001 PREC ( i, j ) = dummy ( i, j ) end do end do print*,'PREC at ',ii,jj,' = ',PREC(ii,jj) ! VarName='RAINC' VarName='CUPREC' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im ! CUPREC ( i, j ) = dummy ( i, j ) * 0.001 CUPREC ( i, j ) = dummy ( i, j ) end do end do print*,'CUPREC at ',ii,jj,' = ',CUPREC(ii,jj) write(6,*) 'getting RAINTOTAL' ! VarName='RAINNC' VarName='ACPREC' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im ACPREC( i, j ) = dummy ( i, j ) ANCPRC ( i, j ) = ACPREC(I,J)-CUPREC(I,J) end do end do print*,'ACPREC at ',ii,jj,' = ',ACPREC(ii,jj) print*,'ANCPRC at ',ii,jj,' = ',ANCPRC(ii,jj) ! ! hoping to read instantanous convective precip rate soon, initialize it to spval ! for now VarName='CPRATE' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im CPRATE(I,J)=dummy(i,j) enddo enddo VarName='CUPPT' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im CUPPT ( i, j ) = dummy ( i, j ) end do end do print*,'maxval CUPPT: ', maxval(DUMMY) ! adding land surface precipitation accumulation for Yin's precip assimilation VarName='LSPA' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im LSPA ( i, j ) = dummy ( i, j ) end do end do print*,'maxval LSPA: ', maxval(DUMMY) VarName='CLDEFI' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im CLDEFI ( i, j ) = dummy ( i, j ) end do end do print*,'maxval CLDEFI: ', maxval(DUMMY) ! ! Very confusing story ... ! ! Retrieve htop and hbot => They are named CNVTOP, CNVBOT in the model and ! with HBOTS,HTOPS (shallow conv) and HBOTD,HTOPD (deep conv) represent ! the 3 sets of convective cloud base/top arrays tied to the frequency ! that history files are written. ! ! IN THE *MODEL*, arrays HBOT,HTOP are similar to CNVTOP,CNVBOT but are ! used in radiation and are tied to the frequency of radiation updates. ! ! For historical reasons model arrays CNVTOP,CNVBOT are renamed HBOT,HTOP ! and manipulated throughout the post. ! VarName='HTOP' VarName='CNVTOP' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im HTOP ( i, j ) = float(LM)-dummy(i,j)+1.0 HTOP ( i, j ) = max(1.0,min(HTOP(I,J),float(LM))) end do end do print*,'maxval HTOP: ', maxval(DUMMY) ! VarName='HBOT' VarName='CNVBOT' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im HBOT ( i, j ) = float(LM)-dummy(i,j)+1.0 HBOT ( i, j ) = max(1.0,min(HBOT(I,J),float(LM))) end do end do print*,'maxval HBOT: ', maxval(DUMMY) VarName='HTOPD' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im HTOPD ( i, j ) = float(LM)-dummy(i,j)+1.0 end do end do print*,'maxval HTOPD: ', maxval(DUMMY) VarName='HBOTD' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im HBOTD ( i, j ) = float(LM)-dummy(i,j)+1.0 end do end do print*,'maxval HBOTD: ', maxval(DUMMY) VarName='HTOPS' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im HTOPS ( i, j ) = float(LM)-dummy(i,j)+1.0 end do end do print*,'maxval HTOPS: ', maxval(DUMMY) VarName='HBOTS' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im HBOTS ( i, j ) = float(LM)-dummy(i,j)+1.0 end do end do print*,'maxval HBOTS: ', maxval(DUMMY) VarName='CLDFRA' call getVariable(fileName,DateStr,DataHandle,VarName,DUM3D, & IM+1,1,JM+1,LM+1,IM,JS,JE,LM) do l = 1, lm do j = jsta_2l, jend_2u do i = 1, im CFR ( i, j, l ) = dum3d ( i, j, l ) end do end do end do VarName='SR' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im SR ( i, j ) = dummy(i,j) end do end do print*,'maxval SR: ', maxval(DUMMY) ! call getVariable(fileName,DateStr,DataHandle,'RAINCV',DUMMY, ! & IM,1,JM,1,IM,JS,JE,1) ! do j = jsta_2l, jend_2u ! do i = 1, im ! CUPPT ( i, j ) = dummy ( i, j )* 0.001 ! end do ! end do ! ! VarName='GSW' VarName='RSWIN' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im RSWIN ( i, j ) = dummy ( i, j ) ! if(abs(dummy(i,j)).gt. 0.0)print*,'rswin=',dummy(i,j) end do end do VarName='RSWINC' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im RSWINC ( i, j ) = dummy ( i, j ) ! if(abs(dummy(i,j)).gt. 0.0)print*,'rswin=',dummy(i,j) end do end do ! read in zenith angle VarName='CZEN' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im CZEN ( i, j ) = dummy ( i, j ) ! if(abs(czen(i,j)).gt. 0.0)print*,'czen=',czen(i,j) end do end do VarName='CZMEAN' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im CZMEAN ( i, j ) = dummy ( i, j ) ! if(abs(dummy(i,j)).gt. 0.0)print*,'czmean=',dummy(i,j) end do end do VarName='RSWOUT' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im RSWOUT ( i, j ) = dummy ( i, j ) ! if(abs(dummy(i,j)).gt. 0.0)print*,'rswout=',dummy(i,j) end do end do ! VarName='GLW' VarName='RLWIN' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im RLWIN ( i, j ) = dummy ( i, j ) end do end do VarName='RLWTOA' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im RLWTOA ( i, j ) = dummy ( i, j ) end do end do VarName='SIGT4' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im SIGT4 ( i, j ) = dummy ( i, j ) end do end do VarName='RADOT' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im RADOT ( i, j ) = dummy ( i, j ) end do end do ! accumulated incoming short wave VarName='ASWIN' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im ASWIN ( i, j ) = dummy ( i, j ) end do end do ! accumulated outgoing short wave VarName='ASWOUT' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im ASWOUT ( i, j ) = dummy ( i, j ) ! if(abs(dummy(i,j)).gt. 0.0)print*,'aswout=',dummy(i,j) end do end do ! shortwave accumulation frequency VarName='NRDSW' call getIVariableN(fileName,DateStr,DataHandle,VarName,NRDSW, & 1,1,1,1,1,1,1,1) print*,'NRDSW in INITPOST_NMM=',NRDSW VarName='ARDSW' call getVariable(fileName,DateStr,DataHandle,VarName,ARDSW, & 1,1,1,1,1,1,1,1) print*,'ARDSW ARDLW in INITPOST_NMM=',ARDSW, ARDLW ! accumulated incoming long wave VarName='ALWIN' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im ALWIN ( i, j ) = dummy ( i, j ) end do end do ! accumulated outgoing long wave VarName='ALWOUT' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im ALWOUT ( i, j ) = dummy ( i, j ) end do end do ! longwave accumulation frequency VarName='NRDLW' call getIVariableN(fileName,DateStr,DataHandle,VarName,NRDLW, & 1,1,1,1,1,1,1,1) print*,'NRDLW= ',NRDLW ! longwave accumulation counts VarName='ARDLW' call getVariable(fileName,DateStr,DataHandle,VarName,ARDLW, & 1,1,1,1,1,1,1,1) ! obtain time averaged radition at the top of atmosphere VarName='ALWTOA' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im ALWTOA ( i, j ) = dummy ( i, j ) end do end do VarName='ASWTOA' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im ASWTOA ( i, j ) = dummy ( i, j ) end do end do ! VarName='TMN' ! VarName='TG' VarName='TGROUND' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im TG ( i, j ) = dummy ( i, j ) end do end do VarName='SOILTB' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im SOILTB ( i, j ) = dummy ( i, j ) end do end do ! sensible heat fluxes VarName='TWBS' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im TWBS ( i, j ) = dummy ( i, j ) end do end do ! accumulated sensible heat fluxes ! VarName='HFX' VarName='SFCSHX' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im SFCSHX ( i, j ) = dummy ( i, j ) end do end do ! fluxes accumulation frequency VarName='NSRFC' call getIVariableN(fileName,DateStr,DataHandle,VarName,NSRFC, & 1,1,1,1,1,1,1,1) print*,'NSRFC= ',NSRFC ! fluxes accumulation counts VarName='ASRFC' call getVariable(fileName,DateStr,DataHandle,VarName,ASRFC, & 1,1,1,1,1,1,1,1) ! instantanous latent heat fluxes VarName='QWBS' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im QWBS ( i, j ) = dummy ( i, j ) end do end do ! accumulated latent heat fluxes ! VarName='QFX' VarName='SFCLHX' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im SFCLHX ( i, j ) = dummy ( i, j ) end do end do ! instantanous ground heat fluxes VarName='GRNFLX' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im GRNFLX ( i, j ) = dummy ( i, j ) end do end do ! accumulated ground heat fluxes VarName='SUBSHX' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im SUBSHX ( i, j ) = dummy ( i, j ) end do end do ! accumulated ground heat fluxes VarName='POTEVP' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im POTEVP ( i, j ) = dummy ( i, j ) end do end do ! VarName='SNOWC' ! VarName='SNO' VarName='WEASD' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) ! do j = jsta_2l, jend_2u ! do i = 1, im ! SNO ( i, j ) = dummy ( i, j ) ! end do ! end do VarName='SNO' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im SNO ( i, j ) = dummy ( i, j ) end do end do VarName='SI' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im SI ( i, j ) = dummy ( i, j ) end do end do ! snow cover VarName='PCTSNO' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im PCTSNO ( i, j ) = dummy ( i, j ) if(dummy(i,j) .gt. 1.0e-5)print*,'nonzero pctsno' end do end do ! GET VEGETATION TYPE ! call getVariable(fileName,DateStr,DataHandle,'IVGTYP',DUMMY ! & ,IM,1,JM,1,IM,JS,JE,1) ! print*,'sample VEG TYPE',DUMMY(20,20) ! XLAND 1 land 2 sea ! VarName='XLAND' VarName='IVGTYP' call getIVariableN(fileName,DateStr,DataHandle,VarName,IDUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im IVGTYP ( i, j ) = idummy ( i, j ) end do end do print*,'MAX IVGTYP=', maxval(idummy) VarName='ISLTYP' call getIVariableN(fileName,DateStr,DataHandle,VarName,IDUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im ISLTYP ( i, j ) = idummy ( i, j ) end do end do print*,'MAX ISLTYP=', maxval(idummy) VarName='ISLOPE' call getIVariableN(fileName,DateStr,DataHandle,VarName,IDUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im ISLOPE( i, j ) = idummy ( i, j ) end do end do VarName='SM' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im SM ( i, j ) = dummy ( i, j ) end do end do VarName='SICE' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im SICE ( i, j ) = dummy ( i, j ) end do end do VarName='ALBEDO' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im ALBEDO( i, j ) = dummy ( i, j ) end do end do VarName='ALBASE' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im ALBASE( i, j ) = dummy ( i, j ) end do end do VarName='MXSNAL' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im MXSNAL( i, j ) = dummy ( i, j ) end do end do VarName='EPSR' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im EPSR( i, j ) = dummy ( i, j ) end do end do ! VarName='XLAT' VarName='GLAT' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im f(i,j) = 1.454441e-4*sin(dummy(i,j)) GDLAT ( i, j ) = dummy ( i, j ) * RTD end do end do ! pos north print*,'GDLAT at ',ii,jj,' = ',GDLAT(ii,jj) print*,'read past GDLAT' ! VarName='XLONG' VarName='GLON' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im GDLON ( i, j ) = dummy ( i, j ) * RTD ! if(j.eq.1 .or. j.eq.jm)print*,'I,J,GDLON,GDLAT= ',i,j ! 1 ,GDLON( i, j ),GDLAT ( i, j ) ! if(abs(GDLAT(i,j)-20.0).lt.0.5 .and. abs(GDLON(I,J) ! 1 +157.0).lt.5.)print* ! 2 ,'Debug:I,J,GDLON,GDLAT,SM,HGT,psfc= ',i,j,GDLON(i,j) ! 3 ,GDLAT(i,j),SM(i,j),FIS(i,j)/G,PINT(I,j,lm+1) end do end do print*,'GDLON at ',ii,jj,' = ',GDLON(ii,jj) print*,'read past GDLON' ! pos east call collect_loc(gdlat,dummy) get_dcenlat: if(me.eq.0)then latstart=nint(dummy(1,1)*1000.) ! lower left latlast=nint(dummy(im,jm)*1000.) ! upper right icen=im/2 !center grid jcen=jm/2 print *, 'dummy(icen,jcen) = ', dummy(icen,jcen) print *, 'dummy(icen-1,jcen) = ', dummy(icen-1,jcen) print *, 'dummy(icen+1,jcen) = ', dummy(icen+1,jcen) ! Grid navigation for copygb - R.Rozumalski latnm = nint(dummy(icen,jm)*1000.) latsm = nint(dummy(icen,1)*1000.) print *, 'latnm, latsm', latnm, latsm ! temporary patch for nmm wrf for moving nest ! cenlat = glat(im/2,jm/2) -Gopal if(mod(im,2).ne.0)then !per Pyle, jm is always odd if(mod(jm+1,4).ne.0)then dcenlat=dummy(icen,jcen) else dcenlat=0.5*(dummy(icen-1,jcen)+dummy(icen,jcen)) end if else if(mod(jm+1,4).ne.0)then dcenlat=0.5*(dummy(icen,jcen)+dummy(icen+1,jcen)) else dcenlat=dummy(icen,jcen) end if end if endif get_dcenlat write(6,*) 'laststart,latlast,dcenlat B calling bcast= ', & latstart,latlast,dcenlat call mpi_bcast(latstart,1,MPI_INTEGER,0,mpi_comm_comp,irtn) call mpi_bcast(latlast,1,MPI_INTEGER,0,mpi_comm_comp,irtn) call mpi_bcast(dcenlat,1,MPI_REAL,0,mpi_comm_comp,irtn) write(6,*) 'laststart,latlast A calling bcast= ',latstart,latlast call collect_loc(gdlon,dummy) get_dcenlon: if(me.eq.0)then lonstart=nint(dummy(1,1)*1000.) lonlast=nint(dummy(im,jm)*1000.) ! icen, jcen set above print *, 'lon dummy(icen,jcen) = ', dummy(icen,jcen) print *, 'lon dummy(icen-1,jcen) = ', dummy(icen-1,jcen) print *, 'lon dummy(icen+1,jcen) = ', dummy(icen+1,jcen) ! Grid navigation for copygb - R.Rozumalski lonem = nint(dummy(icen,jm)*1000.) lonwm = nint(dummy(icen,1)*1000.) if(mod(im,2).ne.0)then !per Pyle, jm is always odd if(mod(jm+1,4).ne.0)then cen1=dummy(icen,jcen) cen2=cen1 else cen1=min(dummy(icen-1,jcen),dummy(icen,jcen)) cen2=max(dummy(icen-1,jcen),dummy(icen,jcen)) end if else if(mod(jm+1,4).ne.0)then cen1=min(dummy(icen+1,jcen),dummy(icen,jcen)) cen2=max(dummy(icen+1,jcen),dummy(icen,jcen)) else cen1=dummy(icen,jcen) cen2=cen1 end if end if ! Trahan fix: Pyle's code broke at the dateline. if(cen2-cen1>180) then ! We're near the dateline dcenlon=mod(0.5*(cen2+cen1+360)+3600+180,360.)-180. else ! We're not near the dateline. Use the original code, ! unmodified, to maintain bitwise identicality. dcenlon=0.5*(cen1+cen2) endif end if get_dcenlon ! rank 0 write(6,*)'lonstart,lonlast,cenlon B calling bcast= ',lonstart, & lonlast,cenlon call mpi_bcast(lonstart,1,MPI_INTEGER,0,mpi_comm_comp,irtn) call mpi_bcast(lonlast,1,MPI_INTEGER,0,mpi_comm_comp,irtn) call mpi_bcast(dcenlon,1,MPI_REAL,0,mpi_comm_comp,irtn) write(6,*)'lonstart,lonlast,cenlon A calling bcast= ',lonstart, & lonlast,cenlon if(me==0) then open(1013,file='this-domain-center.ksh.inc',form='formatted',status='unknown') 1013 format(A,'=',F0.3) write(1013,1013) 'clat',dcenlat write(1013,1013) 'clon',dcenlon endif ! ! OBTAIN DX FOR NMM WRF VarName='DX_NMM' call getVariable(fileName,DateStr,DataHandle,VarName,DUMMY, & IM,1,JM,1,IM,JS,JE,1) do j = jsta_2l, jend_2u do i = 1, im DX ( i, j ) = dummy ( i, j ) if(DX(i,j)<0.1)print*,'zero dx in INIT: I,J,DX= ',i,j & ,DX( i, j ) ! if(j.eq.1 .or. j.eq.jm)print*,'I,J,DX= ',i,j ! 1 ,DX( i, j ) end do end do varname='ETA1' write(6,*) 'call getVariable for : ', VarName call getVariable(fileName,DateStr,DataHandle,VarName,ETA1, & LM,1,1,1,LM,1,1,1) varname='ETA2' write(6,*) 'call getVariable for : ', VarName call getVariable(fileName,DateStr,DataHandle,VarName,ETA2, & LM,1,1,1,LM,1,1,1) open(75,file='ETAPROFILE.txt',form='formatted',status='unknown') DO L=1,lm+1 IF(L .EQ. 1)THEN write(75,1020)L, 0., 0. ELSE write(75,1020)L, ETA1(lm+2-l), ETA2(lm+2-l) END IF ! print*,'L, ETA1, ETA2= ',L, ETA1(l), ETA2(l) END DO 1020 format(I3,2E17.10) close (75) ! physics calling frequency VarName='NPHS0' call getIVariableN(fileName,DateStr,DataHandle,VarName,NPHS, & 1,1,1,1,1,1,1,1) print*,'NPHS= ',NPHS ! physics calling frequency VarName='NCLOD' call getIVariableN(fileName,DateStr,DataHandle,VarName,NCLOD, & 1,1,1,1,1,1,1,1) ! physics calling frequency VarName='NPREC' call getIVariableN(fileName,DateStr,DataHandle,VarName,NPREC, & 1,1,1,1,1,1,1,1) ! physics calling frequency VarName='NHEAT' call getIVariableN(fileName,DateStr,DataHandle,VarName,NHEAT, & 1,1,1,1,1,1,1,1) print*,'NHEAT= ',NHEAT ! ! ncdump -h !! !! !! write(6,*) 'filename in INITPOST=', filename,' is' ! status=nf_open(filename,NF_NOWRITE,ncid) ! write(6,*) 'returned ncid= ', ncid ! status=nf_get_att_real(ncid,varid,'DX',tmp) ! dxval=int(tmp) ! status=nf_get_att_real(ncid,varid,'DY',tmp) ! dyval=int(tmp) ! status=nf_get_att_real(ncid,varid,'CEN_LAT',tmp) ! cenlat=int(1000.*tmp) ! status=nf_get_att_real(ncid,varid,'CEN_LON',tmp) ! cenlon=int(1000.*tmp) ! status=nf_get_att_real(ncid,varid,'TRUELAT1',tmp) ! truelat1=int(1000.*tmp) ! status=nf_get_att_real(ncid,varid,'TRUELAT2',tmp) ! truelat2=int(1000.*tmp) ! status=nf_get_att_real(ncid,varid,'MAP_PROJ',tmp) ! maptype=int(tmp) ! status=nf_close(ncid) ! dxval=30000. ! dyval=30000. ! ! write(6,*) 'dxval= ', dxval ! write(6,*) 'dyval= ', dyval ! write(6,*) 'cenlat= ', cenlat ! write(6,*) 'cenlon= ', cenlon ! write(6,*) 'truelat1= ', truelat1 ! write(6,*) 'truelat2= ', truelat2 ! write(6,*) 'maptype is ', maptype ! call ext_ncd_get_dom_ti_real(DataHandle,'DX',tmp, & 1,ioutcount,istatus) dxval=nint(tmp*1000.) ! E-grid dlamda in degree write(6,*) 'dxval= ', dxval call ext_ncd_get_dom_ti_real(DataHandle,'DY',tmp, & 1,ioutcount,istatus) dyval=nint(tmp*1000.) write(6,*) 'dyval= ', dyval call ext_ncd_get_dom_ti_real(DataHandle,'CEN_LAT',tmp, & 1,ioutcount,istatus) cenlat=nint(tmp*1000.) ! E-grid dlamda in degree write(6,*) 'cenlat= ', cenlat call ext_ncd_get_dom_ti_real(DataHandle,'CEN_LON',tmp, & 1,ioutcount,istatus) cenlon=nint(tmp*1000.) ! E-grid dlamda in degree write(6,*) 'cenlon= ', cenlon ! JW call ext_ncd_get_dom_ti_real(DataHandle,'TRUELAT1',tmp ! JW + ,1,ioutcount,istatus) ! JW truelat1=nint(1000.*tmp) ! JW write(6,*) 'truelat1= ', truelat1 ! JW call ext_ncd_get_dom_ti_real(DataHandle,'TRUELAT2',tmp ! JW + ,1,ioutcount,istatus) ! JW truelat2=nint(1000.*tmp) ! JW write(6,*) 'truelat2= ', truelat2 call ext_ncd_get_dom_ti_integer(DataHandle,'MAP_PROJ',itmp, & 1,ioutcount,istatus) maptype=itmp gridtype = 'E' write(6,*) 'maptype, gridtype ', maptype, gridtype gridtype='E' call ext_ncd_get_dom_ti_integer(DataHandle,'I_PARENT_START',itmp, & 1,ioutcount,istatus) i_parent_start=itmp call ext_ncd_get_dom_ti_integer(DataHandle,'J_PARENT_START',itmp, & 1,ioutcount,istatus) j_parent_start=itmp do j = jsta_2l, jend_2u do i = 1, im ! DX ( i, j ) = dxval DY ( i, j ) = dyval*DTR*ERAD*0.001 end do end do ! generate look up table for lifted parcel calculations THL=210. PLQ=70000. CALL TABLE(PTBL,TTBL,PT, & RDQ,RDTH,RDP,RDTHE,PL,THL,QS0,SQS,STHE,THE0) CALL TABLEQ(TTBLQ,RDPQ,RDTHEQ,PLQ,THL,STHEQ,THE0Q) ! ! IF(ME.EQ.0)THEN WRITE(6,*)' SPL (POSTED PRESSURE LEVELS) BELOW: ' WRITE(6,51) (SPL(L),L=1,LSM) 50 FORMAT(14(F4.1,1X)) 51 FORMAT(8(F8.1,1X)) ENDIF ! ! COMPUTE DERIVED TIME STEPPING CONSTANTS. ! call ext_ncd_get_dom_ti_real(DataHandle,'DT',tmp, & 1,ioutcount,istatus) DT=tmp print*,'DT= ',DT DTQ2 = DT * NPHS TSPH = 3600./DT TSRFC=float(NSRFC)/TSPH IF(NSRFC.EQ.0)TSRFC=float(ifhr) !in case buket does not get emptied TRDLW=float(NRDLW)/TSPH IF(NRDLW.EQ.0)TRDLW=float(ifhr) !in case buket does not get emptied TRDSW=float(NRDSW)/TSPH IF(NRDSW.EQ.0)TRDSW=float(ifhr) !in case buket does not get emptied THEAT=float(NHEAT)/TSPH IF(NHEAT.EQ.0)THEAT=float(ifhr) !in case buket does not get emptied TCLOD=float(NCLOD)/TSPH IF(NCLOD.EQ.0)TCLOD=float(ifhr) !in case buket does not get emptied TPREC=float(NPREC)/TSPH IF(NPREC.EQ.0)TPREC=float(ifhr) !in case buket does not get emptied print*,'TSRFC TRDLW TRDSW= ',TSRFC, TRDLW, TRDSW !how am i going to get this information? ! NPREC = INT(TPREC *TSPH+D50) ! NHEAT = INT(THEAT *TSPH+D50) ! NCLOD = INT(TCLOD *TSPH+D50) ! NRDSW = INT(TRDSW *TSPH+D50) ! NRDLW = INT(TRDLW *TSPH+D50) ! NSRFC = INT(TSRFC *TSPH+D50) !how am i going to get this information? ! ! IF(ME.EQ.0)THEN ! WRITE(6,*)' ' ! WRITE(6,*)'DERIVED TIME STEPPING CONSTANTS' ! WRITE(6,*)' NPREC,NHEAT,NSRFC : ',NPREC,NHEAT,NSRFC ! WRITE(6,*)' NCLOD,NRDSW,NRDLW : ',NCLOD,NRDSW,NRDLW ! ENDIF ! ! COMPUTE DERIVED MAP OUTPUT CONSTANTS. DO L = 1,LSM ALSL(L) = ALOG(SPL(L)) END DO ! if(submodelname == 'NEST') then print *,'NMM NEST mode: use projection center as projection center' else print *,'NMM MOAD mode: use domain center as projection center' CENLAT=NINT(DCENLAT*1000) CENLON=NINT(DCENLON*1000) endif if(me.eq.0)then ! write out copygb_gridnav.txt ! provided by R.Rozumalski - NWS inav=10 TRUELAT1 = CENLAT TRUELAT2 = CENLAT IFDX = NINT (dxval*107.) IFDY = NINT (dyval*110.) open(inav,file='copygb_gridnav.txt',form='formatted', & status='unknown') print *, ' MAPTYPE :',maptype print *, ' IM :',IM*2-1 print *, ' JM :',JM print *, ' LATSTART :',LATSTART print *, ' LONSTART :',LONSTART print *, ' CENLAT :',CENLAT print *, ' CENLON :',CENLON print *, ' TRUELAT2 :',TRUELAT2 print *, ' TRUELAT1 :',TRUELAT1 print *, ' DX :',IFDX*0.001 print *, ' DY :',IFDY*0.001 IF(MAPTYPE.EQ.0 .OR. MAPTYPE.EQ.203)THEN !A STAGGERED E-GRID IMM = 2*IM-1 IDXAVE = ( IFDY + IFDX ) * 0.5 ! If the Center Latitude of the domain is located within 15 degrees ! of the equator then use a a regular Lat/Lon navigation for the ! remapped grid in copygb; otherwise, use a Lambert conformal. Make ! sure to specify the correct pole for the S. Hemisphere (LCC). ! IF ( abs(CENLAT).GT.15000) THEN write(6,*)' Copygb LCC Navigation Information' IF (CENLAT .GT.0) THEN ! Northern Hemisphere write(6,1000) IMM,JM,LATSTART,LONSTART,CENLON, & IFDX,IFDY,CENLAT,CENLAT write(inav,1000) IMM,JM,LATSTART,LONSTART,CENLON, & IFDX,IFDY,CENLAT,CENLAT ELSE ! Southern Hemisphere write(6,1001) IMM,JM,LATSTART,LONSTART,CENLON, & IFDX,IFDY,CENLAT,CENLAT write(inav,1001) IMM,JM,LATSTART,LONSTART,CENLON, & IFDX,IFDY,CENLAT,CENLAT END IF ELSE dlat = (latnm-latsm)/(JM-1) nlat = INT (dlat) if (lonem .lt. 0) lonem = 360000. + lonem if (lonwm .lt. 0) lonwm = 360000. + lonwm dlon = lonem-lonwm if (dlon .lt. 0.) dlon = dlon + 360000. dlon = (dlon)/(IMM-1) nlon = INT (dlon) write(6,*)' Copygb Lat/Lon Navigation Information' write(6,2000) IMM,JM,latsm,lonwm,latnm,lonem,nlon,nlat write(inav,2000) IMM,JM,latsm,lonwm,latnm,lonem,nlon,nlat ENDIF close(inav) 1000 format('255 3 ',2(I3,x),I6,x,I7,x,'8 ',I7,x,2(I6,x),'0 64', & 2(x,I6)) 1001 format('255 3 ',2(I3,x),I6,x,I7,x,'8 ',I7,x,2(I6,x),'128 64', & 2(x,I6),' -90000 0') 2000 format('255 0 ',2(I3,x),2(I7,x),'8 ',2(I7,x),2(I7,x),'64') END IF ! maptype !HC WRITE IGDS OUT FOR WEIGHTMAKER TO READ IN AS KGDSIN igdout=110 if (maptype .eq. 1)THEN ! Lambert conformal WRITE(igdout)3 WRITE(6,*)'igd(1)=',3 WRITE(igdout)im WRITE(igdout)jm WRITE(igdout)LATSTART WRITE(igdout)LONSTART WRITE(igdout)8 WRITE(igdout)CENLON WRITE(igdout)DXVAL WRITE(igdout)DYVAL WRITE(igdout)0 WRITE(igdout)64 ! JW WRITE(igdout)TRUELAT2 ! JW WRITE(igdout)TRUELAT1 WRITE(igdout)255 ELSE IF(MAPTYPE .EQ. 2)THEN !Polar stereographic WRITE(igdout)5 WRITE(igdout)im WRITE(igdout)jm WRITE(igdout)LATSTART WRITE(igdout)LONSTART WRITE(igdout)8 WRITE(igdout)CENLON WRITE(igdout)DXVAL WRITE(igdout)DYVAL WRITE(igdout)0 WRITE(igdout)64 ! JW WRITE(igdout)TRUELAT2 !Assume projection at +-90 ! JW WRITE(igdout)TRUELAT1 WRITE(igdout)255 ! Note: The calculation of the map scale factor at the standard ! lat/lon and the PSMAPF ! Get map factor at 60 degrees (N or S) for PS projection, which will ! be needed to correctly define the DX and DY values in the GRIB GDS if (TRUELAT1 .LT. 0.) THEN LAT = -60. else LAT = 60. end if CALL MSFPS (LAT,TRUELAT1*0.001,PSMAPF) ELSE IF(MAPTYPE .EQ. 3)THEN !Mercator WRITE(igdout)1 WRITE(igdout)im WRITE(igdout)jm WRITE(igdout)LATSTART WRITE(igdout)LONSTART WRITE(igdout)8 WRITE(igdout)latlast WRITE(igdout)lonlast ! JW WRITE(igdout)TRUELAT1 WRITE(igdout)0 WRITE(igdout)64 WRITE(igdout)DXVAL WRITE(igdout)DYVAL WRITE(igdout)255 ELSE IF(MAPTYPE.EQ.0 .OR. MAPTYPE.EQ.203)THEN !A STAGGERED E-GRID WRITE(igdout)203 WRITE(igdout)im WRITE(igdout)jm WRITE(igdout)LATSTART WRITE(igdout)LONSTART WRITE(igdout)136 WRITE(igdout)CENLAT WRITE(igdout)CENLON WRITE(igdout)DXVAL WRITE(igdout)DYVAL WRITE(igdout)64 WRITE(igdout)0 WRITE(igdout)0 WRITE(igdout)0 ! following for hurricane wrf post open(inav,file='copygb_hwrf.txt',form='formatted', & status='unknown') LATEND=LATSTART+(JM-1)*dyval LONEND=LONSTART+(IMM-1)*dxval write(10,1010) IMM,JM,LATSTART,LONSTART,LATEND,LONEND, & dxval,dyval 1010 format('255 0 ',2(I3,x),I6,x,I7,x,'136 ',I6,x,I7,x, & 2(I6,x),'64') close (inav) END IF end if ! ! ! close up shop call ext_ncd_ioclose ( DataHandle, Status ) RETURN END