[R-sig-Geo] Problem with spTransform en +towgs84 ?

[Roger Bivand]

On Tue, 21 Sep 2010, Agustin Lobo wrote:

> Hi!
>
> Might it be that spTransform is not taking the +towgs84 string into account?
> I have:
>> proj4string(polsfl1_1000)
> [1] " +proj=utm +zone=31 +ellps=intl +units=m +no_defs +towgs84=-87,-98,-121"
>
> and do:
>> polsfl1_1000lonlat <- spTransform(polsfl1_1000, CRS("+proj=longlat +datum=WGS84"))
>> writeOGR(polsfl1_1000,  dsn="/media/TRANSCEND/FLUXPYR/FLIGHT201009/polsfl1_1000",layer="polsfl1_1000", driver="ESRI Shapefile")
>> writeOGR(polsfl1_1000lonlat,  dsn="/media/TRANSCEND/FLUXPYR/FLIGHT201009/polsfl1_1000lonlat.kml",layer="polsfl1_1000", driver="KML")
>
> Then, in QGIS, with reprojection on the fly activated, the 2 vectors
> are not coincident.
> Instead, a polsfl1_1000lonlatv2 made from the original polsfl1_1000
> shapefile in QGIS, is coincident.
>
> As there are 2 different programs here, I'm not positive of my diagnostic. But
> I think there is a problem with spTransform

Could you please provide links to screen shots and input and output data 
objects, as well as sessionInfo() output? Which OS is involved, and are 
spTransform() and QGIS using the same PROJ.4 library and metadata. Is your 
input actually:

CRS("+init=epsg:23031 +towgs84=-84,-107,-120")

with the three-parameter datum shift from:

http://www.asprs.org/resources/GRIDS/

of July 2000? Where does yours come from? The way QGIS treats the *.prj 
file is also of interest.

Roger

>
> Thanks
>
> Agus
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at stat.math.ethz.ch
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>

-- 
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