Does a Fourier transformation on a (pseudo-)Riemannian manifold make sense?

There are three main ways of interpreting the Fourier transform.

Decomposition relative to eigenfunctions of the Laplacian

On $\mathbb{R}^n$, the plane waves $E_\xi(x) = \exp( i \xi\cdot x)$ can be interpreted as generalized eigenfunctions of the Laplacian. That is, let

$$ \triangle = \sum_{i=1}^n \left(\frac{\partial}{\partial x_i}\right)^2 $$

then we have

$$ \triangle E_\xi(x) = -|\xi|^2 E_\xi(x) $$

Notice that I say "generalized eigenfunctions". This is because usually it is preferred to do linear algebra/functional analysis on a complete inner product space (a Hilbert space), and in the case of functions on Euclidean space, it is usually preferred to work over the Hilbert space of square integrable functions. The functions $E_\xi(x)$, however, do not belong in the Hilbert space.

But roughly speaking, the Fourier representation of a function

$$ f(x) = \int_{\mathbb{R}^n} \hat{f}(\xi) E_\xi(x) d\xi $$

can be thought of as a decomposition of $f(x)$ in terms of a linear combination of eigenfunctions of the Laplacian.

In this form, as Raskolnikov mentioned in the comment, there is a natural generalization to compact Riemannian manifolds. Let $(M,g)$ be a compact Riemannian manifold, and let $\triangle_g$ denote the associated Laplace-Beltrami operator. One has that $\triangle_g$ is a densely defined, unbounded, self-adjoint operator on $L^2(M)$, the space of square integrable functions on the manifold $M$. As it turns out (using a bit of linear partial differential equations theory and a bit of functional analysis), $\triangle_g$ is invertible if we restrict the co-domain to functions of mean 0 on $M$. The operator $\triangle_g^{-1}$ furthermore turns out to be a compact operator, and thus has discrete eigenvalues with finite dimensional eigen-spaces in $L^2$ by the spectral theorem for compact operators. It is clear by definition that the eigenfunctions of $\triangle_g^{-1}$ will be also eigenfunctions for $\triangle_g$ (just change the eigenvalue $\lambda \to \lambda^{-1}$. Using the spectral theorem we can write an arbitrary square integrable function on $M$ as a unique linear combination of the set of eigenfunctions.

This is one interpretation of the Fourier transform.

You can further refine this theory in two directions. For functions (and to a small extent tensors), instead of explicitly decomposing using an eigenbasis of $L^2$, you can study "smooth frequency cut-offs" by considering the heat-flow on the manifold. That is, given a function $f$, consider the solution to the heat equation $u(t,x)$ defined on $\mathbb{R}_+ \times M$

$$ \partial_t u = \triangle_g u, \quad u(0,x) = f(x) $$

and study, the cut-off of $f$ by a bump-function $\zeta$ in frequencies in the following sense: take a non-negative function $\zeta\in C^\infty_0(\mathbb{R}_+)$ and consider the one-parameter family of functions

$$ f_\lambda(x) = \int_0^\infty \zeta(\lambda t) u(t,x) dt $$

roughly speaking $f_\lambda(x)$ here corresponds to, in the Euclidean space case, the restriction

$$ \int_{\frac12 \lambda <|\xi|< 2\lambda} \hat{f}(\xi) E_\xi(x) d\xi $$

of the function to the spatial frequencies that are close, in absolute value, to $\lambda$. This theory is developed more in Eli Stein's book "Topics in Harmonic Analysis related to Littlewood-Paley Theory".

Another generalization of this decomposition is applicable to not just functions, but also differential forms. In this generalization, instead of looking at the Laplace-Beltrami operator, we consider the Hodge-Laplace operator defined for differential forms. As it turns out, essentially the same argument as for the scalar case give rise to a decomposition of square-integrable differential forms into three categories, the harmonic forms, the exact forms, and the co-exact forms. This goes under the name of Hodge theory, and I refer you to Wikipedia, and also S. Agmon's "Lectures on Elliptic Boundary Value Problems", among many other very good books in this subject.

Fourier transform and representation theory

Another way of treating the Fourier transform on Euclidean space is from the point of view of Euclidean space as an abelian Lie group (in particular a locally compact topological group) which is "self-dual".

Very roughly speaking, for an abelian locally compact topological group $G$, you can consider the set of group homomorphisms from $G\to T$, where $T$ is the one-dimensional circle group. In the case $G = \mathbb{R}^n$, we see that each of the functions $E_\xi(x)$ is one such group homomorphism. This set of group homomorphisms can be made into another topological group, known as the dual group of $G$. In the case of Euclidean space, its dual group has the same structure as another copy of the Euclidean space, hence we say it is self-dual.

Now, given an invariant measure $\mu$ on $G$, in the sense that group translations in $G$ leaves the measure the same (I am being very rough here and sweeping a lot of stuff under the rug), it turns out that we can define the analogue of the Fourier transform by the same formula

$$ \int_G f(x) E_{-\xi}(x) d\mu(x) $$

where $\xi$ is a natural parameter of the dual group of $G$. So in particular, this theory is applicable to the case when we start with a Riemannian manifold $(M,g)$ that is either a Lie group or a symmetric space. This notion is widely applicable in representation theory and analytic number theory, among other things.

See this Wikipedia article, Walter Rudin's "Fourier analysis on groups", S. Helgason's "Differential geometry, Lie groups, and symmetric spaces", and the article by Atle Selberg called "Harmonic analysis and discontinuous groups in weakly symmetric riemannian spaces, with applications to Dirichlet series" for some standard references.

Symmetries on the manifold and partial Fourier transforms

A much more pedestrian way of getting the Fourier transform for functions on manifolds is to just integrate using the usual formula. But as you noted in your question, depending on how you "slice" the manifold, you can get different representations, and whether those representations are meaningful can be questionable.

Fortunately, in one special case we can have a meaningful formula. Let $(M,g)$ be a (pseudo)Riemannian manifold (note that this is the only case I am familiar with that has a natural extension to the case of pseudo-Riemannian manifolds), and assume that it admits a continuous symmetry. That is, there exists a nowhere vanishing vector field $T$ such that the Lie derivative $\mathcal{L}_Tg = 0$; in other words $T$ is a Killing vector. Assume further that the symmetry being considered is a free and proper $\mathbb{R}$ action, and hence that the integral curves of $T$ are diffeomorphic to $\mathbb{R}$ (so in particular $T$ does not have closed orbits. You can also similarly assume the symmetry comes from a free and proper $\mathbb{S}^1$ action, in which case all integral curves are closed, so $T$ does not have open orbits; in this latter case you will deal with Fourier series instead of Fourier transform, I'll leave the obvious changes necessary as an exercise).

Then the action of $T$ defines an equivalence relation on $M$, and we can consider the quotient manifold $Q = M / T$. By the quotient manifold theorem, this quotient is a smooth manifold and also inherits a Riemannian structure, but that is not too important right now.

The point is the following, we can write a point $x\in M$ as $(t,q)\in \mathbb{R}\times Q$, where $t$ is a natural coordinate on $T$ so that $\partial_t$ corresponds to $\mathcal{L}_T$ in $M$. Using this splitting, you can take the Fourier transform in $t$ the usual way, writing $f(x) = f(t,q)$

$$ \hat{f}(\tau,q) = \int f(t,q) \exp( - it\tau) dt $$

That $T$ is a Killing vector is used here so that if you pick a coordinate system adapted to the coordinate $t$, the actual volume form

$$ \sqrt{|g|} dt dq $$

will be independent of $t$, so when you try to do analysis with this notion of the Fourier transform, the formula above will not look like it is missing a factor coming from the volume form (since anything independent of $t$ can be factored out, and taken outside of the integral, so the worst-case scenario is the entire expression is off be a fixed factor depending on $q$, which can always be inserted back later).

This formula can also be easily generalized to the case where $(M,g)$ admits several mutually commuting Killing vector fields $T_1, T_2, \ldots, T_k$. In the case of $\mathbb{R}^n$, the standard coordinate vector fields form a complete set of $n$ mutually commuting Killing vector fields that span the tangent space at every point, and this allows us to take the Fourier transform "in all directions".

This last formulation is useful especially in studying linear partial differential equations on a space-time that is stationary (for example, a stationary black-hole background). See the papers of Dafermos and Rodnianski for example (paper 1, paper 2 ).