C+
C NAME:
C	quaternion_rotate_xyz
C PURPOSE:
C	Convert vector specified in rectangular coordinates from old
C	to new coordinate system.
C CATEGORY:
C	ucsd/gen/for/lib
C CALLING SEQUENCE:
	subroutine quaternion_rotate_xyz(q,x,y,z)
C INPUTS:
C	q(4)	double precision	quaternion by which rotation is to be performed
C	x,y,z	double precision	vector in rectangular coordinates
C OUTPUTS:
C	x,y,z	double precision	rotated vector in rectangular coordinates
C PROCEDURE:
C	Quaternion is converted to rotation matrix.
C SEE ALSO:
C	quaternion_rotate_angles
C MODIFICATION HISTORY:
C	MAR-2003, Aaron Smith (UCSD/CASS)
C	NOV-2003, Paul Hick (UCSD/CASS; pphick@ucsd.edu)
C	    Changed order of quaternion components from q0,qx,qy,qz to qx,qy,qz,q0
C	    Converted to rotate a vector in rectangular coordinates instead of
C	    vector specified as angles phi and lat (see quaternion_rotate_).
C-
	    double precision	q(4)
	    double precision	x
	    double precision	y
	    double precision	z

	double precision	m(3,3)

	double precision	q11
	double precision	q22
	double precision	q33
	double precision	q44

	double precision	q12
	double precision	q13
	double precision	q23

	double precision	q14
	double precision	q24
	double precision	q34

	q11 = q(1)*q(1)		! Set up rotation matrix r(3,3)
	q22 = q(2)*q(2)
	q33 = q(3)*q(3)
	q44 = q(4)*q(4)

	q12 = q(1)*q(2)
	q13 = q(1)*q(3)
	q23 = q(2)*q(3)

	q14 = q(1)*q(4)
	q24 = q(2)*q(4)
	q34 = q(3)*q(4)

	m(1,1) = q44 + q11 - q22 - q33
	m(1,2) = 2.0d0*( q12 + q34 )
	m(1,3) = 2.0d0*( q13 - q24 )
	m(2,1) = 2.0d0*( q12 - q34 )
	m(2,2) = q44 - q11 + q22 - q33
	m(2,3) = 2.0d0*( q23 + q14 )
	m(3,1) = 2.0d0*( q13 + q24 )
	m(3,2) = 2.0d0*( q23 - q14 )
	m(3,3) = q44 - q11 - q22 + q33

	q11 = m(1,1)*x + m(1,2)*y + m(1,3)*z
	q22 = m(2,1)*x + m(2,2)*y + m(2,3)*z
	q33 = m(3,1)*x + m(3,2)*y + m(3,3)*z

	x = q11
	y = q22
	z = q33

	return
	end