Let $f$ be a holomorphic modular form (of given weight and level one). Since it $1$-périodic and of moderate growth, it has a Fourier expansion, but this one is a Fourier expansion in $x$, that is $$f(z) = \sum_{n \geqslant 0} a_n(y) e(nx).$$
I would like to understand why $a_n(y)$ does not depend on $y$. Many sources just state the holomorphicity, but I do not see the relation with this integral? Is that a question of contour integral expression?
On
As you have written it, the coefficients do indeed depend on $y$
$$f(x+iy) = \sum_{n = -m}^{+\infty} a_n(y) \, e^{2i\pi n x}$$
But what remains to be shown is that the coefficients are of the form $a_n(y) = b_n e^{-2\pi n y}$ with $b_n = a_n(0)$ a constant (which will allow us to write $f$ as $f(z) = \sum_{n = -m}^{+\infty} b_n \, e^{2i\pi n z}$). First, recall from Fourier analysis that we have
$$a_n(y) = \int_0^{1} f(t + iy) \, e^{-2 i \pi n t} dt$$
Now comes the complex analysis trick : we can consider $a_n$ as a complex function :
$$a_n(z) = \int_0^{1} f(t + iz) \, e^{-2 i \pi n t} dt$$
And see that $a_n$ is holomorphic because $f$ is. Now, when $z=ix \in i\mathbb{R}$, we can compute
$$a_n(ix) = \int_0^{1} f(t - x) \, e^{-2 i \pi n t} dt \underset{u = t - x}{=} \int_{-x}^{1-x} f(u) \, e^{-2 i \pi n (u+x)} du = a_n(0) \, e^{-2i \pi n x}$$
Since the holomorphic functions $a_n$ and $z \mapsto a_n(0) e^{-2i \pi n z}$ are equal on $i\mathbb{R}$, they are equal everywhere, and in particular for any $y \in \mathbb{R}$ we have
$$a_n(y) = a_n(0) \, e^{-2i \pi n y}$$
Here is a presentation that could shed some light on a more general feature about Fourier coefficients of Maass forms, and can be used in a broader setting.
The case of holomorphic forms. If $f$ is $1$-periodic, it has a Fourier expansion $$f(x) = f(x+iy) = \sum_{n \in \mathbb{Z}} a_n(y) e^{2i\pi n x}.$$
This could be restricted to some $n$ depending on the conditions on $f$ (typically, moderate growth implies $n \geqslant 0$ and cuspidality $n>0$). Holomorphicity means by definition $$\frac{\partial f}{\partial \bar z} = \frac{1}{2} \left( \frac{\partial f}{\partial x} + i\frac{\partial f}{\partial y}\right) = 0.$$
Inputing it in the Fourier expansion and identifying the coefficients leads to the differential equation $$a_n'(y) + 2\pi n a_n(y) = 0,$$ so that there is a constant $a_n$ such that $$a_n(y) = a_n e^{-2\pi n y} = a_n e^{2i\pi n (iy)}.$$
We are done. Why is this way interesting? As I announced before, it works in a broader setting.
The case of Maass forms. This is another example of the philosophy above. Maass forms $\phi$ are automorphic functions that satisfy a certain invariance (say by $SL_2(\mathbb{Z})$, so that it is also $1$-periodic) and is an eigenvalue of the weight $k$ Laplacian (see below)
Since it is $1$-periodic, it has a Fourier expansion of the form $$\phi(x) = \phi(x+iy) = \sum_{n \in \mathbb{Z}} a_n(y) e^{2i\pi n x}.$$
More precisely, the differential equation satisfied by $\phi$ is $$\Delta_k \phi = -y^2 \left( \frac{\partial^2 \phi}{\partial x^2} + \frac{\partial^2 \phi}{\partial y^2} \right) + ik y \frac{\partial \phi}{\partial x} = \nu(1-\nu) \phi.$$
Inputing it in the Fourier expansion and identifying the coefficients leads to the differential equation, we are led to $$\Delta_k \left( a_n(y) e^{2i\pi n x}\right) = \nu(1-\nu) a_n(y) e^{2i\pi n x}.$$
This is the so-called Whittaker differential equation. It is a second order differential equation, therefore the space of solutions is two-dimensional. But the space of solutions with moderate growth (also part of the definition of a Maass form) is only one dimensional, and is spanned by the Whittaker function $W_{\mathrm{sgn}(n)k}$. We find finally that $$a_n(y) = a_n W_{\mathrm{sgn}(n)k} (4\pi n y).$$
The idea. As soon as $f$ (is $1$-periodic and) satisfies a differential equation, this differential equation translates on the Fourier coefficients $a_n(y)$, that forces $a_n(y)$ to take a certain shape (and automorphic forms, by definition, do satisfy certain partial differential equations).