On Fri, 24 Sep 2010, [email protected] wrote:
Hi again,
However,in my case, if instead of the previously mentioned version in my
linux box, I use my windows version of R:
library(rgdal)
Loading required package: sp
Geospatial Data Abstraction Library extensions to R successfully loaded
Loaded GDAL runtime: GDAL 1.6.3, released 2009/11/19
Path to GDAL shared files: C:/ARCHIV~1/R/R-29~1.2/library/rgdal/gdal
Loaded PROJ.4 runtime: Rel. 4.6.1, 21 August 2008
Path to PROJ.4 shared files: C:/ARCHIV~1/R/R-29~1.2/library/rgdal/proj
The mismatch problem is solved. So, could this be, in this case, a problem
of the PROJ.4 Rel. 4.5.0, 22 Oct 2006, or GDAL 1.4.2.0 versions in the
linux computer?
Older releases of the EPSG database for PROJ.4 and GDAL did not
necessarily include +towgs84= parameters for known +datum= - this has
improved with time. The burden is, though, on users to understand that
both the datum and the ellipsoid must be known, and in the case of
ambiguous "datum" definitions like ED50, what the +towgs84= should be.
See projInfo("datum") for known datum definitions.
Currently, the GDAL version is 1.7.2, and PROJ.4 is 4.7.0 (but declares
itself as 4.7.1). These are included in CRAN binaries for Windows and CRAN
extras binaries for OSX. If you install from source under Linux, these are
also what you should have, and Debian/Fedora packagers for GIS are usually
not far behind, and checking for updates in some fashion should be
possible. In R, do make a habit of running sessionInfo() and recording the
output if you need absolutely reproducible research; then
update.packages() to pick up bug fixes and enhancements. If you need to
step back, all older versions of packages are also available on CRAN.
It is known that datum specifications are a problem, which is why one
needs to be particularly careful about correct +towgs84= parameters.
Hope this helps,
Roger
Thanks,
Javier
---
Hi,
Just to add that I've also detected some problems in this sense. I've have
exported a shapefile from ArcGIS to KML, and also imported it into R, and
from there to KML. The final KML files differ in some 70m of latitude,
while the longitudes match.
# import the shapefile
knapdale.supgeo <- readOGR(dsn="test",layer="knapdale_supgeo")
proj4string(knapdale.supgeo)
[1] " +proj=tmerc +lat_0=49 +lon_0=-2 +k=0.999601 +x_0=400000 +y_0=-100000
+ellps=airy +units=m +no_defs"
# Note, the above is OSGB 1936 / British National Grid, which equals
# CRS("+init=epsg:27700")
# Then:
knapdale.supgeoLL <- spTransform(knapdale.supgeo, CRS("+proj=longlat +
datum=WGS84"))
kmlf <- "knapdale_supgeo.kml"
writeOGR(knapdale.supgeoLL, dsn=kmlf, layer="supgeo", driver="KML")
I'm not saying that the R part is responsible for the difference, but the
mismatch is there. Some details:
ArcGIS 9.3 (Windows XP Pro)
R 2.10.1
library(rgdal)
Loading required package: sp
Geospatial Data Abstraction Library extensions to R successfully loaded
Loaded GDAL runtime: GDAL 1.4.2.0, released 2007/06/27
Path to GDAL shared files: /usr/local/share/gdal
Loaded PROJ.4 runtime: Rel. 4.5.0, 22 Oct 2006
Path to PROJ.4 shared files: (autodetected)
Best regards,
Javier
Roger,
I apologize for not including version details etc, here they are:
Loading required package: rgdal
Geospatial Data Abstraction Library extensions to R successfully loaded
Loaded GDAL runtime: GDAL 1.7.2, released 2010/04/23
Path to GDAL shared files: /usr/share/gdal17
Loaded PROJ.4 runtime: Rel. 4.7.1, 23 September 2009
Path to PROJ.4 shared files: (autodetected)
version
_
platform i486-pc-linux-gnu
arch i486
os linux-gnu
system i486, linux-gnu
status
major 2
minor 11.1
year 2010
month 05
day 31
svn rev 52157
language R
version.string R version 2.11.1 (2010-05-31)
Also, I've put SPDF polsfl1_1000 here:
https://sites.google.com/site/openfiles2/home/polsfl1_1000.rda?attredirects=0&d=1
Finally, the towgs84 parameters come from the last paragraph in the pdf
for The Kingdom of Spain (July) listed in the page you mention. (BTW,
I recall having seen
another list with slightly different parameters but cannot find it
anymore)
I had to put
the same 3 parameters for the definition of the projection in QGIS
(which qgis does not
do by default, had to define a custom projection to get the datum
shift working). The towgs
parameters are not in the proj file of the shapefile set, at least not
in the ones created
by writeOGR.
Hope this helps to clarify the problem. I'm still not certain that
spTransform is doing anything wrong.
Thanks for your help,
Agus
2010/9/21 Roger Bivand <[email protected]>:
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
[email protected]
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: [email protected]
_______________________________________________
R-sig-Geo mailing list
[email protected]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
_______________________________________________
R-sig-Geo mailing list
[email protected]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
_______________________________________________
R-sig-Geo mailing list
[email protected]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
_______________________________________________
R-sig-Geo mailing list
[email protected]
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: [email protected]
_______________________________________________
R-sig-Geo mailing list
[email protected]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
