This paper considers multivariate Gaussian fields with their associated matrix valued covariance functions. In particular, we characterize the class of stationary-isotropic matrix valued covariance functions on $d$-dimensional Euclidean spaces, as being the scale mixture of the characteristic function of a $d$ dimensional random vector being uniformly distributed on the spherical shell of $\mathbb{R}^{d}$, with a uniquely determined matrix valued and signed measure. This result is the analogue of celebrated Schoenberg theorem, which characterizes stationary and isotropic covariance functions associated to an univariate Gaussian fields. The elements $\mathbf{C}$, being matrix valued, radially symmetric and positive definite on $\mathbb{R}^{d}$, have a matrix valued generator $\mathbf{{\varphi}}$ such that $\mathbf{C}(\boldsymbol{\tau} )=\boldsymbol{\varphi} (\Vert \boldsymbol{\tau} \Vert)$, $\forall\boldsymbol{\tau} \in\mathbb{R}^{d}$, and where $\Vert \cdot\Vert $ is the Euclidean norm. This fact is the crux, together with our analogue of Schoenberg’s theorem, to show the existence of operators that, applied to the generators $\mathbf{{\varphi}}$ of a matrix valued mapping $\mathbf{C}$ being positive definite on $\mathbb{R}^{d}$, allow to obtain generators associated to other matrix valued mappings, say $\tilde{\mathbf{C}}$, being positive definite on Euclidean spaces of different dimensions.