After some time spent with the literature, I've decided to write an answer. Hopefully, this will sharpen my own thoughts and potentially help the others.

Currently a detailed calculation from the first principles of QCD is not available for the hadronic processes. Instead, one has to resort to some effective description. A common tool is the chiral perturbation theory (http://en.wikipedia.org/wiki/Chiral_perturbation_theory) -- effective theory constructed by taking the (approximate) chiral symmetry as well as other symmetries of QCD into account.

The standard chiral Lagrangian is symmetric w.r.t. replacement $\pi\to-\pi$ and therefore forbids processes with the odd number of $\pi$-mesons, like $\pi^0\to\gamma\gamma$ or the process of interest $\gamma\pi^+\to\pi^0\pi^+$. However, this $\pi\to-\pi$ symmetry is not a legitimate symmetry of QCD but only of the simple effective description. One can add a contribution to the chiral Lagrangian known as the Wess-Zumino term (http://www.sciencedirect.com/science/article/pii/037026937190582X) which violates this symmetry and permits the corresponding processes. This anomalous term is very special for it is not a local Lagrangian in 4d. Another side of this peculiarity is that in contrast to the other terms in the effective Lagrangian the anomalous term comes with the predefined coefficient (Witten http://www.sciencedirect.com/science/article/pii/0550321383900639).

Effective theory amplitudes are naturally expanded in powers of $E/\Lambda\ll1$ with the $E$ the energy scale of the process and $\Lambda$ the characteristic scale up to which the effective description holds. The lowest-energy amplitude for reaction $\gamma\pi^+\to\pi^0\pi^+$ is produced by a single tree diagram coming from the Wess-Zumino term and it is elementary to evaluate (for example see http://www.sciencedirect.com/science/article/pii/037026938490100X).

Unlike the case of $\pi^0\to\gamma\gamma$ the agreement between the theoretical prediction and the experimental value is not excellent: $A_{theor}=9.70\pm0.05 GeV^{-3},\quad A_{exp}=12.9\pm0.9\pm0.5 GeV^{-3}$. The discrepancy is believed to be due to disregarding the higher-energy corrections. However, the task of finding/estimating the corrections to the lowest-order result is not at all simple. One can adopt a variety of approaches on this way. Here my understanding falls rapidly and I will content myself with a reference to paper http://arxiv.org/abs/hep-ph/0107127 where combining a number of higher-energy effects (including loop corrections and previously overlooked electromagnetic corrections due to a single-photon exchange) the authors were able to bring the theory and experiment into agreement.

This summarizes my understanding. Any clarifications and corrections are welcome.