

























































































































 


































































































































































 






















































































































































































 






















































 










 








































































































































































































































 










                                
                                

                                
                                

                                
                                

                                





































module rotation
  use constants
  implicit none
  private
  
  public make_euler
  
contains
  
  subroutine make_euler1 (phi, tt)
    REAL*8, intent(in)  :: phi
    REAL*8, intent(out) :: tt(0:3,0:3)
    
    tt = delta4
    tt(1,1) =  cos(phi)
    tt(1,2) = -sin(phi)
    tt(2,1) =  sin(phi)
    tt(2,2) =  cos(phi)
  end subroutine make_euler1
  
  subroutine make_euler2 (theta, tt)
    REAL*8, intent(in)  :: theta
    REAL*8, intent(out) :: tt(0:3,0:3)
    
    tt = delta4
    tt(2,2) =  cos(theta)
    tt(2,3) = -sin(theta)
    tt(3,2) =  sin(theta)
    tt(3,3) =  cos(theta)
  end subroutine make_euler2
  
  subroutine make_euler3 (psi, tt)
    REAL*8, intent(in)  :: psi
    REAL*8, intent(out) :: tt(0:3,0:3)
    
    tt = delta4
    tt(1,1) =  cos(psi)
    tt(1,2) = -sin(psi)
    tt(2,1) =  sin(psi)
    tt(2,2) =  cos(psi)
  end subroutine make_euler3
  
  subroutine make_euler (phi, theta, psi, tt)
    REAL*8, intent(in)  :: phi, theta, psi
    REAL*8, intent(out) :: tt(0:3,0:3)
    REAL*8 :: tt1(0:3,0:3), tt2(0:3,0:3), tt3(0:3,0:3), tt4(0:3,0:3)
    call make_euler1 (phi, tt1)
    call make_euler2 (theta, tt2)
    call make_euler3 (psi, tt3)
    tt4 = matmul(tt2, tt1)
    tt = matmul(tt3, tt4)
  end subroutine make_euler
  
end module rotation
