$\DeclareMathOperator\sgn{sgn}$For modular forms for $\mathrm{SL}_2(\mathbb{Z})$ there is an easy argument why the (essentially only) Eisenstein series has to be an eigenform of the Hecke operators. What I am now trying to see is that the same assertion is true for Hecke operators acting on modular forms for $\Gamma(N)$ and the Eisenstein series

$$E^{(0, 1)} (\tau) = \frac{1}{2} \sum_{\gcd(c, d)=1 ~\text{and}~ (c,d) \equiv (0,1) \mod N} {(c \tau + d)^{-k}}$$

(see for example Diamond/Shurman, A First course in Modular Forms, p.111).

So the question is: is this function an eigenform for the Hecke operators? (I have the feeling that this might only be true for some of the Hecke operators (see below))

Further remarks:

One can show that $E^{(0, 1)}$ is a modular form for $\Gamma(N)$ and there are the Hecke operators

$$T_p(f) = \sum_{j=0}^{p-1} f\big|_{\begin{pmatrix} 1 & jN \\ 0 & p \end{pmatrix}} + f\big|_{R_p \begin{pmatrix} p & 0 \\ 0 & 1 \end{pmatrix}}$$

(up to a constant) where $R_p$ is an arbitrary matrix in $\mathrm{SL}_2(\mathbb{Z})$ such that $R_p \equiv \begin{pmatrix} p^{-1} & 0 \\ 0 & p \end{pmatrix} \mod N$ (we always assume that $p$ is a prime with $\gcd(p,N)=1$).

One can show that the Fourier expansion of $E^{(0, 1)}$ is given by

$a_n(E^{(0, 1)}) = \sum_{a \in \mathbb{Z}_N^\times} \zeta^{a^{-1}} \sum_{m \in \mathbb{Z},~ m|n,~ n/m \equiv 0 \mod N} \sgn(m) m^{k-1} e(am)$

where $e(z) = e^{2 \pi i z}$, and $\zeta^{a} = \sum_{m \in \mathbb{Z}, m \equiv a \mod N} {\frac{\mu(m)}{m^k}}$ is the restricted zeta function associated to the Moebius function $\mu$. Let us assume for a moment that $p$ has a root $x$ mod N. Let us further assume that $\gcd(n,p)=1$, then the effect on the Fourier coefficients of the $p$-th Hecke operator is

$$a_n(T_p E^{(0,1)}) = a_{np}(E^{(0,1)}) = \sum_{a \in \mathbb{Z}_N^\times} \zeta^{a^{-1}} \sum_{m|np} (...)$$ Now we separate into cases $m|n$ and $m = pm'$ for some $m'|n$ giving $$\sum_{a \in \mathbb{Z}_N^\times} \zeta^{a^{-1}} \sum_{m|n} (\sgn(m)m^{k-1} e(am)) + \sum_{a \in \mathbb{Z}_N^\times} \zeta^{a^{-1}} \sum_{m|n} (\sgn(m)(mp)^{k-1} e(amp))$$

We substitute $a \mapsto ax$ in the first sum and $a\mapsto a/x$ in the second sum yielding

$$\sum_{a \in \mathbb{Z}_N^\times} (\zeta^{a^{-1}x^{-1}} + p^{k-1} \zeta^{a^{-1} x}) \sum_{m|n} (\sgn(m)m^{k-1} e(am))$$

Let us assume that $x \equiv 1 \mod N$ (i.e. also $p \equiv 1 \mod N$). Since the zeta function does only depend on the value modulo $N$, we have $\zeta^{a^{-1}x^{-1}} = \zeta^{a^{-1} x} = \zeta^{a^{-1}}$ so that the above is nothing else than $(1 + p^{k-1}) a_n(E^{(0,1)})$. Since the case $\gcd(n,p)\neq 1$ works similarly, this really proves that $E^{(0,1)}$ is an eigenform whenever $p \equiv 1 \mod N$. In the case I am interested in, $p$ still has a root $x$ modulo $N$ but it is not congruent to one. Does anybody see/know how to handle this case? (Is it even true?)