A semi-Lagrangian Characteristic Mapping method for the solution of the tracer transport equations on the sphere is presented. The method solves for the solution operator of the equations by approximating the inverse of the diffeomorphism generated by a given velocity field. The evolution of any tracer and mass density can then be computed via pullback with this map. We present a spatial discretization of the manifold-valued map using a projection-based approach with spherical spline interpolation. The numerical scheme yields C1 continuity for the map and global second-order accuracy for the solution of the tracer transport equations. Error estimates are provided and supported by convergence tests involving solid body rotation, moving vortices, deformational, and compressible flows. Additionally, we illustrate some features of computing the solution operator using a numerical mixing test and the transport of a fractal set in a complex flow environment.