* PACKAGE P2RADI   !" 物理過程(2) 放射(入射)
*
*"履歴 90/10/09 沼口  敦
*"     92/04/17 沼口  敦  
*"     92/08/27 竹広真一
*
***********************************************************************
      SUBROUTINE SHTINS            !" 太陽入射
     O         ( SINS  , SECZET ,
     I           IDATE ,
     C           ALON  , ALAT    )
*
*"   理想的季節・日変化(オプション)
*
*   [PARAM] 
#ifdef SYS_IBMS
      INCLUDE  (ZCDIM)
#else
#include       "zcdim.F"
#endif
*
*   [OUTPUT] 
      REAL       SINS  ( IDIM, JDIM )        !" 入射フラックス
      REAL       SECZET( IDIM, JDIM )        !" sec(入射角)
*
*   [INPUT] 
      INTEGER    IDATE ( 6 )                 !" 日付
      REAL       ALON  ( IDIM )              !" 経度
      REAL       ALAT  ( JDIM )              !" 緯度
*
*   [INTERNAL PARM] 
      REAL       SOLCON
      LOGICAL    OYRAVR, ODYAVR
      REAL       ATALB
      REAL       EPSORB, EQNORB
      REAL       AINS, BINS
      REAL       AZET, BZET
      DATA       SOLCON     / 1380.   /      !" 太陽定数
      DATA       OYRAVR     / .TRUE. /       !" 年平均日射フラグ
      DATA       ODYAVR     / .TRUE. /       !" 日平均日射フラグ
      DATA       ATALB      / 0.2     /      !" 大気アルベド
      DATA       EPSORB     / 23.5    /      !" 赤道傾斜角
      DATA       EQNORB     / 110.    /      !" 昇降点黄経
      DATA       AINS       / 0.127   /      !" 年平均入射の係数
      DATA       BINS       / 0.183   /      !" 年平均入射の係数
      DATA       AZET       / 0.410   /      !" 年平均入射角の係数
      DATA       BZET       / 0.590   /      !" 年平均入射角の係数
*
      NAMELIST  /NMRADI/ SOLCON, OYRAVR, ODYAVR, ATALB ,
     &                   EPSORB, EQNORB , 
     &                   AINS  , BINS  , AZET  , BZET
      SAVE
*
      LOGICAL    OFIRST
      DATA       OFIRST / .TRUE. /
      SAVE       OFIRST
*
*   [INTERNAL WORK] 
      REAL       RDATE ( 2 )                 !" 季節, 時刻
*
      INTEGER    I, J
      INTEGER    IFPAR, JFPAR
      REAL       PI, COSZET, ANGHR, SINDEL
*      
*   [ONCE] 
      IF ( OFIRST ) THEN
         WRITE ( 6,* ) ' SOLAR INSOLATION  P2.2 DATE=92/04/17'
         OFIRST = .FALSE.
*
         CALL   REWNML ( IFPAR , JFPAR )
         READ   ( IFPAR, NMRADI, END=190 )
  190    WRITE  ( JFPAR, NMRADI )
      ENDIF
*
*"          < 0. 季節, 時刻を小数に >
*
      CALL RDATES
     O      ( RDATE ,
     I        IDATE  )
*
      PI = ATAN( 1. ) *4.
*
      IF ( OYRAVR .AND. ODYAVR ) THEN
*
*"           < 1. 年・日平均日射 >
*
         DO 1100 I = 1, IDIM
            DO 1110 J = 1, JDIM
               SINS( I,J ) = - SOLCON * ( 1. - ATALB )
     &                         * ( AINS + BINS*COS( ALAT(J) )**2 )
*
               IF ( SINS( I,J ) .LT. 0. ) THEN
                  SECZET( I,J ) = 1.  /( AZET + BZET*COS( ALAT(J) )**2 )
               ELSE
                  SINS  ( I,J ) = 0.  
                  SECZET( I,J ) = 0.  
               ENDIF
 1110       CONTINUE
 1100    CONTINUE
*
      ELSE IF ( (.NOT.OYRAVR).AND.(.NOT.ODYAVR) ) THEN
*
*"          < 2. 季節・日変化あり >  
*"                   92/08/27 bug-fixed by takepiro <-- ここから

         SINDEL =  SIN( EPSORB/180.* PI  )
     &           * SIN( ( RDATE( 1 ) - EQNORB/360. )*2.*PI ) 
*
         DO 2100 I = 1, IDIM
            DO 2110 J = 1, JDIM
               ANGHR  = RDATE( 2 )*2.*PI - PI + ALON( I ) 
               COSZET = SIN( ALAT(J) )*SINDEL
     &                + COS( ALAT(J) )*SQRT(1.-SINDEL**2)*COS( ANGHR )
*
               IF ( COSZET .GT. 0. ) THEN
                  SINS  ( I,J ) = - SOLCON * ( 1. - ATALB ) * COSZET
                  SECZET( I,J ) = 1. / COSZET
               ELSE
                  SINS  ( I,J ) = 0.  
                  SECZET( I,J ) = 0.  
               ENDIF
 2110       CONTINUE
 2100    CONTINUE
*
*"                                                  <--- ここまで
*
      ELSE
*
         CALL MSGDMP( 'E', 'SHTINS',
     &                'THIS TYPE OF INSOLATION IS NOT IMPREMENTED' ) 
*
      ENDIF
*
      RETURN
      END


