As a complement to the other answer here is a solution that uses Gauss-Bonnet:
Let $S_{r}=\{(x,y,z) \mid x^2+y^2=z\leq r^2\}$. By symmetry $\int_{M_{r}}K\;dA=\frac 14\int_{S_{r}}K\;dA$.
Since $S_{r}$ is a compact two-dimensional Riemannian manifold by Gauss-Bonnet
$$\int_{S_{r}} K\;dA+\int_{\partial S_{r}}k_g\;ds=2\pi\chi(S_{r})$$
As $S_{r}$ is homeomorphic to a disc, $\chi(S_{r})=1$. The boundary $\partial S_{r}$ can be parametrized by the curve $\gamma(t)= (r\cos t,r\sin t,r^2)$.
The unit tangent vector is $T=(-\sin t,\cos t,0)$ and the inward-pointing unit normal to the boundary $\partial S_{r}$ on the surface $S_r$ at $\gamma(t)$ is $N=\frac {-1}{\sqrt{1+4r^2}}(\cos t,\sin t,2r)$.
Hence
$$\int_{\partial S_{r}}k_g\;ds=\int_{0}^{2\pi}\langle T',N\rangle \,dt=\int_{0}^{2\pi}\frac 1{\sqrt{1+4r^2}}\,dt=\frac {2\pi}{\sqrt{1+4r^2}}$$
which implies
$$\int_{S_{r}} K\;dA=2\pi\left(1-\frac{1}{\sqrt{1+4r^2}}\right)\;,\;\text{so}\;\int_{M_{r}}K\;dA=\frac{\pi}{2}\left(1-\frac{1}{\sqrt{1+4r^2}}\right)$$