Let $h(x+iy)$ be a harmonic function in the open neighbourhood of the closed unit disc $\overline\Delta(0;1)$ of $\mathbb{C}.$ Then it can be presented by Poisson integral formula in the following way: $$h(x+iy)=\frac{1}{2\pi}\int_{\partial\Delta(0;1)}\frac{1-x^2-y^2}{||x-iy-\zeta||^2}h(\zeta)d\sigma(\zeta),$$ where $\sigma$ is a uniformly distributed measure on the boundary of the disc.
I need partial derivatives of $h,$ so I will have to differentiate under the integral sign. How can I justify this operation?