Relation between BPHZ renormalization and on-shell renormalization

I am following Sidney Coleman’s lecture of Quantum field theory where in ch 25 he talked BPHZ renormalization and in ch 23 he talked about on-shell renormalization of spin-1/2 theory.

Let us take a theory of Fermi field coupled to a real pseudoscalar field through Yukawa interaction

$mathcal{L} = mathcal{L}_{free} + bar{psi} ig gamma_{5} phi psi$ where $bar{psi}$ and $psi$ are fermi fields, $phi$ is the scalar field and $mathcal{L}_{free}$ is the free part of the Lagrangian.

First, if we try to do on-shell renormalization on the fermion, we will need to impose conditions such as $Sigma(m) = 0$ and $d Sigma(p) / d p = 0$ when $p = m$. $Sigma(p)$ is the self-energy of the fermions, and $m$ is the "physical" mass of the particle. This allows us to compute the counterterms order by order. The reason we impose these conditions is to fix the position of the pole at $p = m$ of the exact propagator, and the residue of the exact propagator also at $p=m$. Therefore, on-shell renormalization makes reference to experimental results, such as physical masses. We take them as granted, and use them to fix our counterterm.

Next, if we try to use BHPZ renormalization, we will make no reference to experimental results such as physical masses of the particles. What we care about in this scheme is to just add whatever counterterms we need to cancel one-particle-irreducible (1PI) diagrams that are divergent. A theory is renormalizable if we can add finite number of terms in the Lagrangian to cancel any divergent 1PI diagrams. A theory is strictly renormalizable if all the counterterms we add to cancel the divergence take the same form of the terms that are already in the Lagrangian.

Now, these two schemes both serve as the method to eliminate the divergence that might come from large momentum when we do momentum integrals. But I am confused about how these two schemes compared with each other. As a concrete question, one thing that we will be interested in is the S-matrix element, which by LSZ reduction formula is related to the Green’s function of the renormalized fields with on-shell external momenta, to wit

$mathcal{A}(p_{1},p_{2} to p_{1}’,p_{2}’) = (-i)^{4} prod_{r=1}^{4} (k_{r}^{2} – mu_{r}^{2}) G^{(4)}(-k_{3},k_{4},k_{1},k_{2})$

where $mu_{r}$ is the mass of the $r^{th}$ particle and $G^{(4)}(-k_{3},k_{4},k_{1},k_{2})$ is the 4 point Green function for the "renormalized" fields.

Now, I can understand how the on-shell scheme can get me to compute this S-matrix element. However, I am not quite sure how BPHZ scheme can be consistent with the on-shell scheme. Since we just add term to eliminate divergence in a (seemingly) brute force way. And I wonder how does this scheme can (or can not) be used to compute any things that can be measured experimentally. Actually I don’t quite know how I can use BPHZ to compute things like S-matrix element.

I will be happy to have any suggestions, including perhaps some books that deal with this carefully and provide some working examples to let us see what is the difference between these two schemes. Thanks!

P.S. A quick question: is "minimal subtraction scheme" the same thing as "BPHZ scheme"?