I want to solve the integral
\begin{align} I=\int_{-\infty}^{\infty}\frac{1}{e^{\frac{x-\mu}{T}}+1}\frac{\gamma}{(x-x_{0})^{2}+\frac{\gamma^{2}}{4}}dx \end{align}
where $T,\mu,\gamma,x_{0}\in\mathbb{R}_{>0}$. Using partial fraction expansion for the Cauchy distribution \begin{align} \frac{\gamma}{(x-x_{0})^{2}+\frac{\gamma^{2}}{4}}=i\left[\frac{1}{x-(x_{0}-i\frac{\gamma}{2})}-\frac{1}{x-(x_{0}+i\frac{\gamma}{2})}\right] \end{align} those are simple poles and easy to deal. But I don't have very clear how to proceed with the poles coming from the Fermi distribution \begin{align} e^{\frac{x-\mu}{T}}+1=0\iff x=in\pi T+\mu\quad n\in\mathbb{Z}. \end{align}
How should I continue to solve the integral?
Using the substitution $x = \mu + T u$ and the abbreviations $u_0 = \frac{x_0 - \mu}{T}$ and $\beta = \frac{\gamma}{2T}$, we can write your integral as $$ f \colon (0,\infty) \times \mathbb{R} \, , \, f (\beta, u_0) = 2 \beta \int \limits_{-\infty}^\infty \frac{\mathrm{d} u}{(\mathrm{e}^{u} + 1)[(u-u_0)^2 + \beta^2]} \equiv 2 \beta \int \limits_{-\infty}^\infty g_{\beta,u_0}(u) \, \mathrm{d} u \, . $$ $g_{\beta,u_0}$ has simple poles at $u_0 \pm \mathrm{i} \beta$ with residues $$ \operatorname{Res}(g_{\beta,u_0},u_0 \pm \mathrm{i} \beta) = \pm \frac{1}{2 \mathrm{i} \beta (\mathrm{e}^{u_0 \pm \mathrm{i} \beta} + 1)}$$ and at $\pm (2n+1) \pi \mathrm{i}$ (so only at the odd-integer multiples of $\pi \mathrm{i}$ !) with residues $$ \operatorname{Res}(g_{\beta,u_0},\pm (2n+1) \pi \mathrm{i}) = - \frac{1}{[(2n+1) \pi \mathrm{i} \mp u_0]^2 + \beta^2}$$ for $n \in \mathbb{N}_0$ . The integrals of $g_{\beta,u_0}$ along semi-circles (avoiding the poles on the imaginary axis) vanish in the limit of large radii, so we can use the residue theorem to evaluate the integral.
Closing the contour in the upper half-plane (the lower half-plane works just as well) yields $$ f(\beta,u_0) = 2 \beta \, 2 \pi \mathrm{i} \left[\operatorname{Res}(g_{\beta,u_0},u_0 + \mathrm{i} \beta) + \sum \limits_{n=0}^\infty \operatorname{Res}(g_{\beta,u_0},(2n+1) \pi \mathrm{i}) \right] \, .$$ Now we plug in the values of the residues, perform a partial fraction decomposition in the infinite sum and use the series formula for the digamma function $\psi$ to find $$ f(\beta,u_0) = 2 \pi \left[\frac{1}{\mathrm{e}^{u_0 + \mathrm{i} \beta} + 1} - \frac{1}{2 \pi \mathrm{i}} \left(\psi \left(\frac{1}{2} + \frac{\beta + \mathrm{i} u_0}{2 \pi}\right) - \psi \left(\frac{1}{2} + \frac{-\beta + \mathrm{i} u_0}{2 \pi}\right) \right)\right] \, .$$ Finally, we apply the reflection formula to the second digamma function and simplify the result: \begin{align} f(\beta,u_0) &= 2 \pi \left[\frac{1}{\mathrm{e}^{u_0 + \mathrm{i} \beta} + 1} - \frac{1}{2 \mathrm{i}} \tan\left(\frac{\beta - \mathrm{i} u_0}{2}\right)\right. \\ &\phantom{222222}- \left.\frac{1}{2 \pi \mathrm{i}} \left(\psi \left(\frac{1}{2} + \frac{\beta + \mathrm{i} u_0}{2 \pi}\right) - \psi \left(\frac{1}{2} + \frac{\beta - \mathrm{i} u_0}{2 \pi}\right) \right)\right] \\ &= 2 \pi \left[\frac{1}{2} - \frac{1}{\pi} \operatorname{Im} \left(\psi \left(\frac{1}{2} + \frac{\beta + \mathrm{i} u_0}{2 \pi}\right) \right) \right] = \pi - 2 \operatorname{Im} \left[\psi \left(\frac{1}{2} + \frac{\beta + \mathrm{i} u_0}{2 \pi}\right) \right] \, . \end{align}
Returning to the original parameters, we end up with $$ \int \limits_{-\infty}^\infty \frac{1}{\mathrm{e}^{(x - \mu)/T} + 1} \frac{\gamma}{(x-x_0)^2 + \gamma^2/4} \, \mathrm{d} x = \pi - 2 \operatorname{Im} \left[\psi \left(\frac{1}{2} + \frac{\frac{\gamma}{2} + \mathrm{i} (x_0 - \mu)}{2 \pi T}\right) \right] \, . $$ Note that in the special case $x_0 = \mu$ the result is simply $\pi$. This can also be shown using elementary methods, so the implicit assumption $u_0 + \mathrm{i} \beta \not\in (2 \mathbb{Z} + 1) \pi \mathrm{i}$ used in the computation of the residues is justified.