[GRASS-user] Correctly Converting LL to LCC

[Glynn Clements]

Rich Shepard wrote:

> > 	242807.02 / 796758.58 = 0.3047435272049408
> > 	397955.05 / 1305792.74 = 0.3047612670905185
> >
> > These look suspiciously close to 0.3048 (foot/metre conversion).
> 
>    You and Hamish point out that the problem is with conversion of the units:
> feet vs. meters. What I'm trying to understand is why there is this problem
> and what to do about it.
> 
>    The imported Long/Lat location is PROJ_UNITS of degrees, which is what I
> would expect. The location with the other maps has PROJ_UNITS of meters.
> Therefore, when I run v.proj to convert the ll map coordinates in degrees to
> lcc coordinates in meters, why do I end up with feet instead?

Do you end up with feet? Or do you end up with metres but you're
comparing against a reference given in feet? Your previous comment
suggests the latter:

>    The middle point should have values of 796758.58 (instead of 242807.02)
> and 130579274 (instead of 397955.05).

Your numbers seem roughly correct for metres (45d19'N is ~3.567 degrees
north of lat_0, which equates to ~397km), while the "desired" numbers
are roughly correct for feet (397000 / 0.3048 ~= 1,302,493 feet).

> I see no length unit option in either v.proj or cs2cs.

v.proj uses the information in the PROJ_INFO and PROJ_UNITS files. 
Whether a location uses metres or feet is a property of the location. 
If you use metres when you create the location, all coordinates will be
given in (or assumed to be in) metres.

cs2cs and proj are units-agnostic. If you give the semi-major axis
(+a=... option) in metres, you get coordinates in metres (or, for
cartographic->lat/lon, the input is treated as metres). If you give
the semi-major axis in furlongs, you get coordinates in furlongs.

All of the built-in ellipsoids which can be selected by the +ellps=... 
option use metres. If you want proj or cs2cs to provide or receive
coordinates in units other than those used when specifying the
semi-major axis, you need to use the +to_meter=... option to specify
the scale factor, e.g. +to_meter=0.3048 if you want units of 0.3048
metres (i.e. feet).

>    Hamish suggested I scale the units using v.transform. I tried this with
> both xscale and yscale at 0.3048 as well as 1.0/0.3048 but neither resulted
> in units within the appropriate range.

What are you comparing against? How far off are the numbers?

>    It may be that all my maps from the state GIS data repository are
> incorrect because soils data from a different agency also cannot be
> overlaid. Since I used the ArcView .prj and Arc/Info .e00 files for the
> projection and boundary information I need to learn what I've done
> incorrectly so I can re-import all these data and be able to move ahead.
> 
>    If someone will work with me off the list to get straightened out I would
> greatly appreciate it. I have a 2-page .pdf showing the Oregon data library
> formats and projection information. Since these data are in International
> feet, I need to understand why the PROJ_UNITS are in meters.

This is determined when you create the location.

-- 
Glynn Clements <glynn at gclements.plus.com>