subroutine print_diags(ipos) 3,17
!@sum print_diag prints out binary and ascii diag output.
!@auth Original Development Team
USE MODEL_COM
, only : itime,itimeI
USE DAGCOM, only : kdiag,keynr,keyct,isccp_diags
IMPLICIT NONE
!@var ipos =1 (after input), =2 (current diags), =3 (end of diag period)
INTEGER, INTENT(IN) :: ipos
IF (KDIAG(1).LT.9) CALL DIAGJ
IF (KDIAG(2).LT.9) CALL DIAGJK
IF (KDIAG(10).LT.9) CALL DIAGIL
IF (KDIAG(7).LT.9) CALL DIAG7P
IF (KDIAG(3).LT.9) CALL DIAGIJ
IF (KDIAG(9).LT.9) CALL DIAGCP
IF (KDIAG(5).LT.9) CALL DIAG5P
IF (ipos.ne.2. .and. KDIAG(6).LT.9) CALL DIAGDD
IF (KDIAG(13).LT.9) CALL DIAGDH
IF (KDIAG(4).LT.9) CALL DIAG4
IF (KDIAG(11).LT.9) CALL diag_RIVER
IF (KDIAG(12).LT.9) CALL diag_OCEAN
IF (KDIAG(12).LT.9) CALL diag_ICEDYN
IF (isccp_diags.eq.1) CALL diag_ISCCP
IF (ipos.ne.2 .or. Itime.LE.ItimeI+1) THEN
CALL DIAGKN
ELSE ! RESET THE UNUSED KEYNUMBERS TO ZERO
KEYNR(1:42,KEYCT)=0
END IF
return
end subroutine print_diags
MODULE BDJ 2
!@sum stores information for outputting composite zonal diagnostics
!@auth M. Kelley
IMPLICIT NONE
SAVE
!@param nj_out number of j-format output fields = 11
integer, parameter :: nj_out=11
!@var units string containing output field units
CHARACTER(LEN=50), DIMENSION(nj_out) :: UNITS_J_O
!@var lname string describing output field
CHARACTER(LEN=50), DIMENSION(nj_out) :: LNAME_J_O
!@var sname string referencing output field in self-desc. output file
CHARACTER(LEN=30), DIMENSION(nj_out) :: NAME_J_O
!@var stitle short title for print out
CHARACTER(LEN=16), DIMENSION(nj_out) :: STITLE_J_O
!@var INUM_J_O,IDEN_J_O numerator and denominator for calculated J diags
INTEGER, DIMENSION(nj_out) :: INUM_J_O, IDEN_J_O
!@var SCALE_J_O scale for calculated J diags
REAL*8, DIMENSION(nj_out) :: SCALE_J_O
END MODULE BDJ
SUBROUTINE J_TITLES 1,2
!@sum J_TITLES calculated zonal diagnostics
!@auth M. Kelley/G. Schmidt
!@ver 1.0
USE DAGCOM, only : j_srincp0,j_srnfp0,j_plavis,j_planir,j_srnfg
* ,j_srincg,j_albvis,j_albnir,j_srrvis,j_srrnir,j_sravis
* ,j_sranir,j_clddep,j_pcldmc
USE BDJ
IMPLICIT NONE
INTEGER :: K
C**** These information are for J zonal budget calulated diagnostics
C**** Note that we assume that they are all ratios of two existing
C**** records.
c
k = 0
c
k = k + 1
name_j_o(k) = 'plan_alb'
lname_j_o(k) = ' TOTAL PLANETARY ALBEDO'
units_j_o(k) = '%'
stitle_j_o(k)= ' PLANETARY ALBDO'
inum_j_o(k) = J_SRNFP0
iden_j_o(k) = J_SRINCP0
scale_j_o(k) = 100.
c
k = k + 1
name_j_o(k) = 'plan_alb_vis'
lname_j_o(k) = 'PLANETARY ALBEDO IN VISUAL'
units_j_o(k) = '%'
stitle_j_o(k)= ' PLAN ALB VISUAL'
inum_j_o(k) = J_PLAVIS
iden_j_o(k) = J_SRINCP0
scale_j_o(k) = 100.
c
k = k + 1
name_j_o(k) = 'plan_alb_nir'
lname_j_o(k) = 'PLANETARY ALBEDO IN NEAR IR'
units_j_o(k) = '%'
stitle_j_o(k)= ' PLAN ALB NEARIR'
inum_j_o(k) = J_PLANIR
iden_j_o(k) = J_SRINCP0
scale_j_o(k) = 100.
c
k = k + 1
name_j_o(k) = 'surf_alb'
lname_j_o(k) = 'GROUND ALBEDO'
units_j_o(k) = '%'
stitle_j_o(k)= ' SURFACE G ALBDO'
inum_j_o(k) = J_SRNFG
iden_j_o(k) = J_SRINCG
scale_j_o(k) = 100.
c
k = k + 1
name_j_o(k) = 'surf_alb_vis'
lname_j_o(k) = 'GROUND ALBEDO IN VISUAL'
units_j_o(k) = '%'
stitle_j_o(k)= ' SURF ALB VISUAL'
inum_j_o(k) = J_ALBVIS
iden_j_o(k) = J_SRINCP0
scale_j_o(k) = 100.
c
k = k + 1
name_j_o(k) = 'surf_alb_nir'
lname_j_o(k) = 'GROUND ALBEDO IN NEAR IR'
units_j_o(k) = '%'
stitle_j_o(k)= ' SURF ALB NEARIR'
inum_j_o(k) = J_ALBNIR
iden_j_o(k) = J_SRINCP0
scale_j_o(k) = 100.
c
k = k + 1
name_j_o(k) = 'atm_alb_vis'
lname_j_o(k) = 'ATMOSPHERIC ALBEDO IN VISUAL'
units_j_o(k) = '%'
stitle_j_o(k)= '0ATMO ALB VISUAL'
inum_j_o(k) = J_SRRVIS
iden_j_o(k) = J_SRINCP0
scale_j_o(k) = 100.
c
k = k + 1
name_j_o(k) = 'atm_alb_nir'
lname_j_o(k) = 'ATMOSPHERIC ALBEDO IN NEAR IR'
units_j_o(k) = '%'
stitle_j_o(k)= ' ATMO ALB NEARIR'
inum_j_o(k) = J_SRRNIR
iden_j_o(k) = J_SRINCP0
scale_j_o(k) = 100.
c
k = k + 1
name_j_o(k) = 'atm_abs_vis'
lname_j_o(k) = 'ATMOSPHERIC ABSORPTION IN VISUAL'
units_j_o(k) = 'W/m**2'
stitle_j_o(k)= ' ATMO ABS VISUAL'
inum_j_o(k) = J_SRAVIS
iden_j_o(k) = J_SRINCP0
scale_j_o(k) = 100.
c
k = k + 1
name_j_o(k) = 'atm_abs_nir'
lname_j_o(k) = 'ATMOSPHERIC ABSORPTION IN NEAR IR'
units_j_o(k) = 'W/m**2'
stitle_j_o(k)= ' ATMO ABS NEARIR'
inum_j_o(k) = J_SRANIR
iden_j_o(k) = J_SRINCP0
scale_j_o(k) = 100.
c
k = k + 1
name_j_o(k) = 'mc_clddp'
lname_j_o(k) = 'MOIST CONVECTIVE CLOUD DEPTH'
units_j_o(k) = 'mb'
stitle_j_o(k)= ' MC CLD DPTH(MB)'
inum_j_o(k) = J_CLDDEP
iden_j_o(k) = J_PCLDMC
scale_j_o(k) = 1.
RETURN
END SUBROUTINE J_TITLES
SUBROUTINE DIAGJ 1,15
!@sum DIAGJ produces area weighted statistics of zonal budget diags
!@+ based on settings and quantities found in j_defs
!@auth G. Schmidt/R. Reto/G. Russell
use filemanager
USE CONSTANT, only : teeny
USE MODEL_COM, only : im,jm,lm,fim,flice,
& dtsrc,fland,idacc,jhour,jhour0,jdate,jdate0,amon,amon0,
& jyear,jyear0,ls1,sige,itime,itime0,nday,xlabel,lrunid,ntype
USE GEOM, only : dxyp,lat,lat_dg
USE DAGCOM, only :
& QDIAG,acc_period,aj,areg,jreg,kdiag,namreg,nreg,kaj,ia_j,
& j_srabs,j_srnfp0,j_srnfg,j_trnfp0,j_hsurf,j_trhdt,j_trnfp1,
* j_hatm,j_rnfp0,j_rnfp1,j_srnfp1,j_rhdt,j_hz1,j_prcp,j_prcpss,
* j_prcpmc,j_hz0,j_hmelt,j_implh,j_shdt,j_evhdt,j_eprcp,j_erun,
* j_hz2,j_type,j_ervr,scale_j,stitle_j,lname_j,name_j,units_j,
* k_j_out,ia_srf,ia_src,ia_rad,j_h2och4
USE BDJ
IMPLICIT NONE
REAL*8, DIMENSION(JM), SAVE ::S1
REAL*8, DIMENSION(NREG), SAVE :: SAREA
REAL*8, DIMENSION(JM) :: FLAT
REAL*8, DIMENSION(NTYPE,JM) :: SPTYPE
REAL*8, DIMENSION(2) :: FHEM
INTEGER, DIMENSION(JM) :: MLAT
LOGICAL QALB,qIbp
INTEGER, PARAMETER :: INC=1+(JM-1)/24,JMHALF=JM/2
!@param NTYPE_OUT number of output budgets pages
INTEGER, PARAMETER :: NTYPE_OUT=NTYPE+3 ! to include comp/regio
C**** Expanded version of surfaces (including composites)
!@var TERRAIN name of surface type
CHARACTER*16, DIMENSION(0:NTYPE_OUT), PARAMETER :: TERRAIN = (/
* ' (GLOBAL)','(OPEN OCEAN)',' (OCEAN ICE)',' (OCEAN)',
* ' (LAND)',' (LAND ICE)',' (OPEN LAKE)',' (LAKE ICE)',
* ' (LAKES)',' (REGIONS)'/)
C**** Arrays needed for full output
REAL*8, DIMENSION(JM+3,KAJ) :: BUDG
CHARACTER*16, DIMENSION(KAJ) :: TITLEO
CHARACTER*50, DIMENSION(KAJ) :: LNAMEO
CHARACTER*30, DIMENSION(KAJ) :: SNAMEO
CHARACTER*50, DIMENSION(KAJ) :: UNITSO
C**** weighting functions for surface types
REAL*8, DIMENSION(0:NTYPE_OUT-1,NTYPE), PARAMETER ::
* WT=RESHAPE( ! separate types + composites
* (/1.,1.,0.,1.,0.,0.,0.,0.,0., 1.,0.,1.,1.,0.,0.,0.,0.,0.,
* 1.,0.,0.,0.,1.,0.,0.,0.,0., 1.,0.,0.,0.,0.,1.,0.,0.,0.,
* 1.,0.,0.,0.,0.,0.,1.,0.,1., 1.,0.,0.,0.,0.,0.,0.,1.,1./),
* (/NTYPE_OUT,NTYPE/) )
!@var DERPOS character array that determines where derived arrays go
!@var NDERN how many of the derived arrays go in
!@var NDMAX max number of derived array place holders
INTEGER, PARAMETER :: NDMAX=2
INTEGER, DIMENSION(NDMAX), PARAMETER :: ! currently only 2 points
* NDERN = (/10, 1/) ! 10 rad/alb diags and 1 cld diag
CHARACTER*20, DIMENSION(NDMAX), PARAMETER ::
* DERPOS = (/'inc_sw','totcld'/)
REAL*8 :: A1BYA2,A2BYA1,BYA1,BYIACC,FGLOB,GSUM,GSUM2,GWT
* ,HSUM,HSUM2,HWT,QDEN,QJ,QNUM,DAYS,WTX
INTEGER :: I,IACC,J,JH,JHEMI,JR,K,KA,M,MD,N,ND,NN,IT,NDER,KDER
* ,iu_Ibp
character*80 line
logical, save :: Qbp(0:NTYPE_OUT)
INTEGER, SAVE :: IFIRST = 1
integer, external :: NINTlimit
IF (IFIRST.EQ.1) THEN
IFIRST=0
C**** INITIALIZE CERTAIN QUANTITIES
call j_titles
SAREA=0.
DO J=1,JM
S1(J)=IM
DO I=1,IM
JR=JREG(I,J)
SAREA(JR)=SAREA(JR)+DXYP(J)
END DO
END DO
S1(1)=1.
S1(JM)=1.
inquire(file='Ibp',exist=qIbp)
Qbp=.true.
if(.not.qIbp) then
call openunit('Ibp',iu_Ibp,.false.,.false.)
write (iu_Ibp,'(a)') 'List of budget-pages'
do m = 0,ntype_out
write (iu_Ibp,'(i3,1x,a)') m,terrain(m)
end do
else if(kdiag(1).gt.0) then
Qbp=.false.
call openunit('Ibp',iu_Ibp,.false.,.true.)
read (iu_Ibp,'(a)',end=20) line
10 read (iu_Ibp,'(a)',end=20) line
read(line,'(i3)') m
Qbp(m)=.true.
go to 10
20 continue
end if
END IF
C**** OPEN PLOTTABLE OUTPUT FILE IF DESIRED
IF (QDIAG) ! excl. regions, but types dimensioned 0:ntype_out
& call open_j(trim(acc_period)//'.j'//XLABEL(1:LRUNID)
* ,ntype_out,jm,lat_dg)
C**** CALCULATE THE DERIVED QUANTTIES
BYA1=1./(IDACC(ia_srf)+teeny)
A2BYA1=FLOAT(IDACC(ia_rad))/FLOAT(IDACC(ia_src))
A1BYA2=IDACC(ia_src)/(IDACC(ia_rad)+teeny)
DO JR=1,23 ! only 23 will fit on a green sheet
AREG(JR,J_TRNFP0)=AREG(JR,J_HSURF)+A2BYA1*AREG(JR,J_TRHDT)/DTSRC
AREG(JR,J_TRNFP1)=AREG(JR,J_HATM)+A2BYA1*AREG(JR,J_TRHDT)/DTSRC
AREG(JR,J_SRABS) =AREG(JR,J_SRNFP0)-AREG(JR,J_SRNFG)
AREG(JR,J_RNFP0) =AREG(JR,J_SRNFP0)+AREG(JR,J_TRNFP0)
AREG(JR,J_RNFP1) =AREG(JR,J_SRNFP1)+AREG(JR,J_TRNFP1)
AREG(JR,J_RHDT) =A1BYA2*AREG(JR,J_SRNFG)*DTSRC+AREG(JR,J_TRHDT)
AREG(JR,J_PRCP) =AREG(JR,J_PRCPSS)+AREG(JR,J_PRCPMC)
AREG(JR,J_HZ0)=AREG(JR,J_RHDT)+AREG(JR,J_SHDT)+
* AREG(JR,J_EVHDT)+AREG(JR,J_EPRCP)
AREG(JR,J_HZ1)=AREG(JR,J_HZ0)+AREG(JR,J_ERVR)
AREG(JR,J_HZ2)=AREG(JR,J_HZ1)-AREG(JR,J_ERUN)-AREG(JR,J_IMPLH)
END DO
DO J=1,JM
DO IT=1,NTYPE
SPTYPE(IT,J) =AJ(J,J_TYPE,IT)*BYA1
AJ(J,J_TRNFP0,IT)=AJ(J,J_HSURF,IT)+A2BYA1*AJ(J,J_TRHDT,IT)/DTSRC
AJ(J,J_TRNFP1,IT)=AJ(J,J_HATM,IT) +A2BYA1*AJ(J,J_TRHDT,IT)/DTSRC
AJ(J,J_SRABS ,IT)=AJ(J,J_SRNFP0,IT)-AJ(J,J_SRNFG,IT)
AJ(J,J_RNFP0 ,IT)=AJ(J,J_SRNFP0,IT)+AJ(J,J_TRNFP0,IT)
AJ(J,J_RNFP1 ,IT)=AJ(J,J_SRNFP1,IT)+AJ(J,J_TRNFP1,IT)
AJ(J,J_RHDT ,IT)=A1BYA2*AJ(J,J_SRNFG,IT)*DTSRC+AJ(J,J_TRHDT,IT)
AJ(J,J_PRCP ,IT)=AJ(J,J_PRCPSS,IT)+AJ(J,J_PRCPMC,IT)
AJ(J,J_HZ0,IT)=AJ(J,J_RHDT,IT)+AJ(J,J_SHDT,IT)+
* AJ(J,J_EVHDT,IT)+AJ(J,J_EPRCP,IT)
AJ(J,J_HZ1,IT)=AJ(J,J_HZ0,IT)+AJ(J,J_ERVR,IT)
AJ(J,J_HZ2,IT)=AJ(J,J_HZ1,IT)-AJ(J,J_ERUN,IT)-AJ(J,J_IMPLH,IT)
END DO
END DO
DAYS=(Itime-Itime0)/FLOAT(nday)
C****
C**** LOOP OVER SURFACE TYPES: 1 TO NTYPE
C****
IF (KDIAG(1).GT.7) GO TO 510
DO M=0,NTYPE_OUT-1
if(.not.Qbp(M)) cycle
WRITE (6,901) XLABEL
WRITE (6,902) TERRAIN(M),JYEAR0,AMON0,JDATE0,JHOUR0,
* JYEAR,AMON,JDATE,JHOUR,ITIME,DAYS
WRITE (6,903) (NINT(LAT_DG(J,1)),J=JM,INC,-INC)
WRITE (6,905)
NDER=1
KDER=1
DO N=1,k_j_out
IACC=IDACC(IA_J(N))
C**** set weighting for denominator (different only for J_TYPE)
MD=M
IF (name_j(N).eq.'J_surf_type_frac') MD=0
GSUM=0.
GWT=0.
DO JHEMI=1,2
HSUM=0.
HWT=0.
DO JH=1,JMHALF
J=(JHEMI-1)*JMHALF+JH
C**** Sum over types
QJ=0
WTX=0
DO IT=1,NTYPE
QJ =QJ +WT(M ,IT)*AJ(J,N,IT)
WTX=WTX+WT(MD,IT)*SPTYPE(IT,J)
END DO
QJ=QJ*SCALE_J(N)
WTX=WTX*IACC
FLAT(J)=QJ/(WTX+teeny)
MLAT(J)=NINTlimit( FLAT(J) )
HSUM=HSUM+QJ*DXYP(J)*(FIM+1.-S1(J))
HWT=HWT+WTX*DXYP(J)*(FIM+1.-S1(J))
END DO
FHEM(JHEMI)=HSUM/(HWT+teeny)
GSUM=GSUM+HSUM
GWT=GWT+HWT
END DO
FGLOB=GSUM/(GWT+teeny)
IF (M.EQ.0) CALL KEYDJ (N,FGLOB,FHEM(2))
C**** Save BUDG for full output
BUDG(1:JM,N)=FLAT(1:JM)
BUDG(JM+1,N)=FHEM(1)
BUDG(JM+2,N)=FHEM(2)
BUDG(JM+3,N)=FGLOB
TITLEO(N)=STITLE_J(N)
LNAMEO(N)=LNAME_J(N)
SNAMEO(N)=NAME_J(N)
UNITSO(N)=UNITS_J(N)
C**** select output format depending on field name
SELECT CASE (name_j(N)(3:len_trim(name_j(N))))
CASE ('sstab_trop')
WRITE (6,906) STITLE_J(N),FGLOB,FHEM(2),FHEM(1),
* (FLAT(J),J=JM,INC,-INC)
CASE ('evap','prec','ross_num_strat','ross_num_trop'
* ,'ross_radius_strat','ross_radius_trop','ht_runoff'
* ,'river_discharge','ice_melt','impl_m_flux','ht_rvr_disch',
* 'wat_runoff','ssprec','mcprec','h2o_from_ch4'
* ,'lapse_rate','lapse_rate_m','lapse_rate_c'
* ,'ht_thermocline','salt_runoff','s_ice_melt')
WRITE (6,911) STITLE_J(N),FGLOB,FHEM(2),FHEM(1),
* (FLAT(J),J=JM,INC,-INC)
CASE ('ocn_ht_trans','prec_ht_flx','ht_ice_melt')
WRITE (6,912) STITLE_J(N),FGLOB,FHEM(2),FHEM(1),
* (MLAT(J),J=JM,INC,-INC)
CASE DEFAULT
WRITE (6,907) STITLE_J(N),FGLOB,FHEM(2),FHEM(1),
* (MLAT(J),J=JM,INC,-INC)
END SELECT
IF (NDER.le.NDMAX) THEN ! needed to avoid out of bounds address
if (name_j(N)(3:len_trim(name_j(N))).EQ.DERPOS(NDER)) THEN
C**** CALCULATE AND PRINT DERIVED RATIOS
DO KA=KDER,KDER+NDERN(NDER)-1
NN=INUM_J_O(KA)
ND=IDEN_J_O(KA)
C**** differentiate normal ratios from albedo calculations
QALB=(name_j_o(ka).eq.'plan_alb'.or.name_j_o(ka).eq.'surf_alb')
GSUM=0.
GSUM2=0.
DO JHEMI=1,2
HSUM=0.
HSUM2=0.
DO JH=1,JMHALF
J=(JHEMI-1)*JMHALF+JH
C**** Sum over types
QNUM=0
QDEN=0
DO IT=1,NTYPE
QNUM=QNUM+WT(M,IT)*AJ(J,NN,IT)
QDEN=QDEN+WT(M,IT)*AJ(J,ND,IT)
END DO
QNUM=QNUM*SCALE_J_O(KA)
FLAT(J)=QNUM/(QDEN+teeny)
if (QALB) FLAT(J)=100.-FLAT(J)
MLAT(J)=FLAT(J)+.5
HSUM=HSUM+QNUM*DXYP(J)*(FIM+1.-S1(J))
HSUM2=HSUM2+QDEN*DXYP(J)*(FIM+1.-S1(J))
END DO
FHEM(JHEMI)=HSUM/(HSUM2+teeny)
if (QALB) FHEM(JHEMI)=100.-FHEM(JHEMI)
GSUM=GSUM+HSUM
GSUM2=GSUM2+HSUM2
END DO
FGLOB=GSUM/(GSUM2+teeny)
if (QALB) FGLOB=100.-FGLOB
IF (M.EQ.0.AND.name_j_o(ka).eq.'plan_alb') CALL KEYDJA (FGLOB)
C**** Save BUDG for full output
BUDG(1:JM,KA+k_j_out)=FLAT(1:JM)
BUDG(JM+1,KA+k_j_out)=FHEM(1)
BUDG(JM+2,KA+k_j_out)=FHEM(2)
BUDG(JM+3,KA+k_j_out)=FGLOB
TITLEO(KA+k_j_out)=STITLE_J_O(KA)
LNAMEO(KA+k_j_out)=LNAME_J_O(KA)
SNAMEO(KA+k_j_out)=NAME_J_O(KA)
UNITSO(KA+k_j_out)=UNITS_J_O(KA)
C****
SELECT CASE (name_j_o(ka))
CASE ('mc_clddp')
WRITE (6,907) STITLE_J_O(KA),FGLOB,FHEM(2),FHEM(1),
* (MLAT(J),J=JM,INC,-INC)
CASE DEFAULT
WRITE (6,912) STITLE_J_O(KA),FGLOB,FHEM(2),FHEM(1),
* (MLAT(J),J=JM,INC,-INC)
END SELECT
END DO
KDER=KDER+NDERN(NDER)
NDER=NDER+1
END IF
END IF
END DO
WRITE (6,903) (NINT(LAT_DG(J,1)),J=JM,INC,-INC)
WRITE (6,905)
IF (QDIAG) CALL POUT_J(TITLEO,SNAMEO,LNAMEO,UNITSO,BUDG,k_j_out
* +nj_out,TERRAIN(M),M+1) ! the +1 is because M starts at 0
END DO
if(qdiag) call close_j
IF (.not.Qbp(ntype_out)) RETURN
510 CONTINUE
C****
C**** PRODUCE REGIONAL STATISTICS
C****
WRITE (6,901) XLABEL
WRITE (6,902) ' (REGIONS) ',JYEAR0,AMON0,JDATE0,JHOUR0,
* JYEAR,AMON,JDATE,JHOUR,ITIME,DAYS
WRITE(6,918)(NAMREG(1,K),K=1,23),(NAMREG(2,K),K=1,23)
NDER=1
KDER=1
DO N=1,k_j_out
BYIACC=1./(IDACC(IA_J(N))+teeny)
DO JR=1,23
FLAT(JR)=AREG(JR,N)*SCALE_J(N)*BYIACC/SAREA(JR)
MLAT(JR)=NINT(FLAT(JR))
END DO
C**** select output format based on field name
SELECT CASE (name_j(N)(3:len_trim(name_j(N))))
CASE ('evap','prec','ocn_lak_ice_frac','snow_cover'
* ,'ht_ice_melt','impl_m_flux','impl_ht','h2o_from_ch4'
* ,'ice_melt','ht_runoff','wat_g1','river_discharge'
* ,'ht_rvr_disch','ice_g1','snowdp','wat_runoff','ssprec'
* ,'mcprec','atmh2o','ht_thermocline','salt_runoff'
* ,'s_ice_melt')
WRITE (6,910) STITLE_J(N),(FLAT(JR),JR=1,23)
CASE ('sstab_trop','sstab_strat','ross_num_strat'
* ,'ross_num_trop','ross_radius_strat','ross_radius_trop'
* ,'surf_type_frac','lapse_rate','lapse_rate_m'
* ,'lapse_rate_c','rich_num_trop','rich_num_strat'
* ,'dtdlat_strat','dtdlat_trop')
CONTINUE ! no output for not-calculated quantities
CASE DEFAULT
WRITE (6,909) STITLE_J(N),(MLAT(JR),JR=1,23)
END SELECT
IF (NDER.le.NDMAX) THEN ! needed to avoid out of bounds address
IF (name_j(N)(3:len_trim(name_j(N))).EQ.DERPOS(NDER)) THEN
C**** CALCULATE AND PRINT DERIVED RATIOS FOR REGIONAL STATISTICS
DO KA=KDER,KDER+NDERN(NDER)-1
NN=INUM_J_O(KA)
ND=IDEN_J_O(KA)
C**** differentiate normal ratios from albedo calculations
QALB=(name_j_o(ka).eq.'plan_alb'.or.name_j_o(ka).eq.'surf_alb')
DO JR=1,23
FLAT(JR)=SCALE_J_O(KA)*AREG(JR,NN)/(AREG(JR,ND)+teeny)
IF (QALB) FLAT(JR)=100.-FLAT(JR)
MLAT(JR)=FLAT(JR)+.5
END DO
WRITE (6,909) STITLE_J_O(KA),(MLAT(JR),JR=1,23)
END DO
KDER=KDER+NDERN(NDER)
NDER=NDER+1
END IF
END IF
END DO
WRITE (6,905)
WRITE(6,918)(NAMREG(1,K),K=1,23),(NAMREG(2,K),K=1,23)
RETURN
C****
901 FORMAT ('1',A)
902 FORMAT ('0** BUDGETS ',A16,'** From:',I6,A6,I2,', Hr',I3,
* 6X,'To:',I6,A6,I2,', Hr',I3,' Model-Time:',I9,5X,
* 'Dif:',F7.2,' Days')
903 FORMAT ('0',131('-')/20X,'G NH SH ',24I4)
905 FORMAT (1X,131('-'))
906 FORMAT (A16,3F7.2,2X,24F4.1)
907 FORMAT (A16,3F7.2,2X,24I4)
909 FORMAT (A16,1X,23I5)
910 FORMAT (A16,1X,23F5.1)
911 FORMAT (A16,3F7.3,2X,24F4.1)
912 FORMAT (A16,3F7.3,2X,24I4)
918 FORMAT ('0',16X,23(1X,A4)/17X,23(1X,A4)/1X,131('-'))
END SUBROUTINE DIAGJ
MODULE BDjkjl 2
!@sum stores information for outputting lat-sigma/pressure diagnostics
!@auth M. Kelley
IMPLICIT NONE
SAVE
!@param names of derived jk/jl output fields
INTEGER :: jl_rad_cool,jk_dudt_econv,jl_nt_lh_e,jl_vt_lh_e,
* jk_psi_cp,jk_dudt_epdiv,jk_stdev_dp,
* jk_dtempdt_econv,jl_phi_amp_wave1,jl_phi_phase_wave1,
* jl_epflx_div,jk_vt_dse_e,jk_vt_lh_eddy,jk_vt_se_eddy,
* jk_tot_vt_se,jk_psi_tem,jk_epflx_v,
* jk_nt_eqgpv,jk_dyn_conv_eddy_geop,jk_nt_sheat_e,
* jk_dyn_conv_dse,jk_seke,jk_eke,
* jk_nt_dse_se,jk_nt_dse_e,jk_tot_nt_dse,
* jk_nt_lh_e,jk_nt_see,jk_tot_nt_se,
* jk_nt_am_stand_eddy,jk_nt_am_eddy,jk_tot_nt_am,
* jk_we_flx_nor,jk_we_flx_div,jk_refr_ind_wave1,
* jk_del_qgpv,jk_nt_lh_se,jk_wstar,jk_vstar,
* jl_mcdrgpm10,jl_mcdrgpm40,jl_mcdrgpm20,jl_sumdrg
END MODULE BDjkjl
SUBROUTINE JKJL_TITLEX 1,11
!@sum JKJL_TITLEX titles etc for composite jl, jk output
!@auth G. Schmidt/M. Kelley/J. Lerner
!@ver 1.0
use filemanager
USE CONSTANT, only : sday,bygrav,sha,lhe
USE MODEL_COM, only : byim,DTsrc
USE BDjkjl
USE DAGCOM
IMPLICIT NONE
INTEGER :: K,kk,iu_Ijk
LOGICAL qIjk,Ql(KAJLx),Qk(KAJKx)
character*80 line
c
c derived JL-arrays
c
k = kajl
c
k = k + 1
jl_rad_cool = k ; jgrid_jl(k) = 1
sname_jl(k) = 'rad_cool'
lname_jl(k) = 'TOTAL RADIATION COOLING RATE'
units_jl(k) = 'W/(m^2*mb)'
ia_jl(k) = ia_rad
pow_jl(k) = -2
k = k + 1
jl_phi_amp_wave1 = k ; jgrid_jl(k) = 1
sname_jl(k) = 'phi_amp_wave1'
lname_jl(k) ='AMPLITUDE OF GEOPOTENTIAL HEIGHT FOR WAVE NUMBER 1'
units_jl(k) = 'METERS'
k = k + 1
sname_jl(k) = 'phi_amp_wave2' ; jgrid_jl(k) = 1
lname_jl(k) ='AMPLITUDE OF GEOPOTENTIAL HEIGHT FOR WAVE NUMBER 2'
units_jl(k) = 'METERS'
k = k + 1
sname_jl(k) = 'phi_amp_wave3' ; jgrid_jl(k) = 1
lname_jl(k) ='AMPLITUDE OF GEOPOTENTIAL HEIGHT FOR WAVE NUMBER 3'
units_jl(k) = 'METERS'
k = k + 1
sname_jl(k) = 'phi_amp_wave4' ; jgrid_jl(k) = 1
lname_jl(k) ='AMPLITUDE OF GEOPOTENTIAL HEIGHT FOR WAVE NUMBER 4'
units_jl(k) = 'METERS'
k = k + 1
jl_phi_phase_wave1 = k ; jgrid_jl(k) = 1
sname_jl(k) = 'phi_phase_wave1'
lname_jl(k) = 'PHASE OF GEOPOTENTIAL HEIGHT FOR WAVE NUMBER 1'
units_jl(k) = 'DEG WEST LONG'
k = k + 1
sname_jl(k) = 'phi_phase_wave2' ; jgrid_jl(k) = 1
lname_jl(k) = 'PHASE OF GEOPOTENTIAL HEIGHT FOR WAVE NUMBER 2'
units_jl(k) = 'DEG WEST LONG'
k = k + 1
sname_jl(k) = 'phi_phase_wave3' ; jgrid_jl(k) = 1
lname_jl(k) = 'PHASE OF GEOPOTENTIAL HEIGHT FOR WAVE NUMBER 3'
units_jl(k) = 'DEG WEST LONG'
k = k + 1
sname_jl(k) = 'phi_phase_wave4' ; jgrid_jl(k) = 1
lname_jl(k) = 'PHASE OF GEOPOTENTIAL HEIGHT FOR WAVE NUMBER 4'
units_jl(k) = 'DEG WEST LONG'
k = k + 1
jl_epflx_div = k ; jgrid_jl(k) = 1
sname_jl(k) = 'epflx_div'
lname_jl(k) = 'DIVERGENCE OF THE ELIASSEN-PALM FLUX'
units_jl(k) = 'm/s^2'
scale_jl(k) = 1.
pow_jl(k) = -6
k = k + 1
jl_mcdrgpm10 = k ; jgrid_jl(k) = 2
sname_jl(k) = 'dudt_mcdrgpm10'
lname_jl(k) = 'DU/DT BY STRAT. MC DRAG C=+/-10R'
units_jl(k) = 'm/s^2'
pow_jl(k) = -6
k = k + 1
jl_mcdrgpm40 = k ; jgrid_jl(k) = 2
sname_jl(k) = 'dudt_mcdrgpm40' !AJL24+25
lname_jl(k) = 'DU/DT BY STRAT. MC DRAG C=+/-40R'
units_jl(k) = 'm/s^2'
pow_jl(k) = -6
k = k + 1
jl_mcdrgpm20 = k ; jgrid_jl(k) = 2
sname_jl(k) = 'dudt_mcdrgpm20' !AJL26+27
lname_jl(k) = 'DU/DT BY STRAT MC DRAG C=+/-20R'
units_jl(k) = 'm/s^2'
pow_jl(k) = -6
k = k + 1
jl_sumdrg = k ; jgrid_jl(k) = 2
sname_jl(k) = 'dudt_sumdrg' !AJL(18+20-27)
lname_jl(k) = 'ZONAL WIND CHANGE BY MTN+DEFORM+SHR+MC DRAG'
units_jl(k) = 'm/s^2'
pow_jl(k) = -6
k = k + 1
jl_nt_lh_e = k
sname_jl(k) = 'nt_lh_eddy' ; jgrid_jl(k) = 2
lname_jl(k) = 'N. TRANSPORT OF LATENT HEAT BY EDDIES'
units_jl(k) = 'W/mb'
scale_jl(k) = 100.*bygrav*xwon*lhe*fim/DTsrc
pow_jl(k) = 10
ia_jl(k) = ia_src
k = k + 1
jl_vt_lh_e = k
sname_jl(k) = 'vt_lh_eddy1' ; jgrid_jl(k) = 1
lname_jl(k) = 'V. TRANSPORT OF LATENT HEAT BY EDDIES'
units_jl(k) = 'W/m^2'
scale_jl(k) = 100.*bygrav*xwon*lhe*byim/DTsrc
pow_jl(k) = 0
ia_jl(k) = ia_src
c Check the count
if (k .gt. KAJLx) then
write (6,*) 'Increase KAJLx=',KAJLx,' to at least ',k
call stop_model('JL_TITLES: KAJLx too small',255)
end if
inquire(file='Ijk',exist=qIjk)
if(.not.qIjk) then
call openunit('Ijk',iu_Ijk,.false.,.false.)
write (iu_Ijk,'(a)') 'List of JL-fields'
do kk = 1,k
write (iu_Ijk,'(i3,1x,a)') kk,lname_jl(kk)
end do
else if(kdiag(2).gt.0) then
Ql=.false.
call openunit('Ijk',iu_Ijk,.false.,.true.)
read (iu_Ijk,'(a)',end=20) line
10 read (iu_Ijk,'(a)',end=20) line
if(line(1:1).eq.'l') go to 20
read(line,'(i3)') kk
Ql(kk)=.true.
go to 10
20 continue
do kk=1,KAJLx
if(.not.Ql(kk)) sname_jl(kk)='skip'
end do
end if
c
c derived JK-arrays
c
k = KAJK
c
k = k + 1
jk_dudt_econv = k ; jgrid_jk(k) = 2
sname_jk(k) = 'dudt_eddy_conv'
lname_jk(k) = 'DU/DT BY EDDY CONVERGENCE (CP)'
units_jk(k) = '10**-6 m/s^2'
scale_jk(k)= 1.D6
k = k + 1
jk_psi_cp = k ; jgrid_jk(k) = 2
sname_jk(k) = 'psi_cp'
lname_jk(k) = 'STREAM FUNCTION (CP)'
units_jk(k) = 'kg/s' !'10**9 KILOGRAMS/SECOND'
scale_jk(k) = 100.*BYGRAV
ia_jk(k) = ia_dga
pow_jk(k) = 9
k = k + 1
jk_dudt_epdiv = k ; jgrid_jk(k) = 2
sname_jk(k) = 'dudt_epdiv'
lname_jk(k) = 'DU/DT BY ELIASSEN-PALM DIVERGENCE (CP)'
units_jk(k) = 'm/s^2'
pow_jk(k) = -6
k = k + 1
jk_stdev_dp = k ; jgrid_jk(k) = 2
sname_jk(k) = 'stdev_dp'
lname_jk(k) = 'STANDARD DEVIATION OF PRESSURE DIFFERENCES'
units_jk(k) = 'MB'
scale_jk(k) = 1.
k = k + 1
jk_dtempdt_econv = k ; jgrid_jk(k) = 1
sname_jk(k) = 'dtempdt_eddy_conv'
lname_jk(k) = 'DTEMP/DT BY EDDY CONVERGENCE (CP)'
units_jk(k) = 'K/DAY'
scale_jk(k) = SDAY
pow_jk(k) = -1
k = k + 1
jk_vt_dse_e = k ; jgrid_jk(k) = 1
sname_jk(k) = 'vt_dse_e'
lname_jk(k) = 'VERT. TRANS. OF DRY STATIC ENERGY BY EDDIES (CP)'
units_jk(k) = 'W/m^2'
scale_jk(k) = -100.*BYGRAV*BYIM
pow_jk(k) = -1
ia_jk(k) = ia_dga
k = k + 1
jk_vt_lh_eddy = k ; jgrid_jk(k) = 1
sname_jk(k) = 'vt_lh_eddy'
lname_jk(k) = 'VERTICAL TRANSPORT OF LATENT HEAT BY EDDIES (CP)'
units_jk(k) = 'W/m^2'
scale_jk(k) = -100.*BYGRAV*BYIM*LHE
ia_jk(k) = ia_dga
k = k + 1
jk_vt_se_eddy = k ; jgrid_jk(k) = 1
sname_jk(k) = 'vt_se_eddy'
lname_jk(k) ='VERTICAL TRANSPORT OF STATIC ENERGY BY EDDIES (CP)'
units_jk(k) = 'W/m^2'
scale_jk(k) = -100.*BYGRAV*BYIM
ia_jk(k) = ia_dga
k = k + 1
jk_tot_vt_se = k ; jgrid_jk(k) = 1
sname_jk(k) = 'tot_vt_se'
lname_jk(k) =
& 'TOTAL LARGE SCALE VERT. TRANS. OF STATIC ENRG (CP)'
units_jk(k) = 'W/m^2'
scale_jk(k) = -100.*BYGRAV*BYIM
pow_jk(k) = 1
ia_jk(k) = ia_dga
k = k + 1
jk_psi_tem = k ; jgrid_jk(k) = 2
sname_jk(k) = 'psi_tem'
lname_jk(k) = 'TRANSFORMED STREAM FUNCTION (CP)'
units_jk(k) = 'kg/s'
scale_jk(k) = 100.*BYGRAV*XWON
ia_jk(k) = ia_dga
pow_jk(k) = 9
k = k + 1
jk_epflx_v = k ; jgrid_jk(k) = 1
sname_jk(k) = 'epflx_vert_cp'
lname_jk(k) = 'VERTICAL ELIASSEN-PALM FLUX (CP)'
units_jk(k) = 'm^2/s^2'
scale_jk(k) = -100.*BYGRAV*BYIM
ia_jk(k) = ia_dga
pow_jk(k) = -2
k = k + 1
jk_nt_eqgpv = k ; jgrid_jk(k) = 1
sname_jk(k) = 'nt_eddy_qgpv'
lname_jk(k) = 'NORTH. TRANS. OF EDDY Q-G POT. VORTICITY'
units_jk(k) = '10**-6 m/s^2'
scale_jk(k) = 1.D6
k = k + 1
jk_dyn_conv_eddy_geop = k ; jgrid_jk(k) = 1
sname_jk(k) = 'dyn_conv_eddy_geop'
lname_jk(k) = 'DYNAMIC CONVERGENCE OF EDDY GEOPOTENTIAL'
units_jk(k) = 'W/(m^2*mb)'
scale_jk(k) = 1d2*BYGRAV
pow_jk(k) = -4
k = k + 1
jk_nt_sheat_e = k ; jgrid_jk(k) = 2
sname_jk(k) = 'nt_sheat_eddy'
lname_jk(k) = 'NORTH. TRANS. OF SENSIBLE HEAT BY EDDIES'
units_jk(k) = 'W/mb'
scale_jk(k) = .25*SHA*XWON*FIM*1d2*BYGRAV
pow_jk(k) = 11
k = k + 1
jk_dyn_conv_dse = k ; jgrid_jk(k) = 1
sname_jk(k) = 'dyn_conv_dse'
lname_jk(k) = 'DYNAMIC CONVERGENCE OF DRY STATIC ENERGY'
units_jk(k) = 'W/(m^2*mb)'
scale_jk(k) = 1d2*BYGRAV
pow_jk(k) = -2
k = k + 1
jk_seke = k ; jgrid_jk(k) = 2
sname_jk(k) = 'stand_eddy_ke'
lname_jk(k) = 'STANDING EDDY KINETIC ENERGY'
units_jk(k) = 'm^2/s^2'
scale_jk(k) = .5
k = k + 1
jk_eke = k ; jgrid_jk(k) = 2
sname_jk(k) = 'eddy_ke'
lname_jk(k) = 'EDDY KINETIC ENERGY'
units_jk(k) = 'm^2/s^2'
scale_jk(k) = .5
k = k + 1
jk_nt_dse_se = k ; jgrid_jk(k) = 2
sname_jk(k) = 'nt_dse_stand_eddy'
lname_jk(k) = 'NOR. TRANS. OF DRY STAT. ENERGY BY STAND. EDDIES'
units_jk(k) = 'W/mb'
scale_jk(k) = .25*XWON*FIM*1d2*BYGRAV
pow_jk(k) = 11
k = k + 1
jk_nt_dse_e = k ; jgrid_jk(k) = 2
sname_jk(k) = 'nt_dse_eddy'
lname_jk(k) = 'NORTH. TRANS. OF DRY STATIC ENERGY BY EDDIES'
units_jk(k) = 'W/mb'
scale_jk(k) = .25*XWON*FIM*1d2*BYGRAV
pow_jk(k) = 11
k = k + 1
jk_tot_nt_dse = k ; jgrid_jk(k) = 2
sname_jk(k) = 'tot_nt_dse'
lname_jk(k) = 'TOTAL NORTH. TRANSPORT OF DRY STATIC ENERGY'
units_jk(k) = 'W/mb'
scale_jk(k) = .25*XWON*FIM*1d2*BYGRAV
pow_jk(k) = 12
k = k + 1
jk_nt_lh_e = k ; jgrid_jk(k) = 2
sname_jk(k) = 'nt_lh_e'
lname_jk(k) = 'NORTHWARD TRANSPORT OF LATENT HEAT BY EDDIES'
units_jk(k) = 'W/mb'
scale_jk(k) = .25*lhe*XWON*FIM*1d2*BYGRAV
pow_jk(k) = 10
k = k + 1
jk_nt_lh_se = k
sname_jk(k) = 'nt_lh_stand_eddy' ; jgrid_jk(k) = 2
lname_jk(k) = 'N. TRANSPORT OF LATENT HEAT BY STAND. EDDIES'
units_jk(k) = 'W/mb'
scale_jk(k) = .25*lhe*XWON*FIM*1d2*BYGRAV
pow_jk(k) = 9
k = k + 1
jk_nt_see = k ; jgrid_jk(k) = 2
sname_jk(k) = 'nt_se_eddy'
lname_jk(k) = 'NORTH.TRANSPORT OF STATIC ENERGY BY EDDIES'
units_jk(k) = 'W/mb'
scale_jk(k) = .25*XWON*FIM*1d2*BYGRAV
pow_jk(k) = 11
k = k + 1
jk_tot_nt_se = k ; jgrid_jk(k) = 2
sname_jk(k) = 'tot_nt_se'
lname_jk(k) = 'TOTAL NORTHWARD TRANSPORT OF STATIC ENERGY'
units_jk(k) = 'W/mb'
scale_jk(k) = .25*XWON*FIM*1d2*BYGRAV
pow_jk(k) = 12
k = k + 1
jk_nt_am_stand_eddy = k ; jgrid_jk(k) = 2
sname_jk(k) = 'nt_u_stand_eddy'
lname_jk(k) = 'NORTH. TRANS. ZONAL MOM. BY STAND. EDDIES'
units_jk(k) = 'm^2/s^2'
scale_jk(k) = 1.
k = k + 1
jk_nt_am_eddy = k ; jgrid_jk(k) = 2
sname_jk(k) = 'nt_u_eddy'
lname_jk(k) = 'NORTH. TRANS. ZONAL MOM. BY EDDIES'
units_jk(k) = 'm^2/s^2'
scale_jk(k) = 1.
k = k + 1
jk_tot_nt_am = k ; jgrid_jk(k) = 2
sname_jk(k) = 'tot_nt_u'
lname_jk(k) = 'TOTAL NORTH. TRANS. ZONAL MOM.'
units_jk(k) = 'm^2/s^2'
scale_jk(k) = 1.
pow_jk(k) = 1
k = k + 1
jk_we_flx_nor = k ; jgrid_jk(k) = 2
sname_jk(k) = 'we_flx_nor'
lname_jk(k) = 'NORTHWARD WAVE ENERGY FLUX'
c units_jk(k) = '10**11 JOULES/METER/UNIT SIGMA'
units_jk(k) = 'm^2/s^2'
scale_jk(k) = .25
k = k + 1
jk_we_flx_div = k ; jgrid_jk(k) = 1
sname_jk(k) = 'we_flx_div'
lname_jk(k) = 'DIVERGENCE OF THE WAVE ENERGY FLUX'
units_jk(k) = 'm/s^2'
scale_jk(k) = 1
pow_jk(k) = -6
k = k + 1
jk_refr_ind_wave1 = k !!!!! Refraction Inicies must be in order
jgrid_jk(k) = 2
sname_jk(k) = 'refr_ind_wave1'
lname_jk(k) = 'REFRACTION INDEX FOR WAVE NUMBER 1'
units_jk(k) = '10**-8 m^-2'
k = k + 1
sname_jk(k) = 'refr_ind_wave2'
lname_jk(k) = 'REFRACTION INDEX FOR WAVE NUMBER 2'
units_jk(k) = '10**-8 m^-2'
k = k + 1
sname_jk(k) = 'refr_ind_wave3'
lname_jk(k) = 'REFRACTION INDEX FOR WAVE NUMBER 3'
units_jk(k) = '10**-8 m^-2'
k = k + 1
sname_jk(k) = 'refr_ind_wave6'
lname_jk(k) = 'REFRACTION INDEX FOR WAVE NUMBER 6'
units_jk(k) = '10**-8 m^-2'
k = k + 1
sname_jk(k) = 'refr_ind_wave9'
lname_jk(k) = 'REFRACTION INDEX FOR WAVE NUMBER 9'
units_jk(k) = '10**-8 m^-2'
k = k + 1
jk_del_qgpv = k ; jgrid_jk(k) = 2
sname_jk(k) = 'del_qgpv'
lname_jk(k) = 'Q-G POT. VORTICITY CHANGE OVER LATITUDES'
units_jk(k) = '1/(m*s)'
scale_jk(k) = 1.
pow_jk(k) = -12
k = k + 1
jk_wstar = k ; jgrid_jk(k) = 1
sname_jk(k) = 'wstar'
lname_jk(k) = 'W* RESIDUAL VERTICAL VELOCITY'
units_jk(k) = 'mb/s'
pow_jk(k) = -5
scale_jk(k) = 1
k = k + 1
jk_vstar = k ; jgrid_jk(k) = 2
sname_jk(k) = 'vstar'
lname_jk(k) = 'V* = V - D(V''TH''/DTHDP)/DP'
units_jk(k) = 'm/s'
pow_jk(k) = -2
c Check the count
if (k .gt. KAJKx) then
write (6,*) 'Increase KAJKx=',KAJKx,' to at least ',k
call stop_model('JK_TITLES: KAJKx too small',255)
end if
if(.not.qIjk) then
write (iu_Ijk,'(a)') 'list of JK-fields'
do kk = 1,k
write (iu_Ijk,'(i3,1x,a)') kk,lname_jk(kk)
end do
call closeunit(iu_Ijk)
else if(kdiag(2).gt.0) then
Qk=.false.
30 read (iu_Ijk,'(a)',end=40) line
read(line,'(i3)') kk
Qk(kk)=.true.
go to 30
40 continue
do kk=1,KAJKx
if(.not.Qk(kk)) sname_jk(kk)='skip'
end do
call closeunit(iu_Ijk)
end if
RETURN
END SUBROUTINE JKJL_TITLEX
SUBROUTINE DIAGJK 1,124
USE CONSTANT, only :
& grav,rgas,kapa,sday,lhe,twopi,omega,sha,bygrav,tf,teeny
USE MODEL_COM, only :
& im,jm,lm,fim, xlabel,lrunid,DO_GWDRAG,
& BYIM,DSIG,BYDSIG,DT,DTsrc,IDACC,IMH,LS1,NDAA,nidyn,
& PTOP,PMTOP,PSFMPT,SIG,SIGE,JHOUR
USE GEOM, only : JRANGE_HEMI,
& AREAG,BYDXYP,COSP,COSV,DLON,DXV,DXYP,DXYV,DYP,FCOR,RADIUS,WTJ
& ,BYDXYV,lat_dg
USE DAGCOM
USE BDjkjl
IMPLICIT NONE
REAL*8, DIMENSION(JM) ::
& BYDAPO,COSBYPDA,COSBYPV,DXCOSV,ONESPO
& ,DXYPPO,BYDASQR
REAL*8, DIMENSION(JM+LM) :: ONES
REAL*8, DIMENSION(JM,LM) :: AX,BX,CX,DX,VX,EX
REAL*8, DIMENSION(JM,LM) :: BYPDSIG,BYPVDSIG,BYP,BYPV
REAL*8, DIMENSION(JM,LM_REQ) :: ARQX
REAL*8, DIMENSION(LM_REQ) :: BYDPS,BYPKS
REAL*8, DIMENSION(0:IMH) :: AN,BN
REAL*8, DIMENSION(JM,LM,2) :: DSJK
REAL*8, DIMENSION(2,LM,2) :: DSHEM
REAL*8, DIMENSION(LM,2) :: DSGLOB
COMMON/WORK5/DSJK,DSHEM,DSGLOB
REAL*8, DIMENSION(JM,LM,2) :: DPJK
REAL*8, DIMENSION(2,LM,2) :: DPHEM
REAL*8, DIMENSION(LM,2) :: DPGLOB
COMMON/WORKJK/DPJK,DPHEM,DPGLOB
REAL*8, DIMENSION(JM,LM) :: RHO
REAL*8, DIMENSION(LM) :: PM,PKM,PME
REAL*8, DIMENSION(JM,2) :: PJ
REAL*8, DIMENSION(JM,kgz+1,4) :: AMPLTD,PHASE
INTEGER, PARAMETER, DIMENSION(5) :: MW=(/1,2,3,6,9/)
REAL*8, PARAMETER :: ONE=1.
INTEGER ::
& I,IX,J,J1,JH,JHEMI,K,K1,KDN,KM,KUP,L,LR,M,N
REAL*8 ::
& BDN,BUP,BYDP2,BYDPK,BYFSQ,BYIADA,
& BYIMDA,BYN,BYRCOS,DALPHA,DAM4,DE4TI,
& DP,DPG,DPH,DPTI,DTHETA,EL4QI,ELOFIM,
& GBYRSQ,GSQ,PDN,PIG,PIH,PMK,PUP,
& PUTI,PVTI,SCALES,SCALET,SDDP,SKEI,
& SN,SNAMI,SNDEGI,SNELGI,SQM,SQN,SZNDEG,
& SZNELG,THETA,TX,UDXN,UDXS,UX,WTKP1
C**** avoid printing out diagnostics that are not yet defined
if (IDACC(ia_dga).eq.0) RETURN
C**** OPEN PLOTTABLE OUTPUT FILE IF DESIRED
IF(QDIAG) call open_jl(trim(acc_period)//'.jk'//XLABEL(1:LRUNID)
* ,jm,lm,lm_req,lat_dg)
C**** INITIALIZE CERTAIN QUANTITIES
call JKJL_TITLEX
KM=LM
DO 30 L=1,LM
PKM(L)=PLM(L)**KAPA
PME(L)=PSFMPT*SIGE(L)+PTOP
30 PM(L)=PSFMPT*SIGE(L+1)+PTOP
BYDPS(1)=1./(.5*PMTOP)
BYDPS(2)=1./(.3*PMTOP)
BYDPS(3)=1./(.2*PMTOP)
BYPKS(1)=1./(.75*PMTOP)**KAPA
BYPKS(2)=1./(.35*PMTOP)**KAPA
BYPKS(3)=1./(.1*PMTOP)**KAPA
ONES=1.
DO 40 J=1,JM
DXYPPO(J)=DXYP(J)
ONESPO(J)=1.
c BYDXYP(J)=1./DXYP(J)
BYDAPO(J)=BYDXYP(J)
BYDASQR(J)=BYDXYP(J)*BYDXYP(J)
40 CONTINUE
DXYPPO(JM)=DXYP(JM)*FIM
DXYPPO(1)=DXYP(1)*FIM
ONESPO(1)=FIM
ONESPO(JM)=FIM
BYDAPO(1)=BYDAPO(1)*BYIM
BYDAPO(JM)=BYDAPO(JM)*BYIM
DO 50 J=2,JM
DXCOSV(J)=DXV(J)*COSV(J)
50 CONTINUE
DO J=1,JM
BYP(J,1)=IDACC(ia_dga)/(APJ(J,1)+teeny)
BYPV(J,1)=IDACC(ia_dga)/(.25*APJ(J,2)+teeny)
ENDDO
DO L=2,LS1-1
BYP(:,L) = BYP(:,1)
BYPV(:,L) = BYPV(:,1)
ENDDO
DO L=LS1,LM
BYP(:,L) = ONESPO(:)*BYIM/PSFMPT
BYPV(:,L) = ONES(1:JM)*BYIM/PSFMPT
ENDDO
DO L=1,LM
BYPDSIG(:,L) = BYP(:,L)*BYDSIG(L)
BYPVDSIG(:,L) = BYPV(:,L)*BYDSIG(L)
ENDDO
LINECT=65
WRITE (6,901)
BYIADA=1./(IDACC(ia_dga)+teeny)
BYIMDA=BYIADA*BYIM
DO J=1,JM
PJ(J,1)=0
PJ(J,2)=0
DO K=1,KM
PJ(J,1)=PJ(J,1)+AJK(J,K,JK_DPA)
PJ(J,2)=PJ(J,2)+AJK(J,K,JK_DPB)
END DO
IF (J.eq.1) THEN
COSBYPV(J)=0. ! are these values used?
COSBYPDA(J)=0.
ELSE
COSBYPV(J)=COSV(J)/(PJ(J,2)+teeny)
COSBYPDA(J)=COSV(J)*BYDXYP(J)/(PJ(J,1)+teeny)
END IF
END DO
C****
C**** INITIALIZE DELTA SIGMA IN PRESSURE COORDINATES
C****
DO 60 J1=1,2
DO 60 K=1,KM
DPG=0.
PIG=0.
DO 58 JHEMI=1,2
DPH=0.
PIH=0.
DO 55 J=JRANGE_HEMI(1,JHEMI,J1),JRANGE_HEMI(2,JHEMI,J1)
DPJK(J,K,J1) = AJK(J,K,J1)
DSJK(J,K,J1)=AJK(J,K,J1)/(PJ(J,J1)+teeny)
DPH=DPH+AJK(J,K,J1)*WTJ(J,2,J1)
55 PIH=PIH+PJ(J,J1)*WTJ(J,2,J1)
DPHEM(JHEMI,K,J1)=DPH
DSHEM(JHEMI,K,J1)=DPH/(PIH+teeny)
DPG=DPG+DPH
58 PIG=PIG+PIH
DPGLOB(K,J1)=DPG
60 DSGLOB(K,J1)=DPG/(PIG+teeny)
C****
C**** Calculate a density field on tracer grid, edge pressure
C**** (Check this!)
DO K=1,KM-1
DO J=1,JM
IF (K.eq.1) RHO(J,K)=100.*PME(K)/(RGAS*(tf+AJK(J,K,jk_temp)/
* (AJK(J,K,jk_dpa)+teeny)))
IF (K.gt.1) RHO(J,K)=100.*PME(K)/(RGAS*(tf+AJK(J,K-1,jk_temp)/
* (AJK(J,K-1,jk_dpa)+teeny)+
* (AJK(J,K ,jk_temp)/(AJK(J,K-1,jk_dpa)+teeny)
* -AJK(J,K-1,jk_temp)/(AJK(J,K-1,jk_dpa)+teeny))
* *(PME(K)-PLM(K-1))/(PLM(K)-PLM(K-1))))
IF(RHO(J,K).LE.1.D-10) THEN
c print*,"rho<1d-10",j,k,rho(j,k)
RHO(J,K)=100.*PME(K)/(RGAS*(tf+AJK(J+1,K-1,jk_temp)/
* (AJK(J+1,K-1,jk_dpa)+teeny)+
* (AJK(J+1,K ,jk_temp)/(AJK(J+1,K-1,jk_dpa)+teeny)
* -AJK(J+1,K-1,jk_temp)/(AJK(J+1,K-1,jk_dpa)+teeny))
* *(PME(K)-PLM(K-1))/(PLM(K)-PLM(K-1))))
END IF
END DO
END DO
C****
C**** V-V* IS D/DP(V'TH'/DTH/DP) , DX=4*V'TH'/DTH/DP AT INTERFACES
C****
DO 70 J=2,JM
KDN=1
DO 70 K=1,KM
CX(J,K)=FIM*IDACC(ia_dga)
AX(J,K)=AJK(J,K,JK_SHETH)/(AJK(J,K,JK_DPB)+teeny)
KUP=K+1
IF (K.EQ.KM) KUP=KM
VX(J,K)=0.
IF (AJK(J,K,JK_DPB).EQ.0.) GO TO 70
IF (AJK(J,KDN,JK_DPB).EQ.0.) KDN=KDN+1
VX(J,K)=AJK(J,K,JK_DPB)*(AJK(J,KUP,JK_SHETH)/AJK(J,KUP,JK_DPB)-
* AJK(J,KDN,JK_SHETH)/AJK(J,KDN,JK_DPB))/(PM(KUP)-PM(KDN)+.5*
* (AJK(J,KUP,JK_DPB)/AJK(J,KUP,JK_NPTSAVG)-
& AJK(J,KDN,JK_DPB)/AJK(J,KDN,JK_NPTSAVG)))
70 KDN=K
C DX(J,K)=AX*CX=4*(TRANSFORMED STREAM FUNCTION-STREAM FUNCTION)
90 DO 95 J=2,JM
DX(J,KM)=AX(J,KM)*CX(J,KM)
DO 95 K=1,KM-1
WTKP1=AJK(J,K,JK_DPB)/(AJK(J,K+1,JK_DPB)+AJK(J,K,JK_DPB)+teeny)
95 DX(J,K)=(AX(J,K)*(1.-WTKP1)+AX(J,K+1)*WTKP1)*CX(J,K)
C****
C**** PROGNOSTIC QUANTITIES AT CONSTANT PRESSURE
C****
C**** # OF GRIDPOINTS, DELTA P, S.D. OF DELTA P
n = JK_NPTSAVG
SCALET = scale_jk(n)/idacc(ia_jk(n))
CALL JLMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PLM,AJK(1,1,n),SCALET,ONES,ONES,LS1-1,1,JGRID_JK(n))
n = JK_DPB
SCALET = scale_jk(n)/idacc(ia_jk(n))
CALL JLMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PLM,AJK(1,1,n),SCALET,ONES,ONES,LS1-1,2,JGRID_JK(n))
DO 98 J=2,JM
DO 98 K=1,LS1-1
BYN=1./(AJK(J,K,JK_NPTSAVG)+1.D-10)
AX(J,K)=0.
SDDP=(AJK(J,K,JK_DPSQR)-AJK(J,K,JK_DPB)*AJK(J,K,JK_DPB)*BYN)*BYN
98 IF (SDDP.GT.0.) AX(J,K)=SQRT(SDDP)
n = jk_stdev_dp
SCALET = scale_jk(n)
CALL JLMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
& PLM,AX,SCALET,ONES,ONES,LS1-1,2,JGRID_JK(n))
C**** TEMPERATURE, HEIGHT, SPECIFIC AND RELATIVE HUMIDITY
n = JK_TEMP
SCALET = SCALE_JK(n)
SCALES = scale_sjl(1)/idacc(ia_sjl(1))
CALL JKMAPS(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PLM,AJK(1,1,n),SCALET,ONES,ONES,KM,2,JGRID_JK(n),
* ASJL(1,1,1),SCALES,ONESPO,ONES)
n = JK_HGHT
SCALET = SCALE_JK(n)
SCALES = scale_sjl(2)/idacc(ia_sjl(2))
CALL JKMAPS(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PLM,AJK(1,1,n),SCALET,ONES,ONES,KM,2,JGRID_JK(n),
* ASJL(1,1,2),SCALES,ONESPO,ONES)
n = JK_Q
SCALET = SCALE_JK(n)
CALL JKMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PLM,AJK(1,1,n),SCALET,ONES,ONES,KM,2,JGRID_JK(n))
n = JK_RH
SCALET = SCALE_JK(n)
CALL JKMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PLM,AJK(1,1,n),SCALET,ONES,ONES,KM,2,JGRID_JK(n))
C**** LIQUID WATER CONTENT
n = JK_CLDH2O
SCALET = SCALE_JK(n)
CALL JKMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PLM,AJK(1,1,n),SCALET,ONES,ONES,KM,2,JGRID_JK(n))
C**** U AND V WINDS, STREAM FUNCTION
n = JK_U
SCALET = SCALE_JK(n)
CALL JKMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PLM,AJK(1,1,n),SCALET,ONES,ONES,KM,2,JGRID_JK(n))
n = JK_V
SCALET = SCALE_JK(n)
CALL JKMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PLM,AJK(1,1,n),SCALET,ONES,ONES,KM,2,JGRID_JK(n))
DO 100 K=1,KM
DO 100 J=2,JM
100 AX(J,K)=AJK(J,K,JK_V)-.25*VX(J,K)
n = jk_Vstar
CALL JKMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PLM,AX,SCALET,ONES,ONES,KM,2,JGRID_JK(n))
c
c Obtain the stream function by integrating meridional velocity
c downward from the model top
c
do j=2,jm
ax(j,km) = 0.
do k=km-1,1,-1
ax(j,k)=ax(j,k+1)-ajk(j,k+1,jk_v)
enddo
enddo
n = jk_psi_cp
SCALET = SCALE_JK(n)/IDACC(IA_JK(n))
CALL JLMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
& PM,AX,SCALET,DXV,ONES,KM,2,JGRID_JK(n))
DO 110 K=1,KM
DO 110 J=2,JM
110 BX(J,K)=AX(J,K)+.25*DX(J,K)
n = jk_psi_tem
SCALET = SCALE_JK(n)/IDACC(IA_JK(n))
CALL JLMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
& PM,BX,SCALET,DXV,ONES,KM,2,JGRID_JK(n))
C**** RESIDUAL VERTICAL VELOCITY (W*)
DO K=2,KM-1
DO J=2,JM-1
c BX(J,K)=-(BX(J+1,K)-BX(J,K))*
c & (100.*XWON*BYIADA*BYGRAV)*DXV(jmby2)/(RHO(J,K)*XWON*FIM*DXYV(J))
BX(J,K)=-(BX(J+1,K)-BX(J,K))*BYIADA*DXV(jmby2)/(FIM*DXYV(J))
END DO
END DO
BX( 1,:) = 0. ; BX(:,KM) = 0.
BX(JM,:) = 0. ; BX(:, 1) = 0.
n = jk_Wstar
SCALET=SCALE_JK(n)
CALL JLMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PM,BX(1,2),SCALET,ONES,ONES,KM-1,2,JGRID_JK(n))
C**** VERTICAL WINDS
n = JK_VVEL
SCALET = SCALE_JK(n)/IDACC(IA_JK(n))
CALL JLMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PM,AJK(1,2,n),SCALET,BYDXYP,ONES,KM-1,2,JGRID_JK(n))
C****
C**** CALCULATIONS FOR STANDING EDDIES
C****
AX=0.
BX=0.
CX=0.
EX=0.
DO 170 J=2,JM
DO 170 K=1,KM
DPTI=0.
PUTI=0.
PVTI=0.
DE4TI=0.
EL4QI=0.
SKEI=0.
SNDEGI=0.
SNELGI=0.
SNAMI=0.
DO 160 I=1,IM
IF (AIJK(I,J,K,IJK_DP).EQ.0.) GO TO 160
DPTI=DPTI+AIJK(I,J,K,IJK_DP)
BYDPK=1./(AIJK(I,J,K,IJK_DP)+teeny)
PUTI=PUTI+AIJK(I,J,K,IJK_U)
PVTI=PVTI+AIJK(I,J,K,IJK_V)
DE4TI=DE4TI+AIJK(I,J,K,IJK_DSE)
EL4QI=EL4QI+AIJK(I,J,K,IJK_Q)
SKEI=SKEI+(AIJK(I,J,K,IJK_U)*AIJK(I,J,K,IJK_U)
* +AIJK(I,J,K,IJK_V)*AIJK(I,J,K,IJK_V))*BYDPK
SNDEGI=SNDEGI+(AIJK(I,J,K,IJK_DSE)*AIJK(I,J,K,IJK_V)*BYDPK)
SNELGI=SNELGI+(AIJK(I,J,K,IJK_Q)*AIJK(I,J,K,IJK_V)*BYDPK)
SNAMI=SNAMI+AIJK(I,J,K,IJK_U)*AIJK(I,J,K,IJK_V)*BYDPK
160 CONTINUE
AX(J,K)=SKEI-(PUTI*PUTI+PVTI*PVTI)/(DPTI+teeny)
SZNDEG=DE4TI*PVTI/(DPTI+teeny)
SZNELG=EL4QI*PVTI/(DPTI+teeny)
EX(J,K)=SNELGI-SZNELG
BX(J,K)=SNDEGI-SZNDEG
CX(J,K)=SNAMI-PUTI*PVTI/(DPTI+teeny)
170 CONTINUE
C**** STANDING EDDY, EDDY AND TOTAL KINETIC ENERGY
n = jk_seke
SCALET = scale_jk(n)
CALL JKMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
& PLM,AX,SCALET,ONES,ONES,KM,2,JGRID_jk(n))
DO 200 K=1,KM
DO 200 J=2,JM
200 AX(J,K)=AJK(J,K,JK_TOTKE)-AJK(J,K,JK_ZMFKE)
n = jk_eke
SCALET = scale_jk(n)
CALL JKMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
& PLM,AX,SCALET,ONES,ONES,KM,2,JGRID_jk(n))
n = jk_totke
SCALET = SCALE_JK(n)
CALL JKMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PLM,AJK(1,1,n),SCALET,ONES,ONES,KM,2,JGRID_JK(n))
C**** POTENTIAL TEMPERATURE, POTENTIAL VORTICITY
DO 205 LR=1,LM_REQ
DO 205 J=1,JM
205 ARQX(J,LR)=ASJL(J,LR,1)*BYIMDA*ONESPO(J)+TF
N = JK_THETA
SCALET = SCALE_JK(n)
SCALES = P1000K
CALL JKMAPS(LNAME_JK(N),SNAME_JK(N),UNITS_JK(N),POW_JK(n),
& PLM,AJK(1,1,N),SCALET,ONES,ONES,KM,2,JGRID_JK(N),
& ARQX,SCALES,ONES,BYPKS)
N = JK_POTVORT
SCALET = SCALE_JK(n)/IDACC(IA_JK(n))
CALL JLMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PLM,AJK(1,1,N),SCALET,BYDXYP,ONES,KM,2,JGRID_JK(n))
C****
C**** NORTHWARD TRANSPORTS AT CONSTANT PRESSURE
C****
C**** NORTHWARD TRANSPORT OF SENSIBLE HEAT BY EDDIES
DO 210 K=1,KM
DO 210 J=2,JM
210 AX(J,K)=AJK(J,K,JK_TOTNTSH)-AJK(J,K,JK_ZMFNTSH)
N = jk_nt_sheat_e
SCALET = scale_jk(n)
CALL JKMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
& PLM,AX,SCALET,DXV,ONES,KM,2,JGRID_JK(n))
C**** NORTHWARD TRANSPORT OF DRY STATIC ENERGY BY STANDING EDDIES,
C**** EDDIES, AND TOTAL
C**** Individual wave transports commented out. (gas - 05/2001)
N = jk_nt_dse_se
SCALET = scale_jk(n)
CALL JKMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
& PLM,BX,SCALET,DXV,ONES,KM,2,JGRID_jk(n))
DO 230 K=1,KM
DO 230 J=2,JM
AX(J,K)=SHA*(AJK(J,K,JK_TOTNTSH)-AJK(J,K,JK_ZMFNTSH))+
& (AJK(J,K,JK_TOTNTGEO)-AJK(J,K,JK_ZMFNTGEO))
230 BX(J,K)=SHA*AJK(J,K,JK_TOTNTSH)+AJK(J,K,JK_TOTNTGEO)
n = jk_nt_dse_e
SCALET = scale_jk(n)
CALL JKMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
& PLM,AX,SCALET,DXV,ONES,KM,2,JGRID_jk(n))
n = jk_tot_nt_dse
SCALET = scale_jk(n)
CALL JKMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
& PLM,BX,SCALET,DXV,ONES,KM,2,JGRID_jk(n))
C**** NORTHWARD TRANSPORT OF LATENT HEAT BY STAND.EDDY, EDDIES AND TOTAL
C**** New way!
n = jl_nt_lh_e
dx = 0.
DX(2:jm,:)=AJL(2:jm,:,Jl_TOTNTLH)-AJL(2:jm,:,Jl_ZMFNTLH)
SCALET = scale_jl(n)/idacc(ia_jl(n))
dx(1:jm,1:lm) = dx(1:jm,1:lm)*bypvdsig(1:jm,1:lm)
CALL jlMAP(LNAME_jl(n),SNAME_jl(n),UNITS_jl(n),POW_jl(n),
& PLM,DX,SCALET,ONES,ONES,lm,2,JGRID_jl(n))
n = jl_totntlh
SCALET = SCALE_jl(n)/idacc(ia_jl(n))
dx(1:jm,1:lm) = ajl(1:jm,1:lm,n)*bypvdsig(1:jm,1:lm)
CALL jlMAP(LNAME_jl(n),SNAME_jl(n),UNITS_jl(n),POW_jl(n),
& PLM,DX,SCALET,ONES,ONES,lm,2,JGRID_jl(n))
C**** NORTHWARD TRANSPORT OF LATENT HEAT BY STAND. EDDY, EDDIES AND TOTA
C**** Old Way! NOTE: AX is needed later
dx=0.
DO 240 K=1,KM
DO 240 J=2,JM
DX(J,K)=AJK(J,K,JK_TOTNTLH)-AJK(J,K,JK_ZMFNTLH)
240 AX(J,K)=AX(J,K)+LHE*DX(J,K)
n = jk_nt_lh_se
SCALET = scale_jk(n)
CALL JKMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
& PLM,EX,SCALET,DXV,ONES,KM,2,JGRID_jk(n))
n = jk_nt_lh_e
SCALET = scale_jk(n)
CALL JKMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
& PLM,DX,SCALET,DXV,ONES,KM,2,JGRID_jk(n))
n = jk_totntlh
SCALET = SCALE_JK(n)
CALL JKMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PLM,AJK(1,1,n),SCALET,DXV,ONES,KM,2,JGRID_JK(n))
C**** NORTHWARD TRANSPORT OF STATIC ENERGY BY EDDIES AND TOTAL
DO 245 K=1,KM
DO 245 J=2,JM
245 DX(J,K)=BX(J,K)+LHE*AJK(J,K,JK_TOTNTLH)
n = jk_nt_see
SCALET = scale_jk(n)
CALL JKMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
& PLM,AX,SCALET,DXV,ONES,KM,2,JGRID_JK(n))
n = jk_tot_nt_se
SCALET = scale_jk(n)
CALL JKMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
& PLM,DX,SCALET,DXV,ONES,KM,2,JGRID_JK(n))
C**** NORTHWARD TRANSPORT OF KINETIC ENERGY
n = jk_totntke
SCALET = SCALE_JK(n)
CALL JKMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PLM,AJK(1,1,n),SCALET,DXV,ONES,KM,2,JGRID_JK(n))
C**** NOR. TRANS. OF MOM, BY STANDING EDDIES, EDDIES, AND TOTAL ANG. MOM
n = jk_nt_am_stand_eddy
SCALET = scale_jk(n)
CALL JKMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
& PLM,CX,SCALET,ONES,ONES,KM,2,JGRID_JK(n))
DO 260 K=1,KM
DO 260 J=2,JM
CX(J,K)=AJK(J,K,JK_TOTNTMOM)-AJK(J,K,JK_ZMFNTMOM)
260 DX(J,K)=AJK(J,K,JK_TOTNTMOM)+RADIUS*OMEGA*COSV(J)*AJK(J,K,JK_V)
n = jk_nt_am_eddy
SCALET = scale_jk(n)
CALL JKMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
& PLM,CX,SCALET,ONES,ONES,KM,2,JGRID_JK(n))
n = jk_tot_nt_am
SCALET = scale_jk(n)
CALL JKMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
& PLM,DX,SCALET,ONES,ONES,KM,2,JGRID_JK(n))
C****
C**** DYNAMIC CONVERGENCE OF ENERGY
C****
DO 370 K=1,KM
C CX(1,K)=-BX(2,K)*.25*DXV(2)
CX(1,K)=0.
CX(JM,K)=BX(JM,K)*.25*DXV(JM)
C DX(1,K)=-(AJK(2,K,JK_TOTNTGEO)-AJK(2,K,JK_ZMFNTGEO))*.25*DXV(2)
DX(1,K)=0.
DX(JM,K)=(AJK(JM,K,JK_TOTNTGEO)-
& AJK(JM,K,JK_ZMFNTGEO))*.25*DXV(JM)
DO 370 J=2,JM-1
CX(J,K)=(BX(J,K)*DXV(J)-BX(J+1,K)*DXV(J+1))*.25
DX(J,K)=((AJK(J,K,JK_TOTNTGEO)-AJK(J,K,JK_ZMFNTGEO))*DXV(J) -
* (AJK(J+1,K,JK_TOTNTGEO)-AJK(J+1,K,JK_ZMFNTGEO))*DXV(J+1))*.25
370 CONTINUE
DO K=1,KM-1
DO J=1,JM
CX(J,K)=CX(J,K)+AJK(J,K,JK_TOTVTDSE)
CX(J,K+1)=CX(J,K+1)-AJK(J,K,JK_TOTVTDSE)
DX(J,K)=DX(J,K)+AJK(J,K,JK_VTGEOEDDY)
DX(J,K+1)=DX(J,K+1)-AJK(J,K,JK_VTGEOEDDY)
END DO
END DO
n = jk_dyn_conv_dse
SCALET = scale_jk(n)
CALL JKMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
& PLM,CX,SCALET,BYDAPO,ONES,KM,2,JGRID_jk(n))
n = jk_dyn_conv_eddy_geop
SCALET = scale_jk(n)
CALL JKMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
& PLM,DX,SCALET,BYDAPO,ONES,KM,2,JGRID_jk(n))
C**** BAROCLINIC EKE GENERATION, P-K BY EDDY PRESSURE GRADIENT FORCE
n = jk_barekegen
SCALET = SCALE_JK(n)
CALL JKMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PLM,AJK(1,1,n),SCALET,BYDXYP,ONES,KM,2,JGRID_JK(n))
n = jk_p2kedpgf
SCALET = SCALE_JK(n)
CALL JKMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PLM,AJK(1,1,n),SCALET,ONES,ONES,KM,2,JGRID_JK(n))
C****
C**** VERTICAL TRANSPORTS
C****
C**** VERTICAL TRANSPORT OF GEOPOTENTIAL ENERGY BY EDDIES
n = jk_vtgeoeddy
SCALET = SCALE_JK(n)/IDACC(IA_JK(n))
CALL JLMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PM,AJK(1,1,n),SCALET,BYDAPO,ONES,KM-1,2,JGRID_JK(n))
C**** VERTICAL TRANSPORT OF DRY STATIC ENERGY BY EDDIES AND TOTAL
DO 390 K=1,KM-1
DO 390 J=1,JM
AX(J,K)=AJK(J,K,JK_TOTVTDSE)-AJK(J,K,JK_ZMFVTDSE)
390 BX(J,K)=AJK(J,K,JK_TOTVTLH)-AJK(J,K,JK_ZMFVTLH)
n = jk_vt_dse_e
SCALET = SCALE_JK(n)/IDACC(IA_JK(n))
CALL JLMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
& PM,AX,SCALET,BYDAPO,ONES,KM-1,2,JGRID_jk(n))
n = jk_totvtdse
SCALET = SCALE_JK(n)/IDACC(IA_JK(n))
CALL JLMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PM,AJK(1,1,n),SCALET,BYDAPO,ONES,KM-1,2,JGRID_JK(n))
C**** VERTICAL TRANSPORT OF LATENT HEAT BY EDDIES AND TOTAL
C**** New way!
n = jl_vt_lh_e
dx = 0.
DX(:,1:lm-1)=AJL(:,1:lm-1,Jl_totvtlh)-AJL(:,1:lm-1,Jl_zmfvtlh)
SCALET = scale_jl(n)/idacc(ia_jl(n))
CALL jlMAP(LNAME_jl(n),SNAME_jl(n),UNITS_jl(n),POW_jl(n),
& PM,DX,SCALET,BYDAPO,ONES,lm-1,2,JGRID_jl(n))
n = jl_totvtlh
SCALET = SCALE_jl(n)/idacc(ia_jl(n))
CALL jlMAP(LNAME_jl(n),SNAME_jl(n),UNITS_jl(n),POW_jl(n),
& PM,Ajl(1,1,n),SCALET,BYDAPO,ONES,lm-1,2,JGRID_jl(n))
C**** VERTICAL TRANSPORT OF LATENT HEAT BY EDDIES AND TOTAL
C**** Old way!
n = jk_vt_lh_eddy
SCALET = SCALE_JK(n)/IDACC(IA_JK(n))
CALL JLMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
& PM,BX,SCALET,BYDAPO,ONES,KM-1,2,JGRID_jk(n))
n = jk_totvtlh
SCALET = SCALE_JK(n)/IDACC(IA_JK(n))
CALL JLMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PM,AJK(1,1,n),SCALET,BYDAPO,ONES,KM-1,2,JGRID_JK(n))
C**** VERTICAL TRANSPORT OF STATIC ENERGY BY EDDIES AND TOTAL
DO 420 K=1,KM-1
DO 420 J=1,JM
AX(J,K)=AX(J,K)+LHE*BX(J,K)
420 BX(J,K)=AJK(J,K,JK_TOTVTDSE)+LHE*AJK(J,K,JK_TOTVTLH)
n = jk_vt_se_eddy
SCALET = SCALE_JK(n)/IDACC(IA_JK(n))
CALL JLMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
& PM,AX,SCALET,BYDAPO,ONES,KM-1,2,JGRID_jk(n))
n = jk_tot_vt_se
SCALET = SCALE_JK(n)/IDACC(IA_JK(n))
CALL JLMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
& PM,BX,SCALET,BYDAPO,ONES,KM-1,2,JGRID_jk(n))
C**** VERTICAL TRANSPORT OF KINETIC ENERGY
n = jk_totvtke
SCALET = SCALE_JK(n)/IDACC(IA_JK(n))
CALL JLMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PM,AJK(1,1,n),SCALET,BYDXYV,ONES,KM-1,2,JGRID_JK(n))
C**** VERTICAL TRANSPORT OF ANGULAR MOMENTUM BY LARGE SCALE MOTIONS
n = jk_vtameddy
SCALET = SCALE_JK(n)/IDACC(IA_JK(n))
AX = 0.
DO K=1,KM-1
DO J=1,JM
AX(J,K)=AJK(J,K,n)/RHO(J,K)
END DO
END DO
CALL JLMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PM,AX,SCALET,BYDXYV,ONES,KM-1,2,JGRID_JK(n))
n = jk_totvtam
SCALET = SCALE_JK(n)/IDACC(IA_JK(n))
DO K=1,KM-1
DO J=1,JM
AX(J,K)=AJK(J,K,n)/RHO(J,K)
END DO
END DO
CALL JLMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PM,AX,SCALET,BYDXYV,ONES,KM-1,2,JGRID_JK(n))
C**** VERTICAL TRANSPORT OF POTENTIAL VORTICITY TOTAL AND BY EDDIES
n = jk_vtpv
SCALET = SCALE_JK(n)/IDACC(IA_JK(n))
CALL JLMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PM,AJK(1,1,n),SCALET,BYDASQR,ONES,KM-1,2,JGRID_JK(n))
n = jk_vtpveddy
SCALET = SCALE_JK(n)/IDACC(IA_JK(n))
CALL JLMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PM,AJK(1,1,N),SCALET,BYDASQR,ONES,KM-1,2,JGRID_JK(n))
C**** NOR. TRANSPORT OF QUASI-GEOSTROPHIC POT. VORTICITY BY EDDIES
DO 490 K=1,KM
AX(1,K)=0.
AX(JM,K)=0.
DX(1,K)=0.
DX(JM,K)=0.
DO 490 J=2,JM-1
AX(J,K)=((AJK(J,K,JK_TOTNTMOM)-AJK(J,K,JK_ZMFNTMOM))*
& DXCOSV(J)/(AJK(J,K,JK_DPB)+teeny)-
* (AJK(J+1,K,JK_TOTNTMOM)-AJK(J+1,K,JK_ZMFNTMOM))*DXCOSV(J+1)/
* (AJK(J+1,K,JK_DPB)+teeny))/COSP(J)
DX(J,K)=FCOR(J)*(VX(J,K)+VX(J+1,K))/
& (AJK(J,K,JK_DPB)+AJK(J+1,K,JK_DPB)+teeny)
490 CONTINUE
DO 500 K=1,KM
DO 500 J=2,JM-1
500 AX(J,K)=AJK(J,K,JK_DPA)*(AX(J,K)+.25*DX(J,K))
n = jk_nt_eqgpv
SCALET = scale_jk(n)
CALL JKMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
& PLM,AX,SCALET,BYDXYP,ONES,KM,2,JGRID_jk(n))
C****
C**** Wave Energy (ELIASSEN PALM) FLUX: NORTHWARD, VERTICAL, DIVERGENCE
C****
c DO 510 K=1,KM
c AX(1,K)=0.
c DO 510 J=2,JM
c UX=AJK(J,K,JK_U)/(AJK(J,K,JK_DPB)+teeny)
c IF (ABS(UX).GE.teeny) GO TO 510
c SN=+1.
c IF (UX.LT.0.) SN=-1.
c UX=SN*teeny
c 510 AX(J,K)=(AJK(J,K,JK_TOTNTGEO)-AJK(J,K,JK_ZMFNTGEO))/UX!*DXV(J)
c n = jk_we_flx_nor
c SCALET = scale_jk(n)
c CALL JKMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
c & PLM,AX,SCALET,ONES,ONES,KM,2,JGRID_jk(n))
c DO 520 K=1,KM-1
c BX(1,K)=0.
c BX(JM,K)=0.
c DO 520 J=1,JM
c IF (J.NE.1.AND.J.NE.JM) GO TO 516 ! corrected 4-25-2000
c IF (J.EQ.1) UX=.5*(AJK(J+1,K,JK_U)+AJK(J+1,K+1,JK_U))/
c * (AJK(J+1,K,JK_DPB)+AJK(J+1,K+1,JK_DPB)+teeny)
c IF (J.EQ.JM) UX=.5*(AJK(J,K,JK_U)+AJK(J,K+1,JK_U))/
c * (AJK(J,K,JK_DPB)+AJK(J,K+1,JK_DPB)+teeny)
c GO TO 518
c 516 UX=(AJK(J,K,JK_U)+AJK(J+1,K,JK_U)
c & +AJK(J,K+1,JK_U)+AJK(J+1,K+1,JK_U))/
c * (AJK(J,K,JK_DPB)+AJK(J+1,K,JK_DPB)+
c & AJK(J,K+1,JK_DPB)+AJK(J+1,K+1,JK_DPB)+teeny)
c 518 IF (ABS(UX).GE.teeny) GO TO 520
c SN=+1.
c IF (UX.LT.0.) SN=-1.
c UX=SN*teeny
c 520 BX(J,K)=AJK(J,K,JK_VTGEOEDDY)/(UX*RHO(J,K))
c n = jk_epflx_v
c SCALET = SCALE_JK(n)/IDACC(IA_JK(n))
c CALL JLMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
c & PM,BX,SCALET,BYDAPO,ONES,KM-1,2,JGRID_jk(n))
c DO 530 K=1,KM
c CX(1,K)=0.
c CX(JM,K)=0.
c DO 530 J=2,JM-1
c 530 CX(J,K)=.25*(AX(J+1,K)-AX(J,K))
c DO 540 K=1,KM-1
c DO 540 J=1,JM
c CX(J,K)=CX(J,K)-BX(J,K)
c 540 CX(J,K+1)=CX(J,K+1)+BX(J,K)
c n = jk_we_flx_div
c SCALET = scale_jk(n)
c CALL JKMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
c & PLM,CX,SCALET,BYDXYP,ONES,KM,2,JGRID_jk(n))
C****
C**** D/DY OF Q-G POTENTIAL VORTICITY AND REFRACTION INDICES
C****
C**** PRELIMINARIES: VERTICAL DERIVATIVES AND N**2
GSQ=GRAV*GRAV
GBYRSQ=GRAV*GRAV/(RGAS*RGAS)
IF (AJK(2,KM,JK_DPA).LT.teeny) GO TO 670 ! ISTART=4,...
DO 600 J=1,JM
K1=1
580 IF (AJK(J,K1,JK_DPA).GT.teeny) GO TO 590
AX(J,K1)=0.
BX(J,K1)=0.
DX(J,K1)=0.
K1=K1+1
GO TO 580
590 KDN=K1
PDN=PM(KDN)+.5*AJK(J,KDN,JK_DPA)/(AJK(J,KDN,JK_NPTSAVG1)+teeny)
DO 600 K=K1,KM
DP=AJK(J,K,JK_DPA)
PMK=PM(K)+.5*AJK(J,K,JK_DPA)/(AJK(J,K,JK_NPTSAVG1)+teeny)
KUP=K+1
IF (K.EQ.KM) KUP=KM
PUP=PM(KUP)+.5*AJK(J,KUP,JK_DPA)/(AJK(J,KUP,JK_NPTSAVG1)+teeny)
DALPHA=(AJK(J,KUP,JK_TEMP)/(AJK(J,KUP,JK_DPA)+teeny)+TF)/PUP-
* (AJK(J,KDN,JK_TEMP)/(AJK(J,KDN,JK_DPA)+teeny)+TF)/PDN
DTHETA=AJK(J,KUP,JK_THETA)/(AJK(J,KUP,JK_DPA)+teeny)-
* AJK(J,KDN,JK_THETA)/(AJK(J,KDN,JK_DPA)+teeny)
THETA=AJK(J,K,JK_THETA)/(AJK(J,K,JK_DPA)+teeny)
TX=AJK(J,K,JK_TEMP)/(AJK(J,K,JK_DPA)+teeny)+TF
IF (ABS(DTHETA).GE.teeny) GO TO 595
SN=+1.
IF (DTHETA.LT.0.) SN=-1.
DTHETA=SN*teeny
595 DX(J,K)=DP*FCOR(J)*PMK*THETA*(PUP-PDN)/(TX*DTHETA*DXYP(J))
AX(J,K)=DALPHA/(PUP-PDN-teeny)
C**** CALCULATE N**2 AT PRESSURE LATITUDES
BX(J,K)=-DP*GSQ*PMK*DTHETA/(RGAS*TX*THETA*(PUP-PDN-teeny))
KDN=K
600 PDN=PMK
C**** CALCULATE Q12 = (D(UDX) + F*DA)/DA
DO 620 K=1,KM
UDXS=0.
DO 610 J=1,JM-1
UDXN=AJK(J+1,K,JK_U)/(AJK(J+1,K,JK_DPB)+teeny)*DXV(J+1)
CX(J,K)=(UDXS-UDXN+FCOR(J))/DXYP(J)
610 UDXS=UDXN
CX(JM,K)=(UDXS+FCOR(JM))/DXYP(JM)
C**** FIND DQ/DY = (Q12(J)-Q12(J-1)+Q3(J)-Q3(J-1))/DY
DO 620 J=JM,2,-1
DP=AJK(J,K,JK_DPB)
AX(J,K)=DP*(CX(J,K)-CX(J-1,K) + (AX(J,K)-AX(J-1,K))*
* (DX(J,K)+DX(J-1,K))/
& (AJK(J,K,JK_DPA)+AJK(J-1,K,JK_DPA)+teeny))/DYP(3)
620 CONTINUE
n = jk_del_qgpv
SCALET = scale_jk(n)
CALL JKMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
& PLM,AX,SCALET,ONES,ONES,KM,2,JGRID_jk(n))
C**** TERMS FOR THE REFRACTION INDEX EXPRESSION
DO 640 J=2,JM
BYFSQ=2.*DXYV(J)*DXYV(J)/(FCOR(J-1)*FCOR(J-1)+FCOR(J)*FCOR(J))
DO 640 K=1,KM
BYDP2=1./(AJK(J-1,K,JK_DPA)+AJK(J,K,JK_DPA)+teeny)
TX=BYDP2*(AJK(J-1,K,JK_TEMP)+AJK(J,K,JK_TEMP))+TF
DX(J,K)=GBYRSQ/(TX*TX)
SQN=BYDP2*(BX(J-1,K)+BX(J,K))
CX(J,K)=SQN*BYFSQ
UX=AJK(J,K,JK_U)
IF (ABS(UX).GE.teeny) GO TO 635
SN=+1.
IF (UX.LT.0.) SN=-1.
UX=SN*teeny
635 AX(J,K)=AX(J,K)/UX
640 CONTINUE
C**** COMBINE TERMS, PRINT OUT REFRACTION INDICES
SCALET = 1.D8
IX = jk_refr_ind_wave1-1
DO 660 M=1,5
SQM=MW(M)*MW(M)
DO 650 J=2,JM
BYRCOS=1./(RADIUS*RADIUS*COSV(J)*COSV(J))
DO 650 K=1,KM
DP=AJK(J,K,JK_DPB)
650 BX(J,K)=DP*(CX(J,K)*(AX(J,K)-SQM*BYRCOS)-.25*DX(J,K))
660 CALL JKMAP(LNAME_jk(M+IX),SNAME_jk(M+IX),UNITS_JK(M+IX),
& POW_JK(M+IX),
& PLM,BX,SCALET,ONES,ONES,KM,2,JGRID_jk(1+IX))
670 CONTINUE
C**** SKIP REMAINING MAPS IF DATA NOT AVAILABLE
IF (AJK(1,1,JK_EPFLXNCP).NE.0.) GO TO 799
C****
C**** CHANGE OF THE MEAN FIELDS OF WIND AND TEMPERATURE
C****
C**** WIND: RATE OF CHANGE, ADVECTION, EDDY CONVERGENCE
IF (IDACC(ia_dga).LE.1) GO TO 730
SCALET = 1./((IDACC(ia_dga)-1)*(DTsrc*NDAA+DT+DT))
n = JK_TOTDUDT
CALL JLMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PLM,AJK(1,1,n),SCALET,ONES,ONES,KM,2,JGRID_JK(n))
730 CONTINUE
C**** Depending on whether EP fluxes have been specially calculated
C**** output full or approximate version
IF (KEP.gt.0) THEN
CALL EPFLXP
ELSE ! these are not very good
AX=0.
BX=0.
DO 720 K=2,KM-1
DO 720 J=2,JM-1
if (AJK(J,K,JK_DPA).gt.0) AX(J,K)=((AJK(J,K,JK_TOTNTMOM)-
& AJK(J,K,JK_ZMFNTMOM))*DXV(J)-(AJK(J+1,K,JK_TOTNTMOM)-
* AJK(J+1,K,JK_ZMFNTMOM))*DXV(J+1))/
& (AJK(J,K,JK_DPA)*DXYP(J))+
* .125*((AJK(J,K,JK_VTAMEDDY)-AJK(J,K-1,JK_VTAMEDDY))/
& (AJK(J,K,JK_DPB)*DXYV(J)+teeny)+
* (AJK(J+1,K,JK_VTAMEDDY)-AJK(J+1,K-1,JK_VTAMEDDY))/
& (AJK(J+1,K,JK_DPB)*DXYV(J+1)+teeny))
BX(J,K)=(AJK(J+1,K,JK_EPFLXNCP)*DXCOSV(J+1)-
& AJK(J,K,JK_EPFLXNCP)*DXCOSV(J))/
* (DXYP(J)*COSP(J))+.5*(AJK(J,K-1,JK_EPFLXVCP)-
& AJK(J,K,JK_EPFLXVCP)+
* AJK(J+1,K-1,JK_EPFLXVCP)-AJK(J+1,K,JK_EPFLXVCP))/(PM(K-1)-PM(K))
720 CONTINUE
n = jk_dudtmadv
SCALET = SCALE_JK(n)/IDACC(IA_JK(n))
CALL JLMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PLM,AJK(1,1,n),SCALET,ONES,ONES,KM,2,JGRID_JK(n))
n = jk_dudt_econv
SCALET = scale_jk(n)
CALL JLMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
& PLM,AX,SCALET,ONES,ONES,KM,2,jgrid_jk(n))
C**** WIND: TRANSFORMED ADVECTION, LAGRANGIAN CONVERGENCE (DEL.F)
n = jk_dudttem
SCALET = SCALE_JK(n)/IDACC(IA_JK(n))
CALL JLMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PLM,AJK(1,1,n),SCALET,ONES,ONES,KM,2,JGRID_JK(n))
n = jk_dudt_epdiv
CALL JLMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
& PLM,BX,SCALET,ONES,ONES,KM,2,jgrid_jk(n))
END IF
C**** WIND: DU/DT BY STRAT. DRAG - MTN, DEFORM., SHEAR ...
if (DO_GWDRAG) then
SCALET = scale_jl(jl_dudfmdrg)/idacc(ia_jl(jl_dudfmdrg))
n = jl_dumtndrg
CALL JLMAP(LNAME_JL(n),SNAME_JL(n),UNITS_JL(n),POW_JL(n),
* PLM,AJL(1,1,n),SCALET,ONES,ONES,LM,2,JGRID_JL(n))
n = jl_dudfmdrg
CALL JLMAP(LNAME_JL(n),SNAME_JL(n),UNITS_JL(n),POW_JL(n),
* PLM,AJL(1,1,n),SCALET,ONES,ONES,LM,2,JGRID_JL(n))
n = jl_dushrdrg
CALL JLMAP(LNAME_JL(n),SNAME_JL(n),UNITS_JL(n),POW_JL(n),
* PLM,AJL(1,1,n),SCALET,ONES,ONES,LM,2,JGRID_JL(n))
DX=0.
DO 740 L=1,LM
DO 740 J=1,JM
AX(J,L)=AJL(J,L,jl_dumcdrgm10)+AJL(J,L,jl_dumcdrgp10)
BX(J,L)=AJL(J,L,jl_dumcdrgm40)+AJL(J,L,jl_dumcdrgp40)
740 DX(J,L)=AJL(J,L,jl_dumcdrgm20)+AJL(J,L,jl_dumcdrgp20)
n = jl_mcdrgpm10
CALL JLMAP(LNAME_jl(n),SNAME_jl(n),UNITS_JL(n),POW_JL(n),
& PLM,AX,SCALET,ONES,ONES,LM,2,JGRID_jl(n))
n = jl_mcdrgpm40
CALL JLMAP(LNAME_jl(n),SNAME_jl(n),UNITS_JL(n),POW_JL(n),
& PLM,BX,SCALET,ONES,ONES,LM,2,JGRID_jl(n))
n = jl_mcdrgpm20
CALL JLMAP(LNAME_jl(n),SNAME_jl(n),UNITS_JL(n),POW_JL(n),
& PLM,DX,SCALET,ONES,ONES,LM,2,JGRID_jl(n))
C**** DU/DT BY STRAT. DRAG - TOTAL
DO 745 L=1,LM
DO 745 J=1,JM
745 AX(J,L)=AJL(J,L,jl_dumtndrg)+AJL(J,L,jl_dushrdrg)+
* (AX(J,L)+BX(J,L)+DX(J,L)) + AJL(J,L,jl_dudfmdrg)
* + AJL(J,L,jl_dudtsdif)
n = jl_sumdrg
CALL JLMAP(LNAME_jl(n),SNAME_jl(n),UNITS_JL(n),POW_JL(n),
* PLM,AX,SCALET,ONES,ONES,LM,2,JGRID_JL(n))
C**** CHANGE IN EAST WIND BY STRATOSPHERIC DIFFUSION
n = jl_dudtsdif
CALL JLMAP(LNAME_jl(n),SNAME_jl(n),UNITS_JL(n),POW_JL(n),
* PLM,AJL(1,1,n),SCALET,ONES,ONES,LM,2,JGRID_JL(n))
C**** CHANGE IN EAST WIND BY VERTICAL DIFFUSION
n = jl_dudtvdif
CALL JLMAP(LNAME_jl(n),SNAME_jl(n),UNITS_JL(n),POW_JL(n),
* PLM,AJL(1,1,n),SCALET,ONES,ONES,LM,2,JGRID_JL(n))
end if
C**** DU/DT BY SDRAG
n = jl_dudtsdrg
SCALET = scale_jl(n)/idacc(ia_jl(n))
CALL JLMAP(LNAME_JL(n),SNAME_JL(n),UNITS_JL(n),POW_JL(n),
& PLM,AJL(1,1,n),SCALET,ONES,ONES,LM,2,JGRID_JL(n))
C**** TEMPERATURE: RATE OF CHANGE, ADVECTION, EDDY CONVERGENCE
IF (IDACC(ia_dga).GT.1) then
SCALET = SDAY/((IDACC(ia_dga)-1)*(DTsrc*NDAA+DT+DT))
n = JK_TOTDTDT
CALL JLMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PLM,AJK(1,1,n),SCALET,ONES,PKM,KM,2,JGRID_JK(n))
end if
n = JK_DTDTMADV
SCALET = SCALE_JK(n)/IDACC(IA_JK(n))
CALL JLMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PLM,AJK(1,1,n),SCALET,ONES,PKM,KM,2,JGRID_JK(n))
cx = 0.
do k=2,km-1
do j=2,jm-1
if (AJK(J,K,JK_DPA).gt.0.) CX(J,K)=.25*(
& (AJK(J, K,JK_TOTNTSH)-AJK(J, K,JK_ZMFNTSH))*DXV(J)-
& (AJK(J+1,K,JK_TOTNTSH)-AJK(J+1,K,JK_ZMFNTSH))*DXV(J+1))/
& (AJK(J,K,JK_DPA)*DXYP(J))+
* (AJK(J,K,JK_EDDVTPT) -AJK(J,K-1,JK_EDDVTPT))/
& (PM(K-1)-PM(K))*BYIADA*PKM(K)
end do
end do
n = jk_dtempdt_econv
SCALET = scale_jk(n)
CALL JLMAP(LNAME_jk(n),SNAME_jk(n),UNITS_JK(n),POW_JK(n),
& PLM,CX,SCALET,ONES,ONES,KM,2,JGRID_JK(n))
C**** TEMPERATURE: TRANSFORMED ADVECTION
n = JK_DTDTTEM
SCALET = SCALE_JK(n)/IDACC(IA_JK(n))
CALL JLMAP(LNAME_JK(n),SNAME_JK(n),UNITS_JK(n),POW_JK(n),
& PLM,AJK(1,1,n),SCALET,ONES,PKM,KM,2,JGRID_JK(n))
C**** CHANGE IN TEMPERATURE BY STRATOSPHERIC DRAG
n = jl_dtdtsdrg
SCALET = scale_jl(n)/idacc(ia_jl(n))
CALL JLMAP(LNAME_JL(n),SNAME_JL(n),UNITS_JL(n),POW_JL(n),
& PLM,AJL(1,1,n),SCALET,ONES,PKM,KM,2,JGRID_JL(n))
C**** CHANGE IN TEMPERATURE BY DYNAMICS
n = JL_DTDYN
SCALET = scale_jl(n)/idacc(ia_jl(n))
CALL JLMAP(LNAME_JL(n),SNAME_JL(n),UNITS_JL(n),POW_JL(n),
& PLM,AJL(1,1,n),SCALET,ONESPO,PKM,KM,2,JGRID_JL(n))
799 CONTINUE
C****
C**** Transplanted from DIAGJL
C****
LINECT=65
C**** MASS FLUX MOIST CONVECTION
n = JL_MCMFLX
SCALET = scale_jl(n)/idacc(ia_jl(n))
CALL JLMAP(LNAME_JL(n),SNAME_JL(n),UNITS_JL(n),POW_JL(n),
& PLE,AJL(1,1,n),SCALET,ONES,ONES,LM-1,2,JGRID_JL(n))
C****
C**** RADIATION, CONDENSATION AND CONVECTION
C****
C**** SOLAR AND THERMAL RADIATION HEATING
n = JL_SRHR
SCALET = scale_jl(n)/idacc(ia_jl(n))
SCALES = scale_sjl(3)/idacc(ia_sjl(3))
ax(1:jm,1:lm) = ajl(1:jm,1:lm,n)*bypdsig(1:jm,1:lm)
CALL JLMAPS(LNAME_JL(n),SNAME_JL(n),UNITS_JL(n),POW_JL(n),
& PLM,AX,SCALET,ONES,ONES,LM,2,JGRID_JL(n),
* ASJL(1,1,3),SCALES,ONESPO,BYDPS)
n = JL_TRCR
SCALET = scale_jl(n)/idacc(ia_jl(n))
SCALES = scale_sjl(4)/idacc(ia_sjl(4))
ax(1:jm,1:lm) = ajl(1:jm,1:lm,n)*bypdsig(1:jm,1:lm)
CALL JLMAPS(LNAME_JL(n),SNAME_JL(n),UNITS_JL(n),POW_JL(n),
& PLM,AX,SCALET,ONES,ONES,LM,2,JGRID_JL(n),
* ASJL(1,1,4),SCALES,ONESPO,BYDPS)
DO J=1,JM
DO LR=1,LM_REQ
ARQX(J,LR)=ASJL(J,LR,3)+ASJL(J,LR,4)
ENDDO
DO L=1,LM
AX(J,L)=AJL(J,L,JL_SRHR)+AJL(J,L,JL_TRCR)
ENDDO
ENDDO
n = jl_rad_cool
SCALET = -1./idacc(ia_jl(n))
SCALES = -1.*1d-2/idacc(ia_sjl(4))
ax(1:jm,1:lm) = ax(1:jm,1:lm)*bypdsig(1:jm,1:lm)
CALL JLMAPS(LNAME_jl(n),SNAME_jl(n),UNITS_JL(n),POW_JL(n),
& PLM,AX,SCALET,ONES,ONES,LM,2,JGRID_JL(n),
& ARQX,SCALES,ONESPO,BYDPS)
C**** TOTAL, SUPER SATURATION, CONVECTIVE CLOUD COVER, EFFECTIVE RH
n = JL_TOTCLD
SCALET = scale_jl(n)/idacc(ia_jl(n))
CALL JLMAP(LNAME_JL(n),SNAME_JL(n),UNITS_JL(n),POW_JL(n),
& PLM,AJL(1,1,n),SCALET,ONESPO,ONES,LM,2,JGRID_JL(n))
n = JL_SSCLD
SCALET = scale_jl(n)/idacc(ia_jl(n))
CALL JLMAP(LNAME_JL(n),SNAME_JL(n),UNITS_JL(n),POW_JL(n),
& PLM,AJL(1,1,n),SCALET,ONESPO,ONES,LM,2,JGRID_JL(n))
n = JL_MCCLD
SCALET = scale_jl(n)/idacc(ia_jl(n))
CALL JLMAP(LNAME_JL(n),SNAME_JL(n),UNITS_JL(n),POW_JL(n),
& PLM,AJL(1,1,n),SCALET,ONESPO,ONES,LM,2,JGRID_JL(n))
n = JL_RHE
SCALET = scale_jl(n)/idacc(ia_jl(n))
CALL JLMAP(LNAME_JL(n),SNAME_JL(n),UNITS_JL(n),POW_JL(n),
& PLM,AJL(1,1,n),SCALET,ONESPO,ONES,LM,2,JGRID_JL(n))
C**** WATER CLOUD COVER AND ICE CLOUD COVER
n = JL_wcld
SCALET = scale_jl(n)/idacc(ia_jl(n))
CALL JLMAP(LNAME_JL(n),SNAME_JL(n),UNITS_JL(n),POW_JL(n),
& PLM,AJL(1,1,n),SCALET,ONESPO,ONES,LM,2,JGRID_JL(n))
n = JL_icld
SCALET = scale_jl(n)/idacc(ia_jl(n))
CALL JLMAP(LNAME_JL(n),SNAME_JL(n),UNITS_JL(n),POW_JL(n),
& PLM,AJL(1,1,n),SCALET,ONESPO