circumsolar_rad.F90

Go to the documentation of this file.

!SFX_LIC Copyright 1994-2014 CNRS, Meteo-France and Universite Paul Sabatier
!SFX_LIC This is part of the SURFEX software governed by the CeCILL-C licence
!SFX_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt  
!SFX_LIC for details. version 1.
!     ###############################################################################
SUBROUTINE circumsolar_rad(PDIR_SW, PSCA_SW, PZENITH, PF1_o_B)
!     ###############################################################################
!
!!****  *CIRCUMSOLAR_RAD * 
!!
!!    PURPOSE
!!    -------
!!
! compute the fraction of the scattered solar radiation that has a
! directionnal behaviour, it is from the circumsolar region
! for an horizontal surface
!!
!!**  METHOD
!!    ------
!!
!!    REFERENCE
!!    ---------
!!      
! Engeneering Reference, EnergyPlus DOE, 2012
! Perez, 1990 
!!
!!    AUTHOR
!!    ------
!!     G. Pigeon 
!!
!!    MODIFICATIONS
!!    -------------
!!      Original    27/11/2012
!!---------------------------------------------------------------
!
!
USE yomhook   ,ONLY : lhook,   dr_hook
USE parkind1  ,ONLY : jprb      
USE modd_csts ,ONLY : xpi

IMPLICIT NONE
!
!*      0.1    declarations of arguments
!
REAL, DIMENSION(:), INTENT(IN) :: pdir_sw ! direct radiation
REAL, DIMENSION(:), INTENT(IN) :: psca_sw ! scattered radiation
REAL, DIMENSION(:), INTENT(IN) :: pzenith ! solar zenith angle
REAL, DIMENSION(:), INTENT(OUT) :: pf1_o_b ! coef
!
!*      0.2    declarations of local variables
!
REAL, DIMENSION(SIZE(PZENITH)) :: ze !sky clearness factor
REAL, DIMENSION(SIZE(PZENITH)) :: zb ! b coef
REAL, DIMENSION(SIZE(PZENITH)) :: zf1 !factors depending on ZE and Z
REAL, DIMENSION(SIZE(PZENITH)) :: zf11, zf13 !factors depending on ZE
REAL, PARAMETER :: zk = 1.041 !coef. from EP Engineering Reference 
REAL, PARAMETER :: zbmin = 0.087 !coef. from EP Engineering Reference 
REAL(KIND=JPRB) :: zhook_handle
REAL, DIMENSION(SIZE(PZENITH)) :: zelev
INTEGER :: jj
!-------------------------------------------------------------------------------
IF (lhook) CALL dr_hook('MODE_CONV_DOE:CHTC_VERT_DOE',0,zhook_handle)

zelev = xpi/2. - pzenith

pf1_o_b = 0.
WHERE(psca_sw > 0)
   ze = ((psca_sw + pdir_sw) / psca_sw + zk * zelev**3) / (1 + zk * zelev**3)
ELSEWHERE
   ze = -1
ENDWHERE

DO jj=1,SIZE(zelev)
   zb(jj) = max(zbmin, cos(zelev(jj)))
ENDDO

 CALL fij(ze, zf11, zf13)
zf1 = zf11 + zf13 * zelev
   
pf1_o_b = zf1 / zb
DO jj=1,SIZE(pf1_o_b)
   pf1_o_b(jj) = min(pf1_o_b(jj), 1.)
   pf1_o_b(jj) = max(pf1_o_b(jj), 0.)
ENDDO

IF (lhook) CALL dr_hook('MODE_CONV_DOE:CHTC_VERT_DOE',1,zhook_handle)
!-------------------------------------------------------------------------------
 CONTAINS
SUBROUTINE fij(PE, PF11, PF13)
REAL, DIMENSION(:), INTENT(IN) :: pe !sky clearness factor
REAL, DIMENSION(:), INTENT(OUT):: pf11
REAL, DIMENSION(:), INTENT(OUT):: pf13

WHERE(pe < 0 )
   pf11 = 0
   pf13 = 0
ELSEWHERE(pe < 1.065)
   pf11 = -0.0083117
   pf13 = -0.0620636
ELSEWHERE(pe >= 1.065 .AND. pe < 1.230 )
   pf11 = 0.1299457
   pf13 = -0.1513752
ELSEWHERE(pe >= 1.230 .AND. pe < 1.500 )
   pf11 = 0.3296958
   pf13 = -0.2210958
ELSEWHERE(pe >= 1.500 .AND. pe < 1.950)
   pf11 = 0.5682053
   pf13 =  -0.2951290
ELSEWHERE(pe >= 1.950 .AND. pe < 2.800 )
   pf11 = 0.8730280
   pf13 = -0.3616149
ELSEWHERE(pe >= 2.800 .AND. pe < 4.500)
   pf11 = 1.1326077
   pf13 = -0.4118494
ELSEWHERE(pe >= 4.500 .AND. pe < 6.200)
   pf11 = 1.0601591
   pf13 = -0.3589221
ELSEWHERE
   pf11 = 0.6777470
   pf13 = -0.2504286
ENDWHERE

END SUBROUTINE fij

END SUBROUTINE circumsolar_rad