Please view this as an supplement of Olivier's answer.
I will derive a sufficient condition on the meromorphic function involved which
allow one to apply a result similar to that in Olivier's answer.
Let $\phi(z)$ be any meromorphic function over $\mathbb{C}$ which
preserve the extended real line $\mathbb{R}^* = \mathbb{R} \cup \{ \infty \}$ in the sense:
$$\begin{cases}\phi(\mathbb{R}) \subset \mathbb{R}^*\\ \phi^{-1}(\mathbb{R}) \subset \mathbb{R}\end{cases}
\quad\implies\quad
P \stackrel{def}{=} \phi^{-1}(\infty) = \big\{\, p \in \mathbb{C} : p \text{ poles of }\phi(z)\,\big\} \subset \mathbb{R}
$$
Split $\mathbb{R} \setminus P$ as a countable union of its connected components $\,\bigcup\limits_{n} ( a_n, b_n )\,$. Each connected component is an open interval $(a_n,b_n)$
and on such an interval, $\phi(z)$ increases from $-\infty$ at $a_n^{+} $ to $\infty$
at $b_n^{-}$.
There exists an ascending chain of Jordan domains $D_1, D_2, \ldots$ that cover $\mathbb{C}$,
$$\{ 0 \} \subset D_1 \subset D_2 \subset \cdots
\quad\text{ with }\quad \bigcup_{k=1}^\infty D_k = \mathbb{C}
$$
whose boundaries $\partial D_k$ are "well behaved", "diverge" to infinity and $| z - \phi(z)|$ is bounded on the boundaries. More precisely, let
$$
\begin{cases}
R_k &\stackrel{def}{=}& \inf \big\{\, |z| : z \in \partial D_k \,\big\}\\
L_k &\stackrel{def}{=}& \int_{\partial D_k} |dz| < \infty\\
M_k &\stackrel{def}{=}& \sup \big\{\, |z - \phi(z)| : z \in \partial D_k \,\big\}
\end{cases}
\quad\text{ and }\quad
\begin{cases}
\lim\limits_{k\to\infty} R_k = \infty\\
\lim\limits_{k\to\infty} \frac{L_k}{R_k^2} = 0\\
M = \sup_k M_k < \infty
\end{cases}
$$
Given such a meromorphic function $\phi(z)$ and any Lebesgue integrable function $f(x)$ on $\mathbb{R}$, we have following identity: $$
\int_{-\infty}^\infty f(\phi(x)) dx = \int_{-\infty}^\infty f(x) dx \tag{*1}
$$
In order to prove this, we split our integral into a sum over the connected components
of $\mathbb{R} \setminus P$.
$$\int_\mathbb{R} f(\phi(x)) dx
= \int_{\mathbb{R} \setminus P} f(\phi(x)) dx
= \sum_n \int_{a_n}^{b_n} f(\phi(x)) dx
$$
For any connected component $( a_n, b_n )$ of $\mathbb{R} \setminus P$ and $y \in \mathbb{R}$, consider the roots of the equation $\phi(x) = y$.
Using properties $(1)$ and $(2)$ of $\phi(z)$, we find there is a unique root
for the equation $y = \phi(x)$ over $( a_n, b_n )$. Let we call this root as $r_n(y)$.
Change variable to $y = \phi(x)$, the integral becomes
$$\sum_n \int_{-\infty}^\infty f(y) \frac{d r_n(y)}{dy} dy
= \int_{-\infty}^\infty f(y) \left(\sum_n \frac{d r_n(y)}{dy}\right) dy
$$
We can use the obvious fact $\frac{d r_n(y)}{dy} \ge 0$ and dominated convergence theorem to justify the switching of order of summation and integral.
This means to prove $(*1)$, one only need to show
$$\sum_n \frac{d r_n(y)}{dy} \stackrel{?}{=} 1\tag{*2}$$
For any $y \in \mathbb{R}$, let $R(y) = \phi^{-1}(y) \subset \mathbb{R}$ be the collection of roots of the equation $\phi(z) = y$.
Over any Jordan domain $D_k$, we have following expansion
$$\frac{\phi'(z)}{\phi(z) - y} = \sum_{r \in R(y) \cap D_k} \frac{1}{z - r} - \sum_{p \in P \cap D_k} \frac{1}{z - p} + \text{something analytic}$$
This leads to
$$\sum_{r \in R(y)\cap D_k} r - \sum_{ p \in P \cap D_k} p
= \frac{1}{2\pi i}\int_{\partial D_k} z \left(\frac{\phi'(z)}{\phi(z) - y}\right) dz$$
As long as $R(y) \cap \partial D_k = \emptyset$, we can differentiate both sides and get
$$\begin{align}
\sum_{r_n(y) \in D_k} \frac{dr_n(y)}{dy}
&=
\frac{1}{2\pi i}\int_{\partial D_k} z \left(\frac{\phi'(z)}{(\phi(z) - y)^2}\right) dz
= -\frac{1}{2\pi i}\int_{\partial D_k} z \frac{d}{dz}\left(\frac{1}{\phi(z)-y}\right) dz\\
&= \frac{1}{2\pi i}\int_{\partial D_k}\frac{dz}{\phi(z) - y}
\end{align}
$$
For those $k$ large enough such that $R_k > 2(M+|y|)$, we can expand the integrand in last line as
$$\frac{1}{\phi(z) - y} = \frac{1}{z - (y + z - \phi(z))}
= \frac{1}{z} + \sum_{j=1}^\infty \frac{(y + z - \phi(z))^j}{z^{j+1}}$$
and obtain a bound
$$\left|\left(\sum_{r_n(y) \in D_k} \frac{dr_n(y)}{dy} \right) - 1\right|
\le \frac{1}{2\pi}\sum_{j=1}^\infty \int_{\partial D_k} \frac{(|y| + |z-\phi(z)|)^j}{|z|^{j+1}} |dz|\\
\le \frac{(M + |y|)L_k}{2\pi R_k^2}\sum_{j=0}^\infty\left(\frac{M+|y|}{R_k}\right)^j
\le \frac{M + |y|}{\pi}\frac{L_k}{R_k^2}
$$
Since $\lim\limits_{k\to\infty} \frac{L_k}{R_k^2} = 0$, this leads to
$$\sum_n \frac{dr_n(y)}{dy} = \lim_{k\to\infty} \sum_{r_n(y) \in D_k} \frac{dr_n(y)}{dy} = 1$$
This justifies $(*2)$ and hence $(*1)$ is proved. Notice all the $\frac{dr_n(y)}{dy}$ are positive, there is no issue in rearranging the order of summation in last line.
Back to the original problem of evaluating
$$\int_{-\infty}^\infty \frac{1}{1+(x+\tan x)^2} dx$$
One can take $\phi(z)$ as $z + \tan z$ and $f(x)$ as $\frac{1}{1+x^2}$.
It is easy to see $\phi(z)$ satisfies:
Condition $(1)$ - For any $y \in \mathbb{R}$ and $u + iv \in \mathbb{C} \setminus \mathbb{R}$, we have
$$\begin{align}
\Im (\phi(u + iv) - y )
&= v + \Im\tan(u+iv) = v + \Im\frac{\tan u + i\tanh v}{1 - i\tan u\tanh v}\\
&= v + \tanh v\frac{1 + \tan^2 u}{1 + \tan^2 u\tanh^2 v} \ne 0
\end{align}$$
Condition $(2)$ - obvious.
Condition $(3)$. - Let $D_k$ to be the square
$$D_k = \big\{\, u + v i \in \mathbb{C} : |u|, |v| \le k \pi \,\big\}$$
It is not hard to show $|z - \phi(z)| = |\tan z|$ is bounded above by $\frac{1}{\tanh k\pi}$ on $\partial D_k$.
Combine these, we can apply $(*1)$ and deduce
$$ \int_{-\infty}^\infty \frac{1}{1+(x+\tan x)^2} dx
= \int_{-\infty}^\infty \frac{1}{1+x^2} dx
= \pi
$$
