(Informal discussion on singularity) The $v^{ij}$ and $w$ should have very specific singularities along $\mathfrak{D}_{1}$, $\mathfrak{D}_{2}$, $\mathfrak{D}_{3}$. Let us focus on what happens around $\mathfrak{D}_{1}$. The delta forcing term appears in the component of (2.4) as

$$\frac{\sqrt{-1}}{4\pi}\left(\frac{\partial^{2}w}{\partial\mu_{1}\partial\mu_{1}}+4\frac{\partial^{2}v^{11}}{\partial\eta\partial\bar{\eta}}\right)d\mu_{1}\wedge d\eta\wedge d\bar{\eta}=\mathfrak{D}_{1}.$$

If we denote the Lebesgue measure $f\mapsto\int_{\mathfrak{D}_{1}}fd\mu_{2}$ as $\int_{\mathfrak{D}_{1}}d\mu_{2}$, we may rewrite this equation as

$$\frac{1}{2\pi}\left(\frac{\partial^{2}w}{\partial\mu_{1}\partial\mu_{1}}+4\frac{\partial^{2}v^{11}}{\partial\eta\partial\bar{\eta}}\right)d\mu_{1}\wedge d\mu_{2}\wedge d\text{Re}\eta\wedge d\text{Im}{\eta}=-\int_{\mathfrak{D}_{1}}d\mu_{2}.$$

Since $v^{12},v^{22}$ do not see the forcing term, our best guess is that they are smooth along $\mathfrak{D}_{1}$. Then modulo smooth terms $w\sim Aa^{11}v^{11}=a_{22}v^{11}$ along $\mathfrak{D}_{1}$, from which the distributional equation gives the singularity structure along $\mathfrak{D}_{1}$:

$$v^{11}\sim\frac{1}{2\sqrt{\mu_{1}^{2}+a_{22}|\eta|^{2}}},\quad w\sim\frac{a_{22}}{2\sqrt{\mu_{1}^{2}+a_{22}|\eta|^{2}}}.$$

The following base metric encoding the Gibbons-Hawking data has the singularity structure along $\mathfrak{D}_{1}$

$$\begin{split}&g_{a}+v^{ij}d\mu_{i}d\mu_{j}+w|d\eta|^{2}\sim g_{a}+\frac{1}{2\sqrt{\mu_{1}^{2}+a_{22}|\eta|^{2}}}(d\mu_{1}^{2}+a_{22}|d\eta|^{2})\\ =&(\frac{1}{2\sqrt{\mu_{1}^{2}+a_{22}|\eta|^{2}}}+\frac{A}{a_{22}})(d\mu_{1}^{2}+a_{22}|d\eta|^{2})+a_{22}(d(\mu_{2}+\frac{a_{12}}{a_{22}}\mu_{1}))^{2}\end{split}$$

from which we recognize the Taub-NUT metric appearing in directions transverse to $\mathfrak{D}_{1}$. See Section 2.3 for further details.