Calculate nearest spherical harmonic for a 3d distribution over a grid? I have a set of solutions to an equation, which are all very similar to spherical harmonics. The solutions are discretised on a regular 3d grid. I would like to label them with which spherical harmonic they are most like, so the $l$ and $m$ values for each one.
Theoretically this is just selecting the largest coefficient from a decomposition of each solution vector onto the spherical harmonic space, or a transform of some kind.
Is there a cheap/simple way to calculate which harmonic each solution is nearest? 
 A: You could do a standard linear least-squares fit to a linear combination of spherical harmonics with variable coefficients:
$$f(\vec x)=\sum_{lm}c_{lm}Y_{lm}(\theta,\phi)$$
$$\sum_{\vec x}\left|f(\vec x)-f_{\vec x}\right|^2\rightarrow \min$$
$$\sum_{\vec x}\left|\sum_{lm}c_{lm}Y_{lm}(\theta,\phi)-f_{\vec x}\right|^2\rightarrow \min$$
$$\sum_{\vec x}Y^*_{lm}(\theta,\phi)\left(\sum_{l'm'}c_{l'm'}Y_{l'm'}(\theta,\phi)-f_{\vec x}\right)=0$$
$$\sum_{l'm'}\left(\sum_{\vec x}Y^*_{lm}(\theta,\phi)Y_{l'm'}(\theta,\phi)\right)c_{lm}=\sum_{\vec x}Y^*_{lm}(\theta,\phi)f_{\vec x}\;,$$
where $\theta$ and $\phi$ are the angular variables corresponding to $\vec x$, $f(\vec x)$ is the fitted function and $f_{\vec x}$ are the function values you have on the grid. This is a system of linear equations for the coefficients $c_{lm}$; you could solve it and say that your function is "nearest" to the spherical harmonic whose coefficient has the highest magnitude.
However, by using only function values within a spherical volume, you can consider the inner sum on the left-hand side as an approximation of an integral:
$$
\begin{eqnarray}
\sum_{\vec x}Y^*_{lm}(\theta,\phi)Y_{l'm'}(\theta,\phi)
&\approx&
\int Y^*_{lm}(\theta,\phi)Y_{l'm'}(\theta,\phi)\mathrm dV
\\
&=&
\iint Y^*_{lm}(\theta,\phi)Y_{l'm'}(\theta,\phi)\mathrm d\Omega r^2\mathrm dr
\\
&=&
\delta_{ll'}\delta_{mm'}\int r^2\mathrm dr
\\
&=&
\frac13\delta_{ll'}\delta_{mm'}R^3\;,
\end{eqnarray}
$$
where $R$ is the radius of the spherical volume. This may actually be quite a good approximation if you have e.g. a cubical grid, since the symmetry of the grid causes many of the coefficients approximated as zero to actually vanish.
This radically simplifies the least-squares solution, since the matrix on the left is now the identity and the coefficients are simply given by the sums on the right.
A: It really depends on what your functions are. Since you stated that your functions are very similar to spherical harmonics, I will assume that your functions are almost completely independent of the radius $r$ in spherical coordinates, and therefore it is well approximated by a function $f(\theta,\phi)$ defined over the unit sphere. 
Now, if we also know, a priori, that the function is supposed to be very close to a fixed spherical harmonic, in the sense that we have a priori knowledge that the spherical harmonic decomposition 
$$ f(\theta,\phi) = \sum c_{lm} Y_{lm}(\theta,\phi) $$
is such that there is a distinguished $l_0,m_0$ such that $c_{l_0m_0}$ is really large and $c_{l'm'}$ very small for $(l',m')\neq (l_0m_0)$ (and our goal is to find this pair $(l_0,m_0)$, you can make a guess at what $l_0$ is by considering that
$$ \int_{\mathbb{S}^2}|\nabla f|^2 dvol = -\int_{\mathbb{S}^2} f \triangle f \sim \int_{\mathbb{S}^2} (l_0+1)l_0|f|^2 dvol $$
So by considering 
$$ \tilde{l} = \left\lfloor \sqrt{\frac{\int |\nabla f|^2}{\int f^2}}\right\rfloor $$
you get a good guess what $l_0$ should be. From this it requires to make a finite number of computations (compute the spherical harmonic coefficients for $l \in [\tilde{l}-2, \tilde{l}+2]$ and their accompanying admissible $m$ values, your need to just compute around $11l_0$ coefficients. 
