I'm afraid I didn't really understand any of your answers, so here is my own attempt. Please let me know whether I'm making sense.
It seems to me that we don't need to talk about rescaling the fields at all. We can produce the counterterms using only the canonically normalized fields throughout, as follows:
First, the ##\delta m## counterterm is an actual modification of the Lagrangian, on physical grounds: the self-interactions of a single particle will naturally add to its mass, but we need ##m## (used in ##H_0##) to equal the physical mass, so the only option is to explicitly subtract the added mass by assuming a new interaction with the form of a mass term (which happens to be infinite but whatever). We find the correct value for the coupling constant ##\delta m## by demanding that the pole in the propagator be at ##p^2=-m^2##.
Next, the LSZ theorem tell us that when using Feynman diagrams to calculate S-matrix elements,, without any field rescaling, we must multiply the result by a factor ##\sqrt{Z_i}## for each incoming and outgoing particle of species ##i##, where each ##\sqrt{Z_i}## is the inner product between a state created from the vacuum by the field operator - an eigenstate of the free Hamiltonian - and the corresponding one-particle state, an eigenstate of the full Hamiltonian. But in order to treat these factors perturbatively, we manipulate them as follows:
Instead of applying the factors to S-matrix elements or diagrams, we multiply each interaction vertex (including the ##\delta m## terms) by ##\sqrt{Z_i}## for every field of species ##i## appearing in the interaction, and divide every internal line of species ##i## by ##Z_i##. Each internal line meets two field operators at vertices, so all these factors will be cancelled, except for those for field operators that attach to the external lines. So the total effect is just to insert the factors for external lines, as needed.
Finally, we expand the ##Z^{-1}_i## factors in the propagators. Using a scalar field for simplicity:
$$\frac 1 Z = \frac 1 {1-(1-Z)} = \Sigma_{n=0}^\infty (1-Z)^n$$. Therefore,$$\frac 1 {Z(p^2+m^2-\epsilon i)}=$$ $$=\frac 1 {(p^2+m^2-\epsilon i)}+\frac 1 {(p^2+m^2-\epsilon i)}(1-Z)(p^2+m^2-\epsilon i)\frac 1 {(p^2+m^2-\epsilon i)}+...$$
So the final effect is the same as creating a new interaction ##(1-Z)(p^2+m^2-\epsilon i)## that can be inserted any number of times in any internal line.
One last detail is that in QED, it is customary to absorb the factor ##\sqrt{Z_3}## for the photon field in the interaction into a redefinition of the coupling constant ##e##, the "renormalized charge".
In summary, all the actual calculations can be justified without ever defining a "renormalized field". In particular, there is no reason to change ##H_0## at all, so my problem disappears.
I still don't see how to make sense of the standard presentations of this material, but if I have one version that works I am happy.