See Michael Lacey's [exposition on Carleson's theorem][1] for some discussion of this, particularly Section 9 there.

In higher dimensions one needs to specify the order of summation. If one sums over the dilation of a fixed polygon (say increasing cubes), then pointwise convergence holds and follows from the one dimensional case. This was initially [proven by C. Fefferman in 1971][2]. If one sums over, say, rectangles of increasing size but with a varying ratio of side lengths, then almost everywhere convergence can fail, as also [observed by C. Fefferman in 1971][3]. 

Arguably the "big open question" here is what happens for spherical summation. This is open even in two dimensions. The problem there seems to require inputs from both the wave packet analysis that was used to prove Carleson's theorem, but also the restriction theory of circles (and, in higher dimensions, spheres). If the function is assumed to be radial, then things simplify somewhat and it can be treated using the one dimensional technology (albeit with additional work and ideas) as done by [Prestini in 1988][4].

As discussed in Lacey's paper, it seems obtaining the analog of the tree lemma used in the Fefferman-Lacey-Thiele-type arguments is a significant obstruction. Very loosely, this has to do with the fact that we understand what taking dyadic projections of a Fourier series does to a one dimensional function. In higher dimensions, these dyadic projections are now taken by restricting the Fourier series to annuli (compared to intervals) which is less understood.


  [1]: https://arxiv.org/abs/math/0307008
  [2]: https://projecteuclid.org/journals/bulletin-of-the-american-mathematical-society-new-series/volume-77/issue-5/On-the-convergence-of-multiple-Fourier-series/bams/1183533033.full
  [3]: https://projecteuclid.org/journals/bulletin-of-the-american-mathematical-society-new-series/volume-77/issue-2/On-the-divergence-of-multiple-Fourier-series/bams/1183532734.full
  [4]: https://link.springer.com/article/10.1007/BF01636929