c
c   Copyright (c) 1990  by  Dwight Meglan
c
      subroutine Mat_Transpose( matrix1, matrix2 )
      implicit none
      
      real*4    matrix1(3,3),matrix2(3,3)

      integer*4 i,j

      do i= 1, 3
        do j= 1, 3
          matrix2(j,i)= matrix1(i,j)
        end do
      end do
      return
      end
c
c   Copyright (c) 1990  by  Dwight Meglan
c
      subroutine Mat_Mult( matrix1,matrix2,matrix3 )
      implicit none
              
      real*4    matrix1(3,3),matrix2(3,3),matrix3(3,3)

      real*4    sum
      integer*4 i,j,k

      do i= 1, 3
        do j= 1, 3
          sum= 0.0
          do k= 1, 3
            sum= sum + matrix1(j,k)*matrix2(k,i)
          end do
          matrix3(j,i)= sum
        end do
      end do
      return
      end
c
c   Copyright (c) 1990  by  Dwight Meglan             
c
      subroutine Vec_Sub( VecA, VecB, VecOut )
      implicit none

      real*4    VecA(3),VecB(3),VecOut(3)

      integer*4 i

      do i= 1,3
        VecOut(i)= VecA(i) - VecB(i)
      end do
      return
      end
c
c   Copyright (c) 1990  by  Dwight Meglan
c
      subroutine Vec_Xprod( VecA, VecB, VecOut )
      implicit none

      real*4    VecA(3),VecB(3),VecOut(3)

      integer*4 x,y,z
      data    x,y,z / 1,2,3 /

      VecOut(x)= VecA(y)*VecB(z) - VecA(z)*VecB(y)
      VecOut(y)= VecA(z)*VecB(x) - VecA(x)*VecB(z)
      VecOut(z)= VecA(x)*VecB(y) - VecA(y)*VecB(x)
      return
      end
c
c   Copyright (c) 1990  by  Dwight Meglan
c
      real*4 function Vec_Dprod( VecA, VecB )
      implicit none
  
      real*4    VecA(3),VecB(3)

      integer*4 i

      Vec_Dprod= 0.0
      do i= 1,3
        Vec_Dprod= Vec_Dprod + VecA(i)*VecB(i)
      end do
      return
      end
c
c   Copyright (c) 1990  by  Dwight Meglan
c
      function vnrm( VecIn, VecOut )
      implicit none

      real*4    VecIn(3), VecOut(3)
      logical*4 vnrm

      real*4    sum
      integer*4 i

      vnrm= .false.
      sum= 0.0
      do i= 1, 3
        sum= sum + VecIn(i)*VecIn(i)
      end do
      sum= sqrt(sum)
      if ( sum .lt. 1.0e-6 ) return
      do i= 1,3
        VecOut(i)= VecIn(i)/sum
      end do
      vnrm= .true.
      return     
      end
