subroutine s2g(ts,t,nc,nc1,nlon,nlath,pln,trigs,& 7,3
factor,ns2,ns3,ncd2)
!$$$ subprogram documentation block
! . . . .
! subprogram: s2g universal routine for spectr. transfer to grid.
! prgmmr: weiyu yang date: 97-11-06
!
! abstract: summation of scalar spherical harmonic series for the
! spectral coefs matrices with vertical distribution format.
!
! program history log:
! 02-02-22 weiyu yang
!
! input argument list on all pes:
! tsx - spectral coefs
! jcap - triangular truncation
! nlon - number of longitudes
! nlath - number of gaussian lats in one hemisphere
! pln - spherical harmonics
! trigs - used by fft
!
! output argument list:
! t - values of desired field on gaussian grid
!
! attributes:
! language: f90
! machine: ibm RS/6000 SP
!
!$$$
use type_kinds
, only: fp_kind
implicit none
#include <machine.h>
integer izero,nlath,nc,nlon,nc1,ncd2,nlla,jj1a,nc11,modn,i2,i3
integer nlath1,i1,i,j,ii,jj,jj1,jj2,jj3,ll,in2,in3,iil,in1,in4
integer iil1,iin1,ii1,ip1,nl1,ns2,ns3,ip5,ip6,ip4,ip2,ip3
real(fp_kind) zero
real(fp_kind),dimension(nc):: ts
real(fp_kind),dimension(2*nlath,nlon+2):: t
#ifdef ibm_sp
real(fp_kind),dimension(nlon,2*nlath)::ttemp
real(fp_kind),dimension(50000+nlon*4):: trigs
#else
real(fp_kind),dimension(nlon+15):: trigs
#endif
real(fp_kind),dimension(ncd2,nlath):: pln
real(fp_kind),dimension(nc1,nlath):: te,to
real(fp_kind),dimension(nc1,nlath):: factor
real(fp_kind),dimension(nlon+2,2*nlath)::tw
!ibm* prefetch_for_store(zero,izero,nlla)
zero=0.0; izero=0
nlla=2*nlath
!ibm* prefetch_for_store(jj1a,nc11,modn,nl1,nlath1)
jj1a=nc1-1; nc11=nc-1
modn=mod(nlath,2); nl1=nlath-modn
nlath1=nlath+1
!ibm* prefetch_for_store(i1)
i1=izero
do i=1,nl1,2
i1=i1+2
do j=1,nc1
!ibm* prefetch_for_store(te(j,i),to(j,i))
te(j,i)=zero
to(j,i)=zero
!ibm* prefetch_for_store(te(j,i1),to(j,i1))
te(j,i1)=zero
to(j,i1)=zero
end do
end do
if(modn==1) then
do j=1,nc1
!ibm* prefetch_for_store(te(j,nlath),to(j,nlath))
te(j,nlath)=zero
to(j,nlath)=zero
end do
end if
if(ns2==1.or.ns2==3.and.ns3==3) then
do j=1,nlath
!ibm* prefetch_for_store(ii,jj,jj1,jj2,jj3)
ii=1; jj=1; jj1=2; jj2=3; jj3=4
do ll=nc1, 4, -4
!ibm* prefetch_for_store(ii1,ip1,ip2,ip3,ip4,ip5,ip6,iil)
ii1=ii+ll-4
ip1=ii-3
ip2=ii-2
ip3=ii-1
ip4=ii1+1
ip5=ii1+2
ip6=ii1+3
iil=ip6
!ibm* prefetch_for_store(in1,in2)
in1=ii/2-1
in2=in1+1
!ibm* prefetch_for_store(in3,in4)
in3=ii1/2+1
in4=in3+1
do i=ii,iil-2,4
ii1=ii1+4
ip1=ip1+4; ip2=ip2+4; ip3=ip3+4
ip4=ip4+4; ip5=ip5+4; ip6=ip6+4
in1=in1+2; in2=in2+2
in3=in3+2; in4=in4+2
te(jj,j)=te(jj,j)+pln(in1,j)*ts(ip1)
te(jj1,j)=te(jj1,j)-pln(in1,j)*ts(i)
to(jj,j)=to(jj,j)+pln(in2,j)*ts(ip3)
to(jj1,j)=to(jj1,j)-pln(in2,j)*ts(ip2)
te(jj2,j)=te(jj2,j)+pln(in3,j)*ts(ip4)
te(jj3,j)=te(jj3,j)-pln(in3,j)*ts(ii1)
to(jj2,j)=to(jj2,j)+pln(in4,j)*ts(ip6)
to(jj3,j)=to(jj3,j)-pln(in4,j)*ts(ip5)
end do
!ibm* prefetch_for_store(iil1,iin1)
iil1=iil-1
iin1=iil/2
te(jj,j)=te(jj,j)+pln(iin1,j)*ts(iil)
te(jj1,j)=te(jj1,j)-pln(iin1,j)*ts(iil1)
ii=ii+ll*2-2; jj=jj+4; jj1=jj1+4; jj2=jj2+4; jj3=jj3+4
end do
te(jj1a,j)=te(jj1a,j)+pln(ncd2,j)*ts(nc)
te(nc1,j)=te(nc1,j)-pln(ncd2,j)*ts(nc11)
end do
else
do j=1,nlath
!ibm* prefetch_for_store(ii,jj,jj1,jj2,jj3)
ii=1; jj=1; jj1=2; jj2=3; jj3=4
do ll=nc1, 4, -4
!ibm* prefetch_for_store(ii1,ip1,ip2,ip3,ip4,ip5,ip6,iil)
ii1=ii+ll-4
ip1=ii-3
ip2=ii-2
ip3=ii-1
ip4=ii1+1
ip5=ii1+2
ip6=ii1+3
iil=ip6
!ibm* prefetch_for_store(in1,in2)
in1=ii/2-1
in2=in1+1
!ibm* prefetch_for_store(in3,in4)
in3=ii1/2+1
in4=in3+1
do i=ii,iil-2,4
ii1=ii1+4
ip1=ip1+4; ip2=ip2+4; ip3=ip3+4
ip4=ip4+4; ip5=ip5+4; ip6=ip6+4
in1=in1+2; in2=in2+2
in3=in3+2; in4=in4+2
te(jj,j)=te(jj,j)+pln(in1,j)*ts(i)
te(jj1,j)=te(jj1,j)+pln(in1,j)*ts(ip1)
to(jj,j)=to(jj,j)+pln(in2,j)*ts(ip2)
to(jj1,j)=to(jj1,j)+pln(in2,j)*ts(ip3)
te(jj2,j)=te(jj2,j)+pln(in3,j)*ts(ii1)
te(jj3,j)=te(jj3,j)+pln(in3,j)*ts(ip4)
to(jj2,j)=to(jj2,j)+pln(in4,j)*ts(ip5)
to(jj3,j)=to(jj3,j)+pln(in4,j)*ts(ip6)
end do
!ibm* prefetch_for_store(iil1,iin1)
iil1=iil-1
iin1=iil/2
te(jj,j)=te(jj,j)+pln(iin1,j)*ts(iil1)
te(jj1,j)=te(jj1,j)+pln(iin1,j)*ts(iil)
ii=ii+ll*2-2; jj=jj+4; jj1=jj1+4; jj2=jj2+4; jj3=jj3+4
end do
te(jj1a,j)=te(jj1a,j)+pln(ncd2,j)*ts(nc11)
te(nc1,j)=te(nc1,j)+pln(ncd2,j)*ts(nc)
end do
end if
!----------
!---------- now combine even and odd parts
!----------
selectcase(ns2)
case(1)
!ibm* prefetch_for_store(i1,i2,i3)
i1=0; i2=nlla+1; i3=i2+1
do i=1,nl1,2
i1=i1+2; i2=i2-2; i3=i3-2
do j=1,nc1
!ibm* prefetch_for_store(tw(j,i),tw(j,i1))
tw(j,i)=te(j,i)+to(j,i)
tw(j,i1)=te(j,i1)+to(j,i1)
!ibm* prefetch_for_store(tw(j,i2),tw(j,i3))
tw(j,i2)=te(j,i1)-to(j,i1)
tw(j,i3)=te(j,i)-to(j,i)
end do
do j=nc1+1,nlon+2
!ibm* prefetch_for_store(tw(j,i),tw(j,i1))
tw(j,i)=zero
tw(j,i1)=zero
!ibm* prefetch_for_store(tw(j,i2),tw(j,i3))
tw(j,i2)=zero
tw(j,i3)=zero
end do
end do
if(modn==1) then
do j=1,nc1
!ibm* prefetch_for_store(tw(j,nlath),tw(j,nlath1))
tw(j,nlath)=te(j,nlath)+to(j,nlath)
tw(j,nlath1)=te(j,nlath)-to(j,nlath)
end do
do j=nc1+1,nlon+2
!ibm* prefetch_for_store(tw(j,nlath),tw(j,nlath1))
tw(j,nlath)=zero
tw(j,nlath1)=zero
end do
end if
case(2)
!ibm* prefetch_for_store(i1,i2,i3)
i1=0; i2=nlla+1; i3=i2+1
do i=1,nl1,2
i1=i1+2; i2=i2-2; i3=i3-2
do j=1,nc1
!ibm* prefetch_for_store(tw(j,i),tw(j,i1))
tw(j,i)=-te(j,i)-to(j,i)
tw(j,i1)=-te(j,i1)-to(j,i1)
!ibm* prefetch_for_store(tw(j,i2),tw(j,i3))
tw(j,i2)=te(j,i1)-to(j,i1)
tw(j,i3)=te(j,i)-to(j,i)
end do
do j=nc1+1,nlon+2
!ibm* prefetch_for_store(tw(j,i),tw(j,i1))
tw(j,i)=zero
tw(j,i1)=zero
!ibm* prefetch_for_store(tw(j,i2),tw(j,i3))
tw(j,i2)=zero
tw(j,i3)=zero
end do
end do
if(modn==1) then
do j=1,nc1
!ibm* prefetch_for_store(tw(j,nlath),tw(j,nlath1))
tw(j,nlath)=-te(j,nlath)-to(j,nlath)
tw(j,nlath1)=te(j,nlath)-to(j,nlath)
end do
do j=nc1+1,nlon+2
!ibm* prefetch_for_store(tw(j,nlath),tw(j,nlath1))
tw(j,nlath)=zero
tw(j,nlath1)=zero
end do
end if
case(3)
!ibm* prefetch_for_store(i1,i2,i3)
i1=0; i2=nlla+1; i3=i2+1
do i=1,nl1,2
i1=i1+2; i2=i2-2; i3=i3-2
do j=1,nc1
!ibm* prefetch_for_store(tw(j,i),tw(j,i1))
tw(j,i)=factor(j,i)*(te(j,i)+to(j,i))
tw(j,i1)=factor(j,i1)*(te(j,i1)+to(j,i1))
!ibm* prefetch_for_store(tw(j,i2),tw(j,i3))
tw(j,i2)=factor(j,i1)*(te(j,i1)-to(j,i1))
tw(j,i3)=factor(j,i)*(te(j,i)-to(j,i))
end do
do j=nc1+1,nlon+2
!ibm* prefetch_for_store(tw(j,i),tw(j,i1))
tw(j,i)=zero
tw(j,i1)=zero
!ibm* prefetch_for_store(tw(j,i2),tw(j,i3))
tw(j,i2)=zero
tw(j,i3)=zero
end do
end do
if(modn==1) then
do j=1,nc1
!ibm* prefetch_for_store(tw(j,nlath),tw(j,nlath1))
tw(j,nlath)=factor(j,nlath)*(te(j,nlath)+to(j,nlath))
tw(j,nlath1)=factor(j,nlath)*(te(j,nlath)-to(j,nlath))
end do
do j=nc1+1,nlon+2
!ibm* prefetch_for_store(tw(j,nlath),tw(j,nlath1))
tw(j,nlath)=zero
tw(j,nlath1)=zero
end do
end if
case(4)
!ibm* prefetch_for_store(i1,i2,i3)
i1=0; i2=nlla+1; i3=i2+1
do i=1,nl1,2
i1=i1+2; i2=i2-2; i3=i3-2
do j=1,nc1
!ibm* prefetch_for_store(tw(j,i),tw(j,i1))
tw(j,i)=-factor(j,i)*(te(j,i)+to(j,i))
tw(j,i1)=-factor(j,i1)*(te(j,i1)+to(j,i1))
!ibm* prefetch_for_store(tw(j,i2),tw(j,i3))
tw(j,i2)=factor(j,i1)*(te(j,i1)-to(j,i1))
tw(j,i3)=factor(j,i)*(te(j,i)-to(j,i))
end do
do j=nc1+1,nlon+2
!ibm* prefetch_for_store(tw(j,i),tw(j,i1))
tw(j,i)=zero
tw(j,i1)=zero
!ibm* prefetch_for_store(tw(j,i2),tw(j,i3))
tw(j,i2)=zero
tw(j,i3)=zero
end do
end do
if(modn==1) then
do j=1,nc1
!ibm* prefetch_for_store(tw(j,nlath),tw(j,nlath1))
tw(j,nlath)=-factor(j,nlath)*(te(j,nlath)+to(j,nlath))
tw(j,nlath1)=factor(j,nlath)*(te(j,nlath)-to(j,nlath))
end do
do j=nc1+1,nlon+2
!ibm* prefetch_for_store(tw(j,nlath),tw(j,nlath1))
tw(j,nlath)=zero
tw(j,nlath1)=zero
end do
end if
end select
#ifdef ibm_sp
call spffte(nlon,nlon/2+1,nlon,nlla,tw,ttemp,1,trigs)
call dgetmo(ttemp,nlon,nlon,nlla,t,nlla)
#else
call zdfftm1du(1,nlon,nlla,tw,1,nlon+2,trigs)
call dgetmo(tw,nlon+2,nlon,nlla,t,nlla)
#endif
return
end