Let me introduce the notation
$$ \langle \psi,\phi\rangle_\gamma = \int_{\mathbb{R}^2} \bar{\psi}(y) \gamma(x,y) \phi(x) ~\mathrm{d}x~\mathrm{d}y $$
so that
$$ \langle \phi,\psi\rangle_\gamma = \overline{\langle \psi,\phi\rangle_\gamma}$$
Firstly, your definition of $S$ I think is not strictly what they have. Your $S$ should be $$S = S_x\cdot S_y = (1 - \partial_x^2)^\frac12 (1-\partial_y^2)^\frac12$$
This has the advantage that
$$ \langle \psi,\phi\rangle_{S\gamma} = \langle \sqrt{1 - \partial^2} \psi,\sqrt{1 - \partial^2} \phi\rangle_{\gamma} $$
this comes from formally integrating by parts or using that $\sqrt{1-\partial^2}$ is self-adjoint on suitable spaces.
Now, fix $(x,y)\in \mathbb{R}^2$. If $x = y$ then the inequality you want to prove is trivial. So we look at $x \neq y$. Then we can choose disjoint neighbourhoods $U_x, U_y$. Let $f_x$ and $f_y$ be smooth, real valued functions supported in $U_x$ and $U_y$ respectively. Assume further that they are non-negative and not identically zero.
Define $f^{\lambda,\theta} = \frac1\lambda f_x + e^{i\theta} \lambda f_y$, where $\lambda\in \mathbb{R}_+$ and $\theta\in [0,2\pi]$.
Consider now
$$ 0 < \langle \sqrt{1-\partial^2}f^{\lambda,\theta},\sqrt{1-\partial^2}f^{\lambda,\theta}\rangle = \langle f^{\lambda,\theta},f^{\lambda,\theta}\rangle_{S\gamma} $$
This implies
$$ 0 < \frac1{\lambda^2} \langle f_x,f_x\rangle_{S\gamma} + \lambda^2 \langle f_y,f_y\rangle_{S\gamma} + 2\Re \langle f_x, e^{i\theta} f_y\rangle_{S\gamma} $$
Choose $$\lambda^2 = \sqrt{ \frac{\langle f_x,f_x\rangle_{S\gamma}}{\langle f_y,f_y\rangle_{S\gamma}}}$$
and $\theta$ such that $\langle f_x,e^{i\theta}f_y\rangle_{S\gamma}$ is purely real. Then we have
$$ \left| \langle f_x,f_y\rangle_{S\gamma}\right|^2 \leq \langle f_x,f_x\rangle_{S\gamma} \langle f_y,f_y\rangle_{S\gamma} \tag{1}$$
Observe (1) holds for every $f_x\in C^\infty_0(U_x;\overline{\mathbb{R}_+})$ and $f_y\in C^\infty_0(U_y;\overline{\mathbb{R}_+})$. Now you just need to take approximations of identities to get the desired pointwise estimate from the integrated estimate (1).
Remark: the two parameter choice $\lambda,\theta$ is sometimes called the "amplification trick" which Terry Tao explained wonderfully on his blog.