First of all, all the Associated Legendre Polynomials (ALP) I'm mentioning below are NORMALISED according to the convention of Spherical Harmonics, and the ALPs can be accessed in Mathematica using the command SphericalHamonicY[ $l,m,\theta,0$ ].
Now suppose if I have a series of DISCRETE data points for a continuous function $f(\cos(\theta))$ defined on $[0,\pi]$. If the values on the N-th order Gaussian-Legendre quadrature nodes happen to be known, the coefficients for the $l$-th degree, $m$-th order ALP can be calculated using numerical methods:
$$
C_l^m = \sum w_N(i) f(\cos(\theta_i))P_l^m(\cos(\theta_i))
$$
which I suppose is just equivalent to the continuous integration of
$$
C_l^m = \int f(\cos(\theta))P_l^m(\cos(\theta)) d\theta
$$
if $f$ is known.
so I was thinking, since I have many more data points than the only a few ones at the Gaussian nodes, I can fit it with a continuous polynomial function first and do the integration using the second equation and presumably obtain higher accuracy than the numerical method. But to my surprise, the coefficients are completely different to the Gaussian ones. This should not be caused by any normalisation convention difference as I was using the same function (SphericalHamonicY[ $l,m,\theta,0$ ]).I would be thankful if anyone could advise me on the reason for it, and it is my first question.
My second question is, after obtaining all the ALP coefficients, if I want to reconstruct the function to see how well these coefficients can represent the original function. I'm aware that the SphericalHamonicY function has a different normalisation convention than the ALPs, and the ALPs are not necessarily polynomials or orthogonal, but is there a way to visulise the fitting results to be directly comparable with the data please? Currently when I want to use the obtained data to construct a curve which fits the data above, it looks like this:
Thanks a lot!