GenericMappingTools / gmt

The Generic Mapping Tools
https://www.generic-mapping-tools.org
Other
847 stars 355 forks source link

Normalization issues with sph2grd #7004

Open chancelorr opened 2 years ago

chancelorr commented 2 years ago

Description of the problem

sph2grd does not appear to distinguish between Schmidt (-Ns) and 4pi normalizations (-Ng). It treats them both as -Ng, which leaves -Ns too large by a factor of $\sqrt{k(2l+1)}$ (Eqn. 10 in Holmes & Featherstone, 2002). I've attached a script which computes both and maps the difference using grdmath. I also have a fortran code which will downward continue an SH model to the CMB, where I implemented the changes that I think fix the problem. The results are all stored in the maps directory.

Inside this zip is everything you would need to produce these results for yourself. The results are pasted in this report though, so you don't need to run anything. sphNormTest.zip

Error No error msg produced by this but I traced the error back to a few sources.

In sph2grd.c: Within the struct N, there is never any explicit use of the attribute N.mode. The only use I found is in specifying bool ortho on L360. As a result, when gmt_plm_bar_all() is called on L428, it does not carry along N.mode, only ortho.

In gmt_stat.c: The function itself (L1098) could be modified to take another boolean input which specifies whether Plm are to be semi-normalized (schmidt) or fully normalized (4pi). For example,

At 1098: void gmt_plm_bar_all (struct GMT_CTRL *GMT, int lmax, double x, bool ortho, bool semi, double *plm)

At 1137: bool csphase = false, semi;

At 1195: r = (2*l+1.0) / (l+m) / (l-m); if (semi) r = 1.0 / (l+m) / (l-m);

Actual outcome incorrectCMB.pdf incorrectSurface.pdf

Expected outcome correctCMB.pdf correctSurface.pdf

The colorscales are the same at each respective radius. On the CMB, you can see the results differ significantly. The effect on the field shape is more noticeable on the surface

System information

welcome[bot] commented 2 years ago

👋 Thanks for opening your first issue here! Please make sure you filled out the template with as much detail as possible. We appreciate that you took the time to contribute!

Please make sure you read our Contributing Guide and abide by our Code of Conduct.