[R-sig-Geo] Transforming from cartographic to arbitrary local coordinate system

[Roger Bivand]

On Mon, 19 Nov 2007, Don MacQueen wrote:

> I have virtually no experience with cartographic coordinate systems,
> and transforming them, so I would like to ask for some help.
>
> Someone has supplied me with a shapefile. When I open it in Qgis,
> Qgis tells me its "spatial reference system" is:
>
>   +proj=lcc +lat_1=35.4666666666667 +lat_2=34.0333333333333
> +lat_0=33.5 +lon_0=-118 +x_0=1999999.9371016 +y_0=499999.9843516001
> +ellps=GRS80 +to_meter=0.3048006 +no_defs
>
> I can read it into R using
>    tmp <- readShapeLines('filename',
>                       proj4string=CRS(pstring))
>
> where pstring is that thing that Qgis gave me.
> plot(tmp) then gives a picture that makes sense.
>
> I would like to give the thing a convenient (to me) local coordinate
> system. I have quite a few such objects (all with the same proj.4
> spatial reference string (or coordinate specification?). Some are
> lines objects, some are polygon objects, and the polygon objects have
> nested polygons ("holes" or "islands").
>
> Pretend for the moment that the lines of my object represent borders
> of a rectangular shaped plot of land . The rectangle is tilted
> relative to north/south.
>
> What I now would like to do is rotate and shift the rectangle so that
> its borders appear vertical and horizontal (my local axes are
> parallel to its borders), and so that it has a local origin near one
> of the corners. I guess this isn't a true cartographic projection (?).

The new elide methods in maptools would get you some of the way there, but 
only do 90 degree rotation, not arbitrary rotation. In general, sp classes 
expect north to be upwards, and no provision for other constructions is 
made. If you can find an appropriate PROJ.4 description for what you want, 
spTransform will do it, but the rotation is going to be the problem.

The easiest way in is almost certainly to customise elide() and adding 
trigenometry to do the arbitrary rotation.

Roger

>
> Is there a way I can do this with spTransform(), by supplying an
> appropriate CRS argument?
>
>
> Thanks
> -Don
>

-- 
Roger Bivand
Economic Geography Section, Department of Economics, Norwegian School of
Economics and Business Administration, Helleveien 30, N-5045 Bergen,
Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43
e-mail: Roger.Bivand at nhh.no