Spherical actuator is an electric device that can achieve multiple degree-of-freedom rotary motions in a single joint. Permanent magnet array is a key factor that influences the output performance of electromagnetic spherical actuators. In this paper, a novel three-dimensional (3D) pole array is proposed to improve the system flux density and thus the output performance. Analysis of magnetic field distribution is extremely important for spherical actuators with 3D magnet arrays. Thus, the investigation of magnetic field is conducted in analytical, numerical, and experimental ways. The general solution of magnetic scalar potential in 3D space is formulated analytically based on Laplace's equations and spherical harmonics, and then specific solutions of the magnetic scalar potential and magnetic flux density are obtained by using boundary conditions. Numerical computation is utilized to validate the analytical model and to facilitate the observation of the magnetic field variation. A research prototype and a testing platform of magnetic field have been developed for experimental study. The testing platform can move the probe to any position around the spherical actuator and measure the flux density automatically. Experiments are conducted to obtain the flux distribution. Both numerical and experimental results validate the analytical model well.