What makes coming up with a mathematically solid, non-shaky relativistic quantum field theory (RQFT) so hard? This is something I know of but I'm not quite sure I understand the details. Particularly, when it comes to interacting RQFTs, such as even QED, where some posts here have pointed out that it cannot make sense even of a hydrogen atom without some further assumptions. To me that isn't a particularly nice situation, especially for a theory that one would believe from all the hype and the casual talk is supposed to be a "well-trodden area" and not something at the frontier like quantum gravity. And something I think makes sense to want to remedy.
But what I am curious about is, why exactly is this so difficult to do? I've heard of Haag's theorem, but only have access to the Wikipedia terse description so I don't quite get 100% what it's after. I also know of how to make the simple and ubiquitous free-field (R)QFT via the Fock space construction and how you can solve for phonons on a crystal lattice, and so what I don't get is this:

Why can't you just tensor together the two Hilbert (Fock) spaces of the free EM (photon) and free charge (electron/positron) fields, and then write a suitable interaction Hamiltonian?

That is, after all, the way you do interacting systems in NRQM. What makes this naive approach fail? Almost surely this was the first thing tried many decades ago, but can a fairly decent summary of the arguments against it be given?
 A: 
Why can't you just tensor together the two Hilbert (Fock) spaces of the free EM (photon) and free charge (electron/positron) fields, and then write a suitable interaction Hamiltonian?

The quick mathematical answer to this is that when you do that, and try to write the relativistic interaction term as an operator, it does not make sense.
It does not even make sense weakly as a quadratic form or on some subset of vectors.
To overcome this problem, one should first regularize the interaction. In doing that however, often the relativistic covariance is broken. The aim is then to remove the regularizations, restoring covariance, and obtaining a well-defined operator. In order to do that, some infinities in the wavefunctions and energy shall be taken care of.
There are several ways to do that, that have been successful to some extent (scalar $\phi^4$ theory can be defined rigorously in 1 and 2 space dimensions for example). However, most of the strategies even "forget" about the Hamiltonian, since in relativistic theories the generator of time translations does not play a preferred role. Quadratic scalar theories and the Yukawa theory in 2 spatial dimensions have been however defined rigorously using the Hamiltonian approach (works by Ginibre and Velo, and Glimm, respectively).
