I think the following articles can give a clue: http://www.tandfonline.com/doi/abs/10.1080/10652460701688125?journalCode=gitr20 (Closed-form evaluation of some families of cotangent and cosecant integrals, by D. Cvijović) and http://www.sciencedirect.com/science/article/pii/S0377042702003588 (Integral representations of the Riemann zeta function for odd-integer arguments, by D. Cvijović and J. Klinowski).
Update. In particular, using $2x^n=E_n(x)+E_n(1+x)$ and $E_n(1+x)=\sum\limits_{k=0}^n\binom{n}{k}E_k(x)$, we get $$f(n)=\frac{\pi^{n+2}}{2}\left[I_n+\sum\limits_{k=0}^n\binom{n}{k}I_k-I_{n+1}-\sum\limits_{k=0}^{n+1}\binom{n+1}{k}I_k\right]=\frac{\pi^{n+2}}{2}\left[I_n-2I_{n+1}-\sum\limits_{k=1}^n\binom{n}{k-1}I_k\right], \tag{1}$$ where $$I_n=\int\limits_0^1E_n(x)\csc{(\pi x)}\,dx.$$ I suspect (however have not proved it) that odd indices in (1) sum to zero because lead to the integrand which is antisymmetric with regard to $x\to 1-x$. For even indices, $I_n$ was calculated in the second paper indicated above as $$I_{2n}=(-1)^n\frac{(4-2^{1-2n})(2n)!}{\pi^{2n+1}}\zeta(2n+1),$$ and we obtain, for example, for even $n=2m$: $$f(2m)=(-1)^m\frac{\pi}{2}(4-2^{1-2m})(2m)!\,\zeta(2m+1)- \sum\limits_{p=1}^m(-1)^p\binom{2m}{2p-1}\frac{\pi^{2(m-p)+1}}{2}(4-2^{1-2p})(2p)!\,\zeta(2p+1).$$