function Interpolate_array00_w(w_data,alon,alat)
!
! 緯度 alat, 経度 alon における関数値を
! その球面調和変換係数 w_data から補間計算する
!
real(8), intent(IN) :: w_data((nm+1)*(nm+1)) ! スペクトルデータ
real(8), intent(IN) :: alon ! 補間する位置(経度)
real(8), intent(IN) :: alat ! 補間する位置(緯度)
real(8) :: Interpolate_array00_w ! 補間した値
real(8) :: mu
real(8) :: y0, y1, y2, AnmPnm
integer :: k,m
mu = sin(alat)
Interpolate_array00_w = 0.0D0
!---- Σa_n^0 P_n^0 の計算
y2 = 0
y1 = 0
do k=nm,1,-1
y0 = alpha(k,0,mu) * y1 + beta(k+1,0)*y2 + w_data(l_nm(k,0))
y2 = y1
y1 = y0
enddo
Interpolate_array00_w = ( beta(1,0) * y2 + mu*sqrt(3.0D0) * y1 + w_data(l_nm(0,0)) ) * Pmm(0,mu)
!---- cos mλ Σa_n^m P_n^m の計算
do m=1,nm
y2 = 0
y1 = 0
do k=nm,m+1,-1
y0 = alpha(k,m,mu) * y1 + beta(k+1,m) * y2 + w_data(l_nm(k,m))
y2 = y1
y1 = y0
enddo
AnmPnm =(w_data(l_nm(m,m)) + beta(m+1,m)*y2 + mu*sqrt(2.0D0*m+3)*y1 ) * Pmm(m,mu)
Interpolate_array00_w = Interpolate_array00_w + AnmPnm*sqrt(2.0D0)*cos(m*alon)
end do
!---- sin λ Σa_n^m P_n^m の計算
do m=1,nm
y2 = 0
y1 = 0
do k=nm,m+1,-1
y0 = alpha(k,m,mu) * y1 + beta(k+1,m)*y2 + w_data(l_nm(k,-m))
y2 = y1
y1 = y0
enddo
AnmPnm =(w_data(l_nm(m,-m)) + beta(m+1,m)*y2 + mu*sqrt(2.0D0*m+3)*y1 ) * Pmm(m,mu)
Interpolate_array00_w = Interpolate_array00_w - AnmPnm*sqrt(2.0D0)*sin(m*alon)
end do
end function Interpolate_array00_w