Scenario
I have a non-homogeneous Poisson point process (PPP) $X\in\mathbb{R}^2$ with intensity function $\lambda(u)$ that is observed over a bounded region $W$. This PPP is modifyed by a dependent thinning with indicator:
$$
I_{x,W} = \begin{cases}
0,\quad \text{if} \quad \exists x'\in X \cap W : \|x-x'\|\le r\ \land \ x\ne x',\\
1, \quad \text{otherwise}
\end{cases}
$$
that is, a point $x\in X$ is deleted if there exists other point $x'\in X$ at a distance less than $r$. In the example of the image below we have $I_{x,W}=0$, and the point $x$ (and $x'$ as well) would be deleted in the thinning.

Calling $Y$ the PPP obtained after applying the $I_{x,W}$ thinning to the PPP $X$, I want to obtain the intensity measure of $Y$ of an area $B \subset W$: $$ \mathbb{E}N(B) $$ where $N(B)$ denotes the number of points of the PPP $Y$ falling in $B$. So this is what I ask for in this question. Below you can see what I've tried without success.
I've tried this without success
For any $B \subset W$ we define: $$ \begin{align*} f:\quad \mathbb{R}^2 \times \Omega & \to \{0,1\}\\ (x, A) & \mapsto \mathbb{1}_B(x)\ \mathbb{1}\left(dist(x, X_W\setminus x) \gt r\right) \end{align*} $$ where $X_W$ are the points lying at region $W$. Then: $$ \nu_Y(B) = \mathbb{E} N(B) = \mathbb{E}\sum_{x \in X_W} f(x,X) \tag 1 $$ gives the intensity measure for $Y$, and using the Campbell-Mecke formula (Theorem 3.2 of [1]) I get: $$ \begin{align} \mathbb{E}\sum_{x \in X_W} f(x,X) =& \int_{\mathbb{R}^2} \mathbb{E}^x[f(x,X)]\nu(dx)\\ =&\int_{\mathbb{R}^2}\mathbb{E}^x[f(x,X)]\lambda(x)\ dx\\ =&\int_{\mathbb{R}^2}\int_\Omega f(x,X) \mathbb{P}^x(A)\ dA\ \lambda(x)\ dx\\ =&\int_{\mathbb{R}^2}\int_\Omega \mathbb{1}_B(x)\ \mathbb{1}\left(dist(x, X_W\setminus x) \gt r\right) \mathbb{P}^x(A)\ dA\ \lambda(x)\ dx\\ =&\int_{\mathbb{R}^2} \mathbb{1}_B(x) \int_\Omega \mathbb{1}\left(dist(x, X_W\setminus x) \gt r\right) \mathbb{P}^x(A)\ dA\ \lambda(x)\ dx\\ =&\int_B \int_\Omega \mathbb{1}\left(dist(x, X_W\setminus x) \gt r\right) \mathbb{P}^x(A)\ dA\ \lambda(x)\ dx\\ =&\int_B \mathbb{P}^x\left(dist(x, X_W\setminus x) \gt r\right)\ \lambda(x) \ dx\\ \tag 2 =& \int_B e^{-\nu(B(x,r))}\ \lambda(x) \ dx \end{align} $$ where $\nu(B(x,r))$ is the intensity measure of the PPP $X$ at the circle of radius $r$ centered at $x$, and hence can be expressed as: $$ \nu(B(x,r)) = \int_{B(x,r)} \lambda(u)\ du $$ So the final expression I get is: $$ \mathbb{E}\sum_{x \in X_W} f(x,X) = \int_B e^{-\int_{B(x,r)}\lambda(u)\ du}\ \lambda(x)\ dx \tag 3 $$
Why it is wrong
I've done some simulations and the results don't even approach to (3). I think it is because of the last step (2). When I derive expression $\mathbb{P}^x\left(dist(x, X_W\setminus x) > r \right)$ I envision it as the empty space function of the PPP $X$, i.e.: $$ F(r)=\mathbb{P}(N(B(x,r))=0) = e^{-\nu(B(x,r))} $$ but don't consider that $X$ is restricted to $W$. That might be the error, although simulations and intuition tell me that I might have more errors in the equalities present in (2).
References:
[1] Stochastic Geometry. (n.d.). Springer Nature., 48-49
Dr. Møller kindly answered my doubt saying the following:
So it seems like the expression I reached was ok.