A convenient choice of Riemann sum in this case is the following: Let $h>0$ and consider the Riemann sum for $(\phi*\psi)(x)$ given by
$$
\tag{1}
R_h(x):=\sum_{k\in\mathbb{Z}^n} \phi(x-hk)\psi(hk)h^n.
$$
As $h\rightarrow 0$, one finds $R_h\rightarrow \phi*\psi$ with respect to most reasonable topologies. I took this parametrization from Lemma 4.1.3 of the classic 'The Analysis of Linear Partial Differential Operators 1' by Hörmander. It looks quite close to what you had in mind!
As for the techniques that go into proving convergence, it is mostly about appropriately utilizing the rapid decay of Schwartz functions. The rapid decay affords us effective estimates to ensure convergence of sums and integrals, and allows us to construct a dominating function to be used in conjunction with the dominated convergence theorem.
First, for fixed $h>0$, to show that the series (1) converges in $C_b(\mathbb{R}^n)$, the space of bounded continuous functions with supremum norm, note that $\lvert\phi(x-hk)\rvert\leqslant \lVert \phi \rVert_{L^\infty(\mathbb{R}^n)}$, while there is a constant $C>0$ such that $\psi(hk)\leqslant C(1+ h^2 \lvert k\rvert^2)^{-n}$. But then
\begin{align*}
\sum_{k\in\mathbb{Z}^n} \sup_{x\in\mathbb{R}^x}\lvert \phi(x-hk)\psi(hk)h^n \rvert & \leqslant C\lVert \phi\rVert_{L^\infty(\mathbb{R}^n)}\sum_{k\in\mathbb{Z}^n} \left(\frac{h}{1+h^2\lvert k \rvert^2}\right)^n \\
& \leqslant C\lVert \phi\rVert_{L^\infty(\mathbb{R}^n)}\left(h+2\sum_{k=1}^\infty \frac{h}{1+ h^2 k^2 }\right)^n \\
& \leqslant C2^n\lVert \phi\rVert_{L^\infty(\mathbb{R}^n)}\left(h + \int_{-\infty} ^{\infty}\frac{dt}{1+ t^2 }\right)^n < \infty,
\end{align*}
which ensures convergence of (1) in $C_b(\mathbb{R}^n)$. Now argue that convergence also happens in $\mathscr{S}(\mathbb{R}^n)$ (still for fixed $h>0$).
Next, we argue that $R_h\rightarrow \phi*\psi$ in $C_b(\mathbb{R}^n)$ as $h\rightarrow 0$. Let
$$
C_{k,h}=\{x\in\mathbb{R}^n \, | \, \forall j\in\{1,\ldots,n\}:\lvert x_j - k_j \rvert \leqslant h/2 \}
$$
be the cube with side length $h>0$ centered at $k\in\mathbb{R}^n$. Let $f^y(x)=f(y-x)$ and let $f_h=\sum_{k\in\mathbb{Z}^n} f(k) \cdot 1_{C_{k,h}}$. Then we see that
$$
R_h(x)=\int (\phi^x)_h \psi_h \, dy.
$$
Noting that
\begin{align*}
R_h(x)-\phi*\psi(x)&=\int (\phi^x)_h \psi_h \, dy - \int \phi^x \psi \, dy\\
&= \int ((\phi^x)_h-\phi^x) \psi_h \, dy + \int \phi^x (\psi_h-\psi) \, dy,
\end{align*}
we have
$$
\lVert R_h - \phi*\psi\rVert_{L^\infty(\mathbb{R}^n)} \leqslant \lVert \phi_h - \phi\rVert_{L^\infty(\mathbb{R}^n)} \lVert \psi_h\rVert_{L^1(\mathbb{R}^n)} + \lVert \phi\rVert_{L^\infty(\mathbb{R}^n)} \lVert \psi_h - \psi\rVert_{L^1(\mathbb{R}^n)}.
$$
Since $\phi$ is uniformly continuous, $\lVert \phi_h - \phi\rVert_{L^\infty(\mathbb{R}^n)}\rightarrow 0$ as $h\rightarrow 0$, so it remains to show that $\psi_h\rightarrow \psi$ in $L^1(\mathbb{R}^n)$. Since $\psi_h\rightarrow \psi$ pointwise, it suffices to find a dominating function $m\in L^1(\mathbb{R}^n)$ with $\lvert\psi_h(x)\rvert\leqslant m(x)$ for, say, $0\leqslant h \leqslant 1$. Consider $m(x)=\max_{k\in C_{x,1}} \lvert \psi(k)\rvert$. Then we clearly have $\lvert \psi_h(x) \rvert \leqslant m(x)$ when $0\leqslant h \leqslant 1$. Furthermore, there is a constant $C'>0$ such that
$$
\lvert\psi(x)\rvert \leqslant C'(1+n+\lvert x\rvert^2)^{-n}.
$$
thus, if $k\in C_{x,1}$, then
$$
\lvert\psi(k)\rvert \leqslant C'(1+n+\lvert k\rvert^2)^{-n} \leqslant C'(1+\lvert x\rvert^2)^{-n},
$$
since $\lvert x- k\rvert^2\leqslant n$. We conclude that
$$
m(x)=\max_{k\in C_{x,1}} \lvert \psi(k)\rvert \leqslant C'(1+\lvert x\rvert^2)^{-n},
$$
and therefore $m\in L^1(\mathbb{R}^n)$.
Finally, finish the argument by showing that $R_h\rightarrow \phi*\psi$ in $\mathscr S(\mathbb{R}^n)$ as well.
Hörmander, Lars, The analysis of linear partial differential operators. I. Distribution theory and Fourier analysis., Grundlehren der Mathematischen Wissenschaften, 256. Berlin etc.: Springer-Verlag. xi, 440 p. DM 69.00/pbk; DM 128.00/hbk (1990). ZBL0712.35001.
