[Proj] Re: Local cartesian to UTM using PROJ.4

[Frode Tennebø]

Greg Troxel <gdt <at> ir.bbn.com> writes:

> Rereading, I think I follow, and it was inappropriate for me to say
> "far more precise" - I apologize for overreacting.

NP. :-)
 
> I think that you mean that both X (east) and Y (in the direction of UTM
> grid north) are both in the plane that's tangent to the ellipsoid, with
> Z being the normal to the ellipsoid. 

Correct.

> So for small X, Y, you could just
> add the UTM values of the origin.  The key difference from UTM is that
> it's truly cartesian from there instead of projected.

This was also my assumption.

> > So, for a trivial example:  Origin is UTM/WGS84 (500000, 6000000, 10) 
position 
> > vector is (1000, 2000, 100) (1000 meters east, 2000 north and 100 up) -> 
> > UTM/WGS84 (501000, 6002000, 110).  Not factoring in scaling.
> 
> I don't think this is the right answer.  By moving out (1000,2000) in
> your system, you move in the plane that's tangent to the ellipsoid at
> the origin.  So your elevation will be more than 100m.  Granted this
> will be small at a few km distance, but consider 1000 km.

Yes, my simple model also compensates for earth curvature (which I did not show 
in the trivial example above) and rotation.  This appears to work well for most 
distances as long as the origin is near the central meridian.  Now, I'd like to 
modify the calculation for scale as this results to 4 meters per 10000 meters.

However, while doing some tests I noticed that as the origin moves away from 
the central meridian, the actual rotation becomes distorted between my simple 
model and a reference program, Geotrans.

I would like to test some inputs with a third party, hopefully proj.4.

> To solve your problem, I think you could take the UTM values, transform
> to ECEF, and then apply translation/rotation, move in your space, and
> then untranslate/rotate and convert back.

Yes, I asume that would give the *correct* values.  However, the work involved 
would be substantial, unless I can find some free Ada routines for the job. :-)

> Surely proj can convert
> geodetic to ECEF, but I don't see it in 'proj -lp'.  These two messages
> may be helpful.
> 
> http://lists.maptools.org/pipermail/proj/2005-February/001515.html
> http://lists.maptools.org/pipermail/proj/2005-February/001519.html

I didn't quite see how......
 
 -Frode