!--
!------------------------------------------------------------------------
! Copyright (c) 2011-2013 SPMODEL Development Group. All rights reserved.
!------------------------------------------------------------------------
!== c2ajcc
!
! Copyright :: 2011-2013 Shin-ichi Takehiro, Youhei SASAKI,
!              SPMODEL Development Group
! License :: MIT/X11. Please see COPYRIGHT[link:../COPYRIGHT]
!
subroutine c2ajcc(lm,km,jm,im,sa,sb,sc,ws,wg,itj,tj,iti,ti)
  use dc_types, only: DP
  implicit none
  integer, intent(in) :: lm
  integer, intent(in) :: km
  integer, intent(in) :: jm
  integer, intent(in) :: im
  real(DP), intent(in) :: sa(-km:km,lm)
  real(DP), intent(in) :: sb(-km:km,0:lm)
  real(DP), intent(inout) :: sc(-km:km,0:lm)
  real(DP) :: ws(-km:km,0:lm)
  real(DP) :: wg((jm+1)*im,3)
  integer :: itj(5)
  real(DP) :: tj(jm*6)
  integer :: iti(5)
  real(DP) :: ti(im*2)
  integer :: ji
  integer :: k
  integer :: l

!! a --> wg(ji,3)
  call c2s2ga(lm,km,jm,im,sa,wg(1,3),wg,itj,tj,iti,ti,3)

!! ∂b/∂y --> wg(ji,2)
  call bsset0(2*km+1,ws)
  do l=1,lm
    do k=-km,km
      ws(k,l)=-l*sb(k,l)
    end do
  end do
  call c2s2ga(lm,km,jm,im,ws(-km,1),wg(1,2),wg,itj,tj,iti,ti,3)

!! a × ∂b/∂y  --> wg(ji,2)
  do ji=1,(jm+1)*im
    wg(ji,2)=wg(ji,3)*wg(ji,2)
  end do
!! a × ∂b/∂y のスペクトル --> ws
  call c2g2sa(lm,km,jm,im,wg(1,2),ws,wg,itj,tj,iti,ti,4)

!! ∂(a × ∂b/∂y)/∂x --> sc
  do l=0,lm
    do k=-km,km
      sc(k,l)=-k*ws(-k,l)
    end do
  end do
!! ∂b/∂x --> wg(ji,2)
  do l=0,lm
    do k=-km,km
      ws(k,l)=-k*sb(-k,l)
    end do
  end do
  call c2s2ga(lm,km,jm,im,ws,wg(1,2),wg,itj,tj,iti,ti,4)

!! a × ∂b/∂x  --> wg(ji,2)
  do ji=1,(jm+1)*im
    wg(ji,2)=wg(ji,3)*wg(ji,2)
  end do

!! a × ∂b/∂x のスペクトル --> ws
  call c2g2sa(lm,km,jm,im,wg(1,2),ws(-km,1),wg,itj,tj,iti,ti,3)

!! finally calculate jacobian
  do l=1,lm
    do k=-km,km
      sc(k,l)=sc(k,l)-l*ws(k,l)
    end do
  end do

end subroutine c2ajcc
