Assume you have a vector (x,y,z) and you want to rotate it to, say x-axis, you can multiply the rotation matrix to the vector.
First, make the vector a column vector and append 1 to the end. It becomes a 4×1 matrix:
v = x y z 1
Then, convert the vector from Cartesian to spherical coordinate:
[theta, phi, r] = cart2sph(x,y,z)
The meaning of theta, phi and r is illustrated by the following diagram:
Then based on theta and phi, do proper rotations. For example, to rotate the vector to x-axis, you might want to rotate the vector to XZ plane first (by rotate along Z axis with angle theta), and then to x-axis by rotate along Y axis with angle phi. You need to construct the rotation matrix (below) and multiply the matrix with the original vector. Note the order matters.
v_final = M2*M1*v
The following material is from http://www.siggraph.org/education/materials/HyperGraph/modeling/mod_tran/3drota.htm
Rotate along z-axis
( cos q sin q 0 0) (-sin q cos q 0 0) ( 0 0 1 0) ( 0 0 0 1)
Along x-axis
(1 0 0 0) (0 cos q sin q 0) (0 -sin q cos q 0) (0 0 0 1)
Along y-axis
(cos q 0 -sin q 0) (0 1 0 0) (sin q 0 cos q 0) (0 0 0 1)
Note:
To construct a unique rotation matrix, you will need 4 vectors (not on the same plane) in 3D.