The Riemann zeta function is given for $Re s>1$
$$\zeta(s) = \prod\limits_{p \; prime} ( 1-p^{-s}) $$
This product converges absolutely in $Re s >1$, hence it does not vanish in $Re s>1$. Actually the product also converges locally uniformly, which implies that $\zeta(s)$ is holomorphic for $Re s>1$.
The functional equation follows from the Poisson summation formula, which is a Fourier theoretic argument. The functional equation is given by
$$ \Lambda(s) = \pi^{-s/2} \Gamma(s/2) \zeta(s) = \Lambda(1-s).$$
This implies that all zeros lie inside $0 \leq Re s \leq 1$. We will need to evaluate a certain integral
$$ * = \int\limits_{Re\; s = 1 + \epsilon/2 \cup Re \; s = -\epsilon/2} G(s) \frac{\Lambda'}{\Lambda} (s) d s$$
We have now for a meromorphic function $F$ and a holomorphic function $G$ in some region $G$ containing the closure of a simply connected region $O$, that the contour integral
$$\frac{1}{2 \pi i} \int\limits_{\partial O} G(z) \frac{F'}{F}(z) d z = \sum\limits_{\rho \; zero \; of \; F \; in \; O} G( \rho) - \sum\limits_{\nu \; pol \; of \; F \; in \; O} G( \nu),$$
where $\partial O$ denotes the boundary of $O$ and is a Jordan curve by assumption. For this identity to be valid $F$ must have no zero and no pol on the boundary of $O$.
This is what I call the weighted argument principle and can be derived in the same lines as the argument principle. It follows for entire functions easily from the Hadamard factorization theorem, but is a purely local property.
Apply this to the function $F = \Lambda$ and $G$ being holomorphic in $ - \epsilon \leq Re\; s \leq 1$ with $\epsilon>0$ and certain restrictions of the growth. We choose the contour $C=C(T)$ being the boundary of $ - \epsilon/2 < Re\;s < 1 + \epsilon/2$ and $| Im \;s | \leq T$, where \zeta does not vanish on $Im \; s = \pm T$.
This give an expression for $*$ for $T \rightarrow \infty$ involving the nontrivial zeros of $\zeta$.
Using the Euler product, we can also derive a nice explicit expression
$$ \frac{\Lambda'}{\Lambda} (s) = (log \Lambda(s))'= -1/2 \log \pi +\frac{1}{2} \frac{\Gamma'}{\Gamma}(s/2) - \sum\limits_{p \; prime} \frac{p^{-s} \log p }{1-p^{-s}}.$$
This gives an expression for $*$ for $T \rightarrow \infty$ involving the primes.
Assuming certain boundedness conditions on $\Lambda(s)$ and $G(s)$ in $-\epsilon \leq Re \; s \leq 1+\epsilon$, we are actually allowed to choose $C(T)$ with $T \rightarrow \infty$ and derive the formula as the limit. The boundedness conditions for $\Lambda$ follow from the Hadamard three lines principle or the Phragmen Lindeloeff principle (this is not entire function theory, but this only a complex analysis argument), then the explicit formula follows by choosing $G$ appropiately.
Remark: Inserting a Gaussian function into the explicit formula allows to derive the functional equation for $\zeta$, hence the functional equation is equivalent to Weil's explicit formula. Actually Samuel Patterson states in his famous book on $\zeta$ that they are actually both equivalent to the Poisson summation formula, but I do not know how? Of course the Poisson summation implies Functional equation of $\zeta$... How to go back?
So the philosophy is: Functional equation + Euler product = explicit formula. Another example for this is the relation of the Selberg Zeta function and Selberg trace formula.
You are right that the entire function theory implies that $\Lambda$ has necessary many zeros, but not too many, since it is of exponential type $1$ because of the factor $\Gamma$. If you want to derive this without using the merophorphicity of $\Lambda$, you might want to try to deduce this without knowledge over the primes and by inserting an approppiate chosen test function in the explicit formula. I have never seen this been worked out, but choosing an appropiate function $g$ being supported in $ - \log 2 < t < \log 2$ (so no contribution by finite primes) etc. should lead to a rough asymptotic of the zeros without any information used about the primes, but possibly a weaker error term than in the classical van Mangoldt estimate, which I expect to be a square root of the actual main term. Look at similiar techniques used in Werner Mueller and Erez Lapid's article Chapter 2 of http://www.math.uni-bonn.de/people/mueller/papers/orbint09.11.pdf for the Weyl law. The Selberg trace formula has many analogies with the explicit formula. In Iwaniec - Spectral methods in automorphic forms, you can find an argument using Tauberian theorems, which is weaker, but you get only the main term and no error term.
One interesting, but technical derivation of the explicit formula using only the languages of the adeles, harmonic analysis and no entire function theory at all was given by Ralf Meyer:
http://arxiv.org/pdf/math/0311468v3
However, the Fourier transform of a function with certain growth properties has some holomorphicity conditions, so you basically will just hide the complex analysis arguments, but be able to deduce the results from Fourier analysis only. If you want to derive the prime number theorem from Fourier analysis, you might want to consider the treatement in Rudin - Functional analysis, which is based on real analysis only. Also there are elementary proofs of the prime number theorem, which I have no idea if they apply to Beurling primes.
See Chapter 15 ("Oscillation of error terms") in Montgomery-Vaughan: Multiplicative number theory I. See especially Theorems 15.2-15.3 and 15.11.
Added by Steven Clark and GH from MO. For convenience, we copied below the relevant theorems and some additional text from the book. As usual,
$$M(x):=\sum\limits_{n\le x}\mu(n)$$
is the Mertens function.
Theorem 15.2. Let $\Theta$ denote the supremum of the real parts of the zeros of the zeta function. Then for every $\varepsilon>0$,
$$\psi(x)-x=\Omega_\pm(x^{\Theta-\varepsilon})\tag{15.1}$$
and
$$\pi(x)-\mathrm{li}(x)=\Omega_\pm(x^{\Theta-\varepsilon})\tag{15.2}$$
as $x\to\infty$.
Theorem 15.3. Suppose that $\Theta$ is the supremum of the real parts of $\zeta(s)$, and there is a zero $\rho$ with $\Re\rho=\Theta$, say $\rho=\Theta+i\gamma$. Then
$$\underset{x\to\infty}{\text{lim sup}}\ \frac{\psi(x)-x}{x^{\Theta}}\geq \frac{1}{| \rho |}\tag{15.4}$$
and
$$\underset{x\to\infty}{\text{lim inf}}\ \frac{\psi(x)-x}{x^{\Theta}}\leq -\frac{1}{| \rho |}.\tag{15.5}$$
Theorem 15.11. As $x\to\infty$,
$$\psi(x) -x = \Omega_{\pm} \bigl(x^{1/2} \log \log \log x\bigr),$$
and
$$\pi(x) - \mathrm{li}(x) = \Omega_{\pm} \bigl(x^{1/2}(\log x)^{-1}\log \log \log x\bigr).$$
Then in the manner of the proof of Theorem 15.3, we find that if $\Theta+i\gamma$ is a zero of $\zeta(s)$, then
$$\underset{x\to\infty}{\text{lim sup}}\ \frac{M(x)}{x^{\Theta}}\geq \frac{1}{| \rho\,\zeta'(\rho) |},\tag{15.11}$$
and
$$\underset{x\to\infty}{\text{lim inf}}\ \frac{M(x)}{x^{\Theta}}\leq -\frac{1}{| \rho\,\zeta'(\rho) |}.\tag{15.12}$$
Best Answer
The analogy between zeta functions and trace formulas goes at least to Selberg, when he proved his famous trace formula for hyperbolic surfaces and the result turned out to resemble Weil's generalization of Riemann's explicit formula quite a bit.
Riemann-Weil explicit formula:
\begin{equation*} \begin{split} \sum_\gamma h(\gamma)=\frac{1}{2\pi}\int_{-\infty}^{+\infty} h(r) \frac{\Gamma'}{\Gamma}\left(\frac{1}{4}+\frac{1}{2}ir \right)dr &+h\left(\frac{i}{2}\right)+h\left(-\frac{i}{2}\right)\\ &-g(0)\ln\pi-2\sum_{n=1}^\infty\frac{\Lambda(n)}{\sqrt{n}}g(\ln n) \end{split} \end{equation*} Selberg trace formula:
\begin{equation*} \begin{split} \sum_{n=0}^\infty h(r_n)=\frac{\mu(F)}{4\pi}\int_{-\infty}^{+\infty} rh(r) \tanh (\pi r)dr &+\sum_n \Lambda(n) g(\ln N(n)) \end{split} \end{equation*}
There's a big literature on this type of question, and the interplay of number theory, spectral analysis, mathematical physics... I recomend section 3 of Lagarias' survey "The Riemann Hypothesis: Arithmetic and Geometry" for references.
Connes approach started with
And was completed (as far as I know) in
The main result says that given a global field $K$ and a character $\alpha=\prod_v \alpha$ of the space of adele classes $A/K$, and any adecuate test function $h$, we have:
$$\underbrace{\widehat{h}(0)+\widehat{h}(1)-\sum \widehat{h}(\mathcal{X},\rho)}_{\text{spectral side}}=\underbrace{\sum_v \int_{K_v^*}' \frac{h(u^{-1})}{|1-u|}d^*u}_{\text{arithmetic side}}$$