  We prove universality of spin correlations in the scaling limit of the planar
Ising model on isoradial graphs with uniformly bounded angles and Z-invariant
weights. Specifically, we show that in the massive scaling limit, i.e., as the
mesh size $\delta$ tends to zero at the same rate as the inverse temperature
goes to the critical one, the two-point spin correlations in the full plane
behave as \[
\delta^{-\frac{1}{4}}\mathbb{E}\left[\sigma_{u_{1}}\sigma_{u_{2}}\right]\ \to\
C_{\sigma}^{2}\cdot\Xi\left(|u_{1}-u_{2}|,m\right)\quad\text{as}\quad\delta\to0,
\] where the universal constant $C_{\sigma}$ and the function
$\Xi(|u_{1}-u_{2}|,m)$ are independent of the lattice. The mass $m$ is defined
by the relation $k'-1\sim 4m\delta$, where $k'$ is the Baxter elliptic
parameter. This includes $m$ of both signs as well as the critical case when
$\Xi(r,0)=r^{-1/4}.$
  These results, together with techniques developed to obtain them, are
sufficient to extend to isoradial graphs the convergence of multi-point spin
correlations in finite planar domains on the square grid, which was established
in a joint work of the first two authors and C. Hongler at criticality, and by
S.C. Park in the sub-critical massive regime. We also give a simple proof of
the fact that the infinite-volume magnetization in the Z-invariant model is
independent of the site and of the lattice.
  As compared to techniques already existing in the literature, we streamline
the analysis of discrete (massive) holomorphic spinors near their ramification
points, relying only upon discrete analogues of the kernel $z^{-1/2}$ for $m=0$
and of $z^{-1/2}e^{\pm 2m|z|}$ for $m\ne 0$. Enabling the generalization to
isoradial graphs and providing a solid ground for further generalizations, our
approach also considerably simplifies the proofs in the square lattice setup.
