srtm_vrtqdr.F90

Go to the documentation of this file.

#ifdef RS6K
@process hot nostrict
#endif
SUBROUTINE srtm_vrtqdr &
 & ( klev , kw,&
 & pref , prefd, ptra , ptrad,&
 & pdbt , prdnd, prup , prupd , ptdbt,&
 & pfd  , pfu  &
 & )  
 
!**** *SRTM_VRTQDR* - VERTICAL QUADRATURE

!     PURPOSE.
!     --------

!          THIS ROUTINE PERFORMS THE VERTICAL INTEGRATION

!**   INTERFACE.
!     ----------

!          *SRTM_VRTQDR* IS CALLED FROM *SRTM_SPCVRT*

!        IMPLICIT ARGUMENTS :
!        --------------------

!     ==== INPUTS ===
!     ==== OUTPUTS ===

!     METHOD.
!     -------

!     EXTERNALS.
!     ----------
!          NONE

!     REFERENCE.
!     ----------

!        SEE RADIATION'S PART OF THE ECMWF RESEARCH DEPARTMENT
!        DOCUMENTATION, AND FOUQUART AND BONNEL (1980)

!     AUTHOR.
!     -------
!        from Howard Barker
!        JEAN-JACQUES MORCRETTE  *ECMWF*

!     MODIFICATIONS.
!     --------------
!        ORIGINAL : 02-10-04
!        M.Hamrud      01-Oct-2003 CY28 Cleaning
!     ------------------------------------------------------------------

USE parkind1  ,ONLY : jpim     ,jprb
USE yomhook   ,ONLY : lhook,   dr_hook

USE parsrtm  , ONLY : jplay, jpgpt 

!USE YOESWN   , ONLY : NDBUG

IMPLICIT NONE

!     ------------------------------------------------------------------

!*       0.1   ARGUMENTS
!              ---------

INTEGER(KIND=JPIM),INTENT(IN)    :: KLEV 
INTEGER(KIND=JPIM),INTENT(IN)    :: KW 
REAL(KIND=JPRB)   ,INTENT(IN)    :: PREF(jplay+1) 
REAL(KIND=JPRB)   ,INTENT(IN)    :: PREFD(jplay+1) 
REAL(KIND=JPRB)   ,INTENT(IN)    :: PTRA(jplay+1) 
REAL(KIND=JPRB)   ,INTENT(IN)    :: PTRAD(jplay+1) 
REAL(KIND=JPRB)   ,INTENT(IN)    :: PDBT(jplay+1) 
REAL(KIND=JPRB)   ,INTENT(OUT)   :: PRDND(jplay+1) 
REAL(KIND=JPRB)   ,INTENT(INOUT) :: PRUP(jplay+1) 
REAL(KIND=JPRB)   ,INTENT(INOUT) :: PRUPD(jplay+1) 
REAL(KIND=JPRB)   ,INTENT(IN)    :: PTDBT(jplay+1) 
REAL(KIND=JPRB)   ,INTENT(INOUT)   :: PFD(jplay+1,jpgpt) 
REAL(KIND=JPRB)   ,INTENT(INOUT)   :: PFU(jplay+1,jpgpt) 
!     ------------------------------------------------------------------

!              ------------

REAL(KIND=JPRB) :: ZTDN(jplay+1)  

INTEGER(KIND=JPIM) :: IKP, IKX, JK, I_NDBUG

REAL(KIND=JPRB) :: ZREFLECT
REAL(KIND=JPRB) :: ZHOOK_HANDLE

!     ------------------------------------------------------------------

! PREF(JK)   direct reflectance
! PREFD(JK)  diffuse reflectance
! PTRA(JK)   direct transmittance
! PTRAD(JK)  diffuse transmittance

! PDBT(JK)   layer mean direct beam transmittance
! PTDBT(JK)  total direct beam transmittance at levels
                   
IF (lhook) CALL dr_hook('SRTM_VRTQDR',0,zhook_handle)
i_ndbug=3
    
!-- link lowest layer with surface
             
zreflect=1.0_jprb / (1.0_jprb -prefd(klev+1)*prefd(klev))
prup(klev)=pref(klev)+(ptrad(klev)* &
 & ((ptra(klev)-pdbt(klev))*prefd(klev+1)+ &
 & pdbt(klev)*pref(klev+1)))*zreflect  
