Re: [R-sig-Geo] [OT] Coordinates and pojection problems

[Paul Hiemstra]

Hi,

You need to find out what the projection of the shapefile is, the 
ogrinfo does not give this information for your file as you probably do 
not have the .prj file. The R function spTransform needs the projection 
info in a so called proj4string. The easiest way to define such a string 
is through the epsg code, this code enables proj (the projection library 
used in R) to find the parameter it needs. For example, the proj4string 
for the Dutch national projection is:

"+proj=sterea +lat_0=52.15616055555555 +lon_0=5.38763888888889 
+k=0.9999079 +x_0=155000 +y_0=463000 +ellps=bessel +units=m 
+towgs84=565.04,49.91,465.84,-0.409394387,0.35970519561,-1.868491,4.0772 
+no_defs"

using epsg codes:

"+init=epsg:28992"

The bottomline is to find the projection of your second shapefile as 
this is the most likely candidate to cause the difference. A good place 
to ask questions regarding projections is the proj4-mailing list 
(http://lists.maptools.org/mailman/listinfo/proj).

hth and cheers,
Paul

Julien Barnier wrote:
Hi,

I'm sorry to post a rather off-topic question here, but as a new sp
user and total novice in coordinates and projections questions, I'm
facing a problem you may help me to solve. Please ignore this post if
it is too much noise.

I've got two map files (a MapInfo file and an ESRI shapefile) from two
different sources but which represent objects from the same french
region (one represent the cities, the other rivers and water
surfaces). 

The problem is that when imported into R and plotted together, they
doesn't appear in the same place.  They seem to have the same
«longitude», but the «latitude» is different. If I plot them, the
first map ranges in the Y axis from 70,000 to 100,000 and the second
one from 2,060,000 to 2,140,000. The ranges along the x axis seem to
be the same, though.

Here is the output of ogrinfo on the two files :

,----
| $ ogrinfo -so 69_iris.mid 69_iris
| Had to open data source read-only.
| INFO: Open of `69_iris.mid'
|       using driver `MapInfo File' successful.
| 
| Layer name: 69_iris
| Geometry: Unknown (any)
| Feature Count: 774
| Extent: (747910.000000, 2053380.000000) - (819740.000000,
| 2147622.300000)
| Layer SRS WKT:
| PROJCS[«unnamed»,
|     GEOGCS[«unnamed»,
|         DATUM[«MIF 9999,6,-168,-60,320,0,0,0,0,0»,
|             SPHEROID[«Clarke 1880»,6378249.145,293.465],
|             TOWGS84[-168,-60,320,-0,-0,-0,0]],
|         PRIMEM[«non-Greenwich»,0],
|         UNIT[«degree»,0.0174532925199433]],
|     PROJECTION[«Lambert_Conformal_Conic_2SP»],
|     PARAMETER[«standard_parallel_1»,45.90287723937],
|     PARAMETER[«standard_parallel_2»,47.69712276063],
|     PARAMETER[«latitude_of_origin»,46.8],
|     PARAMETER[«central_meridian»,2.337229104484],
|     PARAMETER[«false_easting»,600000],
|     PARAMETER[«false_northing»,2200000],
|     UNIT[«Meter»,1]]
| DepCom: String (5.0)
| Nom_Com: String (40.0)
| Iris: String (4.0)
| DComIris: String (9.0)
| Nom_Iris: String (43.0)
| Typ_Iris: String (1.0)
| Indic: String (1.0)
| Origine: String (1.0)
`----

,----
| $ ogrinfo -so Decoupage_Administratif/SURFDEAU.shp SURFDEAU
| INFO: Open of `Decoupage_Administratif/SURFDEAU.shp'
|       using driver `ESRI Shapefile' successful.
| 
| Layer name: SURFDEAU
| Geometry: Polygon
| Feature Count: 168
| Extent: (785619.190000, 68133.620000) - (818193.410000, 107443.340000)
| Layer SRS WKT:
| (unknown)
| NUMOBJ: String (11.0)
`----

I've tried to play with some parameters in readShapePoly and
spTransform, but without success.

Any help would be greatly appreciated.

Thanks in advance,
  


--
Drs. Paul Hiemstra
Department of Physical Geography
Faculty of Geosciences
University of Utrecht
Heidelberglaan 2
P.O. Box 80.115
3508 TC Utrecht
Phone:  +31302535773
Fax:    +31302531145
http://intamap.geo.uu.nl/~paul

_______________________________________________
R-sig-Geo mailing list
R-sig-Geo@stat.math.ethz.ch
https://stat.ethz.ch/mailman/listinfo/r-sig-geo