[Proj] CEA from an ellipse and +R_A

[David James Brockley]

Hi Gerald

Thanks for the frank reply. I'm new at map projections!

What I want to do is : call pj_fwd() with a pair of geodetic  
coordinates on the WGS84 ellipsoid, have the latitude transformed to  
authalic latitude, and get XY in cylindrical equal area with a radius  
equal to the authalic radius for WGS84. As far as I can see from the  
books, this is the normal way to convert from WGS84 lon/lat to  CEA x/ 
y. Please correct me if I'm hopelessly confused already by this stage.

As I understand it, the Authalic sphere is an intermediate step  
between a reference ellipsoid (in my case WGS84) and a map projection  
(cylindrical equal area). It allows one to project from the ellipsoid  
and causes the result to still be equal area, but you have to replace  
geodetic latitude with authalic latitude.

 From Map Projections: Theory and Application p102

R^2_A = a^2(1-e^2)( 1 + 2e^2/3 + 3e^4/5 + ... )

and

sin phi_a = sin phi ( ( 1 + 2e^2 sin phi / 3 + ... ) / ( 1 + 2e^2/3  
+ ... )

So, given a reference ellipsoid with known a, proj works out R_A -  
this works and is as I expected.

What I would now expect is for proj to take the input Geodetic  
latitude and convert it internally into an Authalic latitude before  
doing the coordinate transform i.e.

y = R_A * sin( phi_a ) / k

but what actually happens (tracing with GDB) is

y = R_A * sin( phi ) / k

This is when I initialise with +R_A

Is the convention that when using +R_A you have to do the conversion  
to Authalic latitude outside of proj?

I can see from PJ_cea.c that for the spheroidal case, pj_qsfn() and  
pj_authlat() are used to handle the conversion to and from authalic  
latitude, but in this case proj is using WGS84 a for the radius, not  
R_A. (If I initialise with +proj=cea +long_0=0 +ellps=WGS84)

When I initialise with +R_A, I get the radius I want, but the  
projection switches to using the spherical code so the transformation  
to authalic latitude doesn't happen.

Hope this helps explain my confusion. Thanks for your help.

Brock

On 11 Jul 2008, at 17:56, Gerald I. Evenden wrote:

> On Friday 11 July 2008 5:57:49 am David James Brockley wrote:
>> Hello
>>
>> I'm trying to convert geodetic lon/lat into CEA xy using the proj API
>> initialised with
>>
>> +proj=cea +long_0=0 +ellps=WGS84 +R_A
>>
>> What I expect to happen is
>>
>> y = a * q / 2
>>
>> where a is the authalic radius calculated from the semi-major of  
>> WGS84.
>>
>> q/2 seems to be a common approximation to q/q_p which is sin
>> ( authalic latitude)
>>
>>
>> What actually happens is that proj uses the equation for a sphere
>> rather than an ellipse.
>
> That is exactly what is supposed to happen as implied with the  
> capital R in
> option's acronym,  The R in this case is based upon a sphere of  
> surface area
> equal to the ellipse selected.  *All* of rhe R_* option series  
> create radii
> to be used with the spherical earth projection formulas.
>
> Quite frankly, I do not understand what you are trying to do and your
> operations do not make sense to me.  Authalic projections are based  
> upon the
> math of the transformation and not on simple manipulation of the  
> elliptic
> elements.
>
>> y = a sin( phi )
>>
>> I can get proj to do what I expect by initialising with
>>
>> +proj=cea +long_0=0 +a=6371007 +e=0.0818191
>>
>> but that seems somewhat back to front.
>>
>>
>> Can anyone tell me the correct way to configure proj for what I am
>> trying to do?
>>
>> I'm fairly sure that my expected function for y is correct, since
>> I've both copied it from the literature and derived it by hand from
>> the spherical case by substituting authalic radius for radius and
>> authalic latitude for latitude - I've also realised that the /2 is
>> due to the maximum possible value of q_p being 2 when e=0, and that
>> it doesn't affect the properties of the transform (still equal area).
>>
>> Regards, Brock
>> _______________________________________________
>> Proj mailing list
>> Proj at lists.maptools.org
>> http://lists.maptools.org/mailman/listinfo/proj
>
>
>
> -- 
> The whole religious complexion of the modern world is due
> to the absence from Jerusalem of a lunatic asylum.
> -- Havelock Ellis (1859-1939) British psychologist
> _______________________________________________
> Proj mailing list
> Proj at lists.maptools.org
> http://lists.maptools.org/mailman/listinfo/proj