prupd(klev)=prefd(klev)+ptrad(klev)* &
 & ptrad(klev)*prefd(klev+1)*zreflect  

!IF (NDBUG.LE.1) THEN
!  print 9201,PRUP(KLEV),PRUPD(KLEV)
9201 format(1x,'link surf:',6e13.6)
!  print *,'SRTM_VRTQDR after linking with surface layer'
!END IF
    
!-- pass from bottom to top 

DO jk=1,klev-1
  ikp=klev+1-jk                       
  ikx=ikp-1
!  print 9202,JK,IKP,IKX
  9202  format(1x,'Pass from bottom to top:',3i3)      
  zreflect=1.0_jprb / (1.0_jprb -prupd(ikp)*prefd(ikx))
  prup(ikx)=pref(ikx)+(ptrad(ikx)* &
   & ((ptra(ikx)-pdbt(ikx))*prupd(ikp)+ &
   & pdbt(ikx)*prup(ikp)))*zreflect  
  prupd(ikx)=prefd(ikx)+ptrad(ikx)* &
   & ptrad(ikx)*prupd(ikp)*zreflect  

!  print 9203,PRUP(IKX),PRUPD(IKX)
  9203 format(1x,'bot2top:',6e13.6)
ENDDO
!print *,'SRTM_VRTQDR after passing from bottom to top'
    
!-- upper boundary conditions

ztdn(1)=1.0_jprb
prdnd(1)=0.0_jprb
ztdn(2)=ptra(1)
prdnd(2)=prefd(1)

!IF (NDBUG.LE.1) THEN
!  print 9204,ZTDN(1),PRDND(1),ZTDN(2),PRDND(2)
9204 format(1x,'link upper bound:',6e13.6)
!  print *,'SRTM_VRTQDR after upper boundary conditions'
!END IF
    
!-- pass from top to bottom

DO jk=2,klev
  ikp=jk+1
  zreflect=1.0_jprb / (1.0_jprb -prefd(jk)*prdnd(jk))
  ztdn(ikp)=ptdbt(jk)*ptra(jk)+ &
   & (ptrad(jk)*((ztdn(jk)-ptdbt(jk))+ &
   & ptdbt(jk)*pref(jk)*prdnd(jk))) * zreflect  
  prdnd(ikp)=prefd(jk)+ptrad(jk)*ptrad(jk) &
   & *prdnd(jk)*zreflect  

!  IF (NDBUG.LE.1) THEN
!    print 9205,ZTDN(IKP),PRDND(IKP)
  9205 format(1x,'top2bot2:',6e13.6)
!  END IF

ENDDO
!print *,'SRTM_VRTQDR after passing from top to bottom'
                                              
!-- up and down-welling fluxes at levels

DO jk=1,klev+1
!  IF (NDBUG.LE.1) THEN
!    print 9207,JK,PRDND(JK),PRUPD(JK)
!    print 9208,JK,PTDBT(JK),PRUP(JK),ZTDN(JK)
  9207 format(1x,'A',i3,4e13.6)      
  9208 format(1x,'B',i3,4e13.6)      
!  END IF 

  zreflect=1.0_jprb / (1.0_jprb - prdnd(jk)*prupd(jk))
  pfu(jk,kw)=(ptdbt(jk)*prup(jk) + &
   & (ztdn(jk)-ptdbt(jk))*prupd(jk))*zreflect  
  pfd(jk,kw)=ptdbt(jk) + (ztdn(jk)-ptdbt(jk)+ &
   & ptdbt(jk)*prup(jk)*prdnd(jk))*zreflect  

!  IF (NDBUG.LE.2) THEN
!    print 9206,JK,PFU(JK,KW),PFD(JK,KW)
  9206 format(1x,'fluxes:',i3,6e13.6)
!  END IF

ENDDO
!print *,'SRTM_VRTQDR after up and down flux'
   
!print *,'SRTM_VRTQDR about to come out'
!     ------------------------------------------------------------------

IF (lhook) CALL dr_hook('SRTM_VRTQDR',1,zhook_handle)
END SUBROUTINE srtm_vrtqdr