  We consider the Schr\"odinger operator in two dimensions with a periodic
potential and a strong constant magnetic field perturbed by slowly varying
non-periodic scalar and vector potentials, $\phi(\epsilon x)$ and $A(\epsilon
x)$, for $\epsilon\ll 1$. For each isolated family of magnetic Bloch bands we
derive an effective Hamiltonian that is unitarily equivalent to the restriction
of the Schr\"odinger operator to a corresponding almost invariant subspace. At
leading order, our effective Hamiltonian can be interpreted as the Peierls
substitution Hamiltonian widely used in physics for non-magnetic Bloch bands.
However, while for non-magnetic Bloch bands the corresponding result is well
understood, for magnetic Bloch bands it is not clear how to even define a
Peierls substitution Hamiltonian beyond a formal expression. The source of the
difficulty is a topological obstruction: magnetic Bloch bundles are generically
not trivializable. As a consequence, Peierls substitution Hamiltonians for
magnetic Bloch bands turn out to be pseudodifferential operators acting on
sections of non-trivial vector bundles over a two-torus, the reduced Brillouin
zone. Part of our contribution is the construction of a suitable Weyl calculus
for such pseudos. As an application of our results we construct a new family of
canonical one-band Hamiltonians $H^B_{\theta,q}$ for magnetic Bloch bands with
Chern number $\theta\in \mathbb{Z}$ that generalizes the Hofstadter model
$H^B_{\rm Hof} = H^B_{0,1}$ for a single non-magnetic Bloch band. It turns out
that $H^B_{\theta,q}$ is isospectral to $H^{q^2B}_{\rm Hof}$ for any $\theta$
and all spectra agree with the Hofstadter spectrum depicted in his famous black
and white butterfly. However, the resulting Chern numbers of subbands,
corresponding to Hall conductivities, depend on $\theta$ and $q$, and thus the
models lead to different colored butterflies.
