We calculate correlation functions of two local operators within the nucleon carrying momentum. We resolve their dependence on the spatial distance of the currents. This is carried out for all Wick contractions, taking into account several operator insertion types. The resulting four-point functions can be related to parton distribution functions as well as to Mellin moments of double parton distributions. For the latter, we analyze their quark spin and flavor dependency. In this first study, we employ an $N_F = 2 + 1$ CLS ensemble on a $96 \times 32^3$ lattice with lattice spacing $a = 0.0856\ \mathrm{fm}$ and the pseudoscalar masses $m_\pi = 355\ \mathrm{MeV}$ and $m_K = 441\ \mathrm{MeV}$.