[CIG-MC] Errors in CitcomS Manual Appendix C.6.9 & two comments in Topo_gravity.c

[Joy Hines]

Hello,

There are a couple of deceptive descriptions in the manual and in
Topo_gravity.c that could easily be fixed, which potentially would save
users some time and confusion in the future.

I wrote my own code to visualize the geoid output from CitcomS and found
that my map differed quite a bit from the one produced by the application
included in the visual subdirectory, project_geoid.app.  The file
project_geoid.c was reading the geoid output file in a way that contradicted
the manual appendix, so I decided to dig into the code a bit to make sure
that I was reading the output correctly.  It turns out that I was misled by
the following reinforcing errors.

Appendix C.6.9 states that the third column of the geoid output file
contains the sine coefficients of the total geoid, while the fourth column
lists the cosine coefficients.  The opposite is true.

Following are lines 497 - 500 o f Topo_gravity.c, which incorrectly support
the description in the manual.

    /* sin coeff */

    geoid[0] = (float*)malloc(E->sphere.hindice*sizeof(float));

    /* cos coeff */

    geoid[1] = (float*)malloc(E->sphere.hindice*sizeof(float));

The above float pointers are received by the sphere_expansion function on
line 103 of Sphere_harmonics.c as follows,

    void sphere_expansion(E,TG,sphc,sphs)

and in the following lines of the code are treated opposite the comment
descriptions, i.e.  geoid[0] points to the cosine coefficient and geoid[1]
points to the sine coefficient.  Those values are scaled and stored in the
variables harm_geoid[0] and harm_geoid[1] on lines 543 and 544 of
Topo_gravity.c, respectively, maintaining the same ordering of cosine terms
followed by sine terms.

    harm_geoid[0][p] += con1*cont*geoid[0][p];

    harm_geoid[1][p] += con1*cont*geoid[1][p];

The geoid coefficients are then written out in the order of cosine terms
then sine terms on lines 424-431 of Output.c.

    fprintf(fp1,"%d %d %.4e %.4e %.4e %.4e %.4e %.4e\n",

                        ll, mm,

                        E->sphere.harm_geoid[0][p],

                        E->sphere.harm_geoid[1][p],

                        E->sphere.harm_geoid_from_tpgt[0][p],

                        E->sphere.harm_geoid_from_tpgt[1][p],

                        E->sphere.harm_geoid_from_bncy[0][p],

                        E->sphere.harm_geoid_from_bncy[1][p]);

It would be helpful to future users if the manual were updated to reflect
the correct output order, along with the comments in Topo_gravity.c.
 Another useful addition to the manual would be the formula for normalizing
the Legendre polynomials that is implemented in the function
modified_plgndr_a on lines 51-97 of Sphere_harmonics.c, but that's just a
suggestion and is not immediately necessary.

Please let me know if I have incorrectly interpreted CitcomS, so that I can
appropriately update my own code for geoid coefficient processing.

Thanks,
Joy Hines
-- 
Graduate Student Researcher
Department of Geology
University of California, Davis
2205 Earth & Physical Sciences
Davis, CA 95616
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://geodynamics.org/pipermail/cig-mc/attachments/20100223/86dbe506/attachment.htm