!================================================================
! 一次元放射対流平衡大気計算; 放射計算
!
! 97/02/23 小高正嗣
! 97/03/12 小高正嗣; 修正
! 97/03/17 小高正嗣
!================================================================
MODULE rad_module
USE para_module
IMPLICIT NONE
PUBLIC :: get_urf,get_drf,rdeqiv

	CONTAINS
!----------------------------------------------------------------
SUBROUTINE get_bbf(temp,bf)	! 黒体放射束密度の計算

!	INTEGER,INTENT(in)		:: dim	! 配列次元
	REAL,DIMENSION(:),INTENT(in)	:: temp	! 温度
	REAL,DIMENSION(:),INTENT(out)	:: bf	! 黒体放射束密度

	bf = sig*temp**4

END SUBROUTINE get_bbf
!----------------------------------------------------------------
SUBROUTINE get_urf1(nu,tau,temp,uf)	! 上向き放射束密度の計算

	INTEGER,INTENT(in)		:: nu	! 大気層番号
	REAL,DIMENSION(zdim),INTENT(in) :: tau	! 光学的深さ
	REAL,DIMENSION(zdim),INTENT(in)	:: temp	! 温度
	REAL,DIMENSION(zdim),INTENT(out):: uf	! 上向き放射束密度
	REAL,DIMENSION(:),ALLOCATABLE	:: bf	! 黒体放射束密度
	INTEGER i

	ALLOCATE(bf(zdim))
	CALL get_bbf(temp,bf)

	uf = bf
	DO i = 2,nu
	   uf(nu) = uf(nu) &
	          +(bf(i-1)-bf(i+1))*EXP(-3*(tau(i)-tau(nu))/2)/2
	END DO

	DEALLOCATE(bf)

END SUBROUTINE get_urf1
!----------------------------------------------------------------
SUBROUTINE get_urf(tau,temp,uf)	! 全層での上向き放射束密度の計算

	REAL,DIMENSION(zdim),INTENT(in) :: tau	! 光学的深さ
	REAL,DIMENSION(zdim),INTENT(in)	:: temp	! 温度
	REAL,DIMENSION(zdim),INTENT(out):: uf	! 上向き放射束密度
	REAL,DIMENSION(:),ALLOCATABLE	:: bf	! 黒体放射束密度
	INTEGER i,j

	ALLOCATE(bf(zdim))
	CALL get_bbf(temp,bf)

	uf = bf
	DO i = 2,zdim-1
	   DO j = 2,i
	      uf(i) = uf(i) &
		    +(bf(j-1)-bf(j+1))*EXP(-3*(tau(j)-tau(i))/2)/2
	   END DO
	END DO
	uf(zdim) = uf(zdim-1) + (bf(zdim-1)-bf(zdim))

	DEALLOCATE(bf)

END SUBROUTINE get_urf
!----------------------------------------------------------------
SUBROUTINE get_drf1(nu,tau,temp,df)	! 下向き放射束密度の計算

	INTEGER,INTENT(in)		:: nu	! 大気層番号
	REAL,DIMENSION(zdim),INTENT(in) :: tau	! 光学的深さ
	REAL,DIMENSION(zdim),INTENT(in)	:: temp	! 温度
	REAL,DIMENSION(zdim),INTENT(out):: df	! 下向き放射束密度
	REAL,DIMENSION(:),ALLOCATABLE	:: bf	! 黒体放射束密度
	INTEGER i

	ALLOCATE(bf(zdim))
	CALL get_bbf(temp,bf)

	df = 0.
	df(nu) = bf(nu) - bf(zdim)*EXP(-3*(tau(nu))/2)
	DO i = zdim-1,nu,-1
	   df(nu) = df(nu) &
	          -(bf(i-1)-bf(i+1))*EXP(-3*(tau(nu)-tau(i))/2)/2
	END DO

	DEALLOCATE(bf)

END SUBROUTINE get_drf1
!----------------------------------------------------------------
SUBROUTINE get_drf(tau,temp,df)	! 全層での下向き放射束密度の計算

	REAL,DIMENSION(zdim),INTENT(in) :: tau	! 光学的深さ
	REAL,DIMENSION(zdim),INTENT(in)	:: temp	! 温度
	REAL,DIMENSION(zdim),INTENT(out):: df	! 下向き放射束密度
	REAL,DIMENSION(:),ALLOCATABLE	:: bf	! 黒体放射束密度
	INTEGER i,j

	ALLOCATE(bf(zdim))
	CALL get_bbf(temp,bf)

	df = 0.
	DO i = zdim-1,2,-1
	   df(i) = bf(i) - bf(zdim)*EXP(-3*(tau(i))/2)
	   DO j = zdim-1,i,-1
	      df(i) = df(i) &
		    -(bf(j-1)-bf(j+1))*EXP(-3*(tau(i)-tau(j))/2)/2
	   END DO
	END DO
	df(1) = df(2) + (bf(2)-bf(1))

	DEALLOCATE(bf)

END SUBROUTINE get_drf
!----------------------------------------------------------------
SUBROUTINE get_cvf(uf,df,sf,cf)	! 対流熱流速密度の計算

	REAL,DIMENSION(zdim),INTENT(in) :: uf	! 上向き射束密度
	REAL,DIMENSION(zdim),INTENT(in) :: df	! 下向き射束密度
	REAL,DIMENSION(zdim),INTENT(in) :: sf	! 太陽放射束密度
	REAL,DIMENSION(zdim),INTENT(out):: cf	! 対流熱流束密度

	cf = df + sf - uf

END SUBROUTINE get_cvf
!----------------------------------------------------------------
SUBROUTINE rdeqiv(bn,tau,bf,temp,uf0)	! 放射平衡解の計算

	INTEGER,INTENT(in)		  :: bn	! 下端高度番号
	REAL,DIMENSION(zdim),INTENT(in)   :: tau! 光学的厚さ
	REAL,DIMENSION(zdim),INTENT(inout):: bf	! 黒体放射束密度
	REAL,DIMENSION(zdim),INTENT(inout):: temp! 温度
	REAL,INTENT(out) :: uf0			! 大気上端上向き放射束密度
	INTEGER i

!-大気上端上向き放射束密度の計算

	uf0 = 2*sig*temp(bn)**4/(1.5*tau(bn)+1)

!-成層圏の温度; 放射平衡解

	DO i = bn+1,zdim
  	   bf(i)   = 0.5*uf0*(1.5*tau(i) + 1)
	   temp(i) = (bf(i)/sig)**(0.25)
	END DO

END SUBROUTINE rdeqiv
!----------------------------------------------------------------
END MODULE rad_module