[Proj] "Double Ellipsoid" error, reproduction

[Frank Warmerdam]

Gerald I. Evenden wrote:
> On Sunday 07 December 2008 10:43:25 pm Richard Greenwood wrote:
>> On Sun, Dec 7, 2008 at 7:59 PM, Frank Warmerdam <warmerdam at pobox.com> wrote:
>>> This is why I believe it was suggested we should be able (potentially)
>>> to use different ellipsoid values for the projection calculations and
>>> for the datum shift.
>> Which seems like it would solve all the problems (at least the
>> technical problems) discussed thus far in this thread. But there is no
>> way to do this in the existing implementation of cs2cs - please
>> correct me if I am wrong. I've not seen software that allows different
>> ellipsoids for the projection and datum. Does such software exist?
>> Frank - is this something that you would contemplate adding to cs2cs?

Rich,

Using a different ellipsoid for the datum shift and the projection logic
is not currently supported by PROJ.4.  I would be willing to contemplate
adding it, though I'm concerned it would be messy and might not really be
used a great deal.

One of the things I had hoped to evaluate was just how often the
ellipsoid for projection and datum shifting purposes are not the same. Also
if it turns out that in essentially all such cases the geographic
coordinates are just supposed to be treated as WGS84 we can just
use the +nadgrids=@null trick which essentially means treat the shift
to WGS84 as being a null shift in geographic (lat/long) space.

Gerald I. Evenden wrote:
> This email and my previous one with the flow chart brings up the question Is 
> it possible to code of a datum shift procedure that merely expects geographic 
> on sysin and spews out datum shifted geographic on sysout.
> 
> With that process, Rich can do anything he wants if he is willing to contruct 
> a pipe procedure and learn to directly use proj4 or lproj.

Yes, this is possible now using pipelining something like:

proj4 -I +proj=merc +ellps=sphere
    | cs2cs +proj=latlong +datum=WGS84 +to +proj=latlong +datum=potsdam
    | proj +proj=... +datum=potsdam

In the above scenario the original projection from mercator would be
done using a sphere, while the datum shift is done assuming a WGS84 ellipsoid.

However, it is not practical for most folks using PROJ.4 in an application
to do the transformation with a pipeline ... instead the API is called
directly with a coordinate system defined using pj_init() for the source
and destination of a pj_transform() call.  In general (for stuff like
MapServer, GDAL, GRASS, QGIS) all the user can specify is the textual
PROJ.4 representation of the source and destination coordinate system.

Best regards,
-- 
---------------------------------------+--------------------------------------
I set the clouds in motion - turn up   | Frank Warmerdam, warmerdam at pobox.com
light and sound - activate the windows | http://pobox.com/~warmerdam
and watch the world go round - Rush    | Geospatial Programmer for Rent