Hi,
Forgive me if this should be posted to another list relating to GMT
or the gshhs database. I do not know of any such list at the moment,
and my goal is to use the data in R anyway.
Until now I have made maps with the high-resolution world map that
comes with PBSmapping. As explained in the user guide, these are a
lower resolution version of the gshhs full-resolution database. But
it is not quite good enough for maps of my part of the world. A few
islands are missing that I'd like to have, and the details of the
coastline are not always at a sufficient resolution.
I would have like to make my own coastline database from the original
gshhs database, a version that would only have covered -80 to -20 in
longitude, and 35 to 75 in latitude. I can lower resolution and clip
to smaller areas with functions provided in PBSmapping so I'd be fine
for most maps I have to make.
The problem is that I have next to no experience with compiling
programs and I don't know how to install GMT (which I could use to
extract the coastline) or the smaller programs provided just to
interact with the database.
I thought I did not need to because a shapefile version of the
database exists, and I have managed (with help from Roger in
particular) to import such files in R before. But this one is huge
(154 MB for the shp file alone). The function read.shape refuses to
open it, telling me right way
Erreur dans read.shape("gshhs_land.shp") : unable to open file
I tried with function read.shapefile, which attempts to open the
shapefile but never finishes (I stopped it after 12 h. Thus the
shapefile is either a bit weird (it had more than the usual 3 files
on the web site, I just downloaded the shp, dbf and shx files) or it
is too big, though I have 1 GB of RAM on my Mac.
It appears I'll have to compile and learn to use the programs
provided with that database to extract the data I want to import to R.
Has anyone done this before? Can you point me to a site that explains
how to compile these small c programs on Mac OS X, as a first step?
Are there instructions, somewhere, on how you'd specify the lat long
you want to extract to an ascii file? I did not see any on the GMT site.
Maybe I'm doing all this for nothing and there is already a world
coastline (preferably with main lakes and rivers too) database at
high resolution available for R, beside the ones that come with maps
(which does not allow me enough flexibility) and PBSmapping (not
quite high enough resolution)?
Sincerely,
Denis Chabot
Accessing gshhs database of world coastlines
5 messages · Denis Chabot, Roger Bivand
On Sun, 16 Oct 2005, Denis Chabot wrote:
Hi,
Forgive me if this should be posted to another list relating to GMT
or the gshhs database. I do not know of any such list at the moment,
and my goal is to use the data in R anyway.
Until now I have made maps with the high-resolution world map that
comes with PBSmapping. As explained in the user guide, these are a
lower resolution version of the gshhs full-resolution database. But
it is not quite good enough for maps of my part of the world. A few
islands are missing that I'd like to have, and the details of the
coastline are not always at a sufficient resolution.
I would have like to make my own coastline database from the original
gshhs database, a version that would only have covered -80 to -20 in
longitude, and 35 to 75 in latitude. I can lower resolution and clip
to smaller areas with functions provided in PBSmapping so I'd be fine
for most maps I have to make.
The problem is that I have next to no experience with compiling
programs and I don't know how to install GMT (which I could use to
extract the coastline) or the smaller programs provided just to
interact with the database.
I thought I did not need to because a shapefile version of the
database exists, and I have managed (with help from Roger in
particular) to import such files in R before. But this one is huge
(154 MB for the shp file alone). The function read.shape refuses to
open it, telling me right way
Erreur dans read.shape("gshhs_land.shp") : unable to open file
This means read.shape() cannot see the file. I managed to read gshhs_land.shp into R with current maptools:
try1 <- read.shape("gshhs_land.shp")
Shapefile type: PolyLine, (3), # of Shapes: 191811
object.size(try1)
[1] 395387692 but it is a shapefile of lines, not polygons - the base GSHHS is now polygonal. I have 1GB under linux, so I think memory management is effectively the same as yours. I've put a PDF (33MB) of your area on: http://spatial.nhh.no/misc/nw_atlantic.pdf please save link as, and open locally, otherwise you'll be sitting looking at a blank screen for a long time. I did this so that you can use the Acroread zoom to find out if your islands are there. read.shape() was not too slow, it seemed to go quite well. Is there any use in writing an R function on top of the gshhs program for extracting polygons from the gshhs databases directly into sp class objects - it doesn't look too difficult, and would give direct access without having to use the shapefile? Roger
I tried with function read.shapefile, which attempts to open the shapefile but never finishes (I stopped it after 12 h. Thus the shapefile is either a bit weird (it had more than the usual 3 files on the web site, I just downloaded the shp, dbf and shx files) or it is too big, though I have 1 GB of RAM on my Mac. It appears I'll have to compile and learn to use the programs provided with that database to extract the data I want to import to R. Has anyone done this before? Can you point me to a site that explains how to compile these small c programs on Mac OS X, as a first step? Are there instructions, somewhere, on how you'd specify the lat long you want to extract to an ascii file? I did not see any on the GMT site. Maybe I'm doing all this for nothing and there is already a world coastline (preferably with main lakes and rivers too) database at high resolution available for R, beside the ones that come with maps (which does not allow me enough flexibility) and PBSmapping (not quite high enough resolution)? Sincerely, Denis Chabot
_______________________________________________ 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
Hi Roger, Le 05-10-16 ? 14:09, Roger Bivand a ?crit :
On Sun, 16 Oct 2005, Denis Chabot wrote:
...
I thought I did not need to because a shapefile version of the
database exists, and I have managed (with help from Roger in
particular) to import such files in R before. But this one is huge
(154 MB for the shp file alone). The function read.shape refuses to
open it, telling me right way
Erreur dans read.shape("gshhs_land.shp") : unable to open file
This means read.shape() cannot see the file.
Maybe my copy is correupted? I tried again giving the full pathname
to the file to insure I had not screwed up with defining my working
directory, and got yet the same (or close to):
> test <- read.shape("Macintosh
HD:Users:dchabot:Documents:Cartographie:Shapefiles:gshhs:gshhs_land.shp"
)
Erreur dans read.shape("Macintosh
HD:Users:dchabot:Documents:Cartographie:Shapefiles:gshhs:gshhs_land.shp"
) :
unable to open file
(I had used a utility to give me the full pathname of the shp file).
Alternative explanation: I just upgraded to R2.2, and wiped out my
2.1 install. Because R2.2 is not officially on line at CRAN yet, I
could not download mac binaries for maptools and sp and compiled them
myself. I did get a warning about some name being used twice, and
although the package seems to work (it certainly imports a smaller
shapefile) maybe it is not 100% functional.
I managed to read gshhs_land.shp into R with current maptools:
try1 <- read.shape("gshhs_land.shp")
Shapefile type: PolyLine, (3), # of Shapes: 191811
object.size(try1)
[1] 395387692 but it is a shapefile of lines, not polygons - the base GSHHS is now polygonal.
I don't understand this part. the "gshhs_land" is not "base GSHHS"? I thought I downloaded the most recent version, at the highest resolution. Where did you find the version with polygons?
I have 1GB under linux, so I think memory management is effectively the same as yours. I've put a PDF (33MB) of your area on: http://spatial.nhh.no/misc/nw_atlantic.pdf please save link as, and open locally, otherwise you'll be sitting looking at a blank screen for a long time. I did this so that you can use the Acroread zoom to find out if your islands are there.
Yes the islands are there! This coastline is an improvement on what I have!
read.shape() was not too slow, it seemed to go quite well. Is there any use in writing an R function on top of the gshhs program for extracting polygons from the gshhs databases directly into sp class objects - it doesn't look too difficult, and would give direct access without having to use the shapefile?
I would love to have it. At a minimum, I'd need clear(er) instructions as to how to compile and use the programs provided on GSHHS's site. But until now I valued being able to work with my dataset in R as I am used to, and produce maps for the right subset of data I'm interested in. So I don't plan to use GMT any time soon and would much prefer R functions to access the very nice gshhs databases and take what I need out of them. As an aside, to you know of a way to do map borders in R that look like those in GMT and many paper maps (thick borders, alternating from white to black from one degree of lat or long to the next)? This is a cosmetic issue only but it would be a nice touch to be able to do this in R. Anyway, I'll wait a few days for the full R2.2 Mac to be at CRAN and get fresh binaries of maptools and sp to see if I can open that shapefile, since you were able to. It irritates me to have failed. With many thanks, Denis
On Sun, 16 Oct 2005, Denis Chabot wrote:
Hi Roger, Le 05-10-16 ? 14:09, Roger Bivand a ?crit :
On Sun, 16 Oct 2005, Denis Chabot wrote:
...
I thought I did not need to because a shapefile version of the
database exists, and I have managed (with help from Roger in
particular) to import such files in R before. But this one is huge
(154 MB for the shp file alone). The function read.shape refuses to
open it, telling me right way
Erreur dans read.shape("gshhs_land.shp") : unable to open file
This means read.shape() cannot see the file.
Maybe my copy is correupted? I tried again giving the full pathname to the file to insure I had not screwed up with defining my working directory, and got yet the same (or close to):
> test <- read.shape("Macintosh
HD:Users:dchabot:Documents:Cartographie:Shapefiles:gshhs:gshhs_land.shp"
)
Erreur dans read.shape("Macintosh
HD:Users:dchabot:Documents:Cartographie:Shapefiles:gshhs:gshhs_land.shp"
) :
unable to open file
(I had used a utility to give me the full pathname of the shp file).
I have no access to OSX, I'm afraid. If you do setwd() to the directory
where gshhs_land.shp is placed, and then list.files(pattern="shp"), is the
file visible?`What does getwd() say when you are there? The error message
is from:
hSHP = SHPOpen(CHAR(STRING_ELT(shpnm,0)), "rb" );
if( hSHP == NULL )
error("unable to open file");
in src/Rshapeget.c, where shpnm is the file name given to read.shape(). So
either the compilation failed, or the file isn't there. Try
getinfo.shape("gshhs_land.shp")
the error message there should say "Unable to open: gshhs_land.shp".
Alternative explanation: I just upgraded to R2.2, and wiped out my 2.1 install. Because R2.2 is not officially on line at CRAN yet, I could not download mac binaries for maptools and sp and compiled them myself. I did get a warning about some name being used twice, and although the package seems to work (it certainly imports a smaller shapefile) maybe it is not 100% functional.
I managed to read gshhs_land.shp into R with current maptools:
try1 <- read.shape("gshhs_land.shp")
Shapefile type: PolyLine, (3), # of Shapes: 191811
object.size(try1)
[1] 395387692 but it is a shapefile of lines, not polygons - the base GSHHS is now polygonal.
I don't understand this part. the "gshhs_land" is not "base GSHHS"? I thought I downloaded the most recent version, at the highest resolution. Where did you find the version with polygons?
It is not a shapefile, but in GSHHS format. It is formatted as polygons. The problem with the shapefile as lines is that it cannot be used for filling land or sea, just for drawing coastlines.
I have 1GB under linux, so I think memory management is effectively the same as yours. I've put a PDF (33MB) of your area on: http://spatial.nhh.no/misc/nw_atlantic.pdf please save link as, and open locally, otherwise you'll be sitting looking at a blank screen for a long time. I did this so that you can use the Acroread zoom to find out if your islands are there.
Yes the islands are there! This coastline is an improvement on what I have!
OK, so I'll follow this route as well.
read.shape() was not too slow, it seemed to go quite well. Is there any use in writing an R function on top of the gshhs program for extracting polygons from the gshhs databases directly into sp class objects - it doesn't look too difficult, and would give direct access without having to use the shapefile?
I would love to have it. At a minimum, I'd need clear(er) instructions as to how to compile and use the programs provided on GSHHS's site. But until now I valued being able to work with my dataset in R as I am used to, and produce maps for the right subset of data I'm interested in. So I don't plan to use GMT any time soon and would much prefer R functions to access the very nice gshhs databases and take what I need out of them.
I think it can be done.
As an aside, to you know of a way to do map borders in R that look like those in GMT and many paper maps (thick borders, alternating from white to black from one degree of lat or long to the next)? This is a cosmetic issue only but it would be a nice touch to be able to do this in R.
For unprojected it could be done, for projected like GMT is more trouble. If anyone would like to contribute GMTaxis() to work like GMT, please say, contributions welcome!
Anyway, I'll wait a few days for the full R2.2 Mac to be at CRAN and get fresh binaries of maptools and sp to see if I can open that shapefile, since you were able to. It irritates me to have failed.
Roger
With many thanks, Denis
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
Hi Roger, Le 05-10-17 ? 02:21, Roger Bivand a ?crit :
I have no access to OSX, I'm afraid. If you do setwd() to the directory where gshhs_land.shp is placed, and then list.files(pattern="shp"), is the file visible?`
Yes, I got > list.files(pattern="shp") [1] "gshhs_land.shp"
What does getwd() say when you are there? The error message is from:
> getwd() [1] "/Users/dchabot/Documents/Cartographie/Shapefiles/gshhs"
hSHP = SHPOpen(CHAR(STRING_ELT(shpnm,0)), "rb" );
if( hSHP == NULL )
error("unable to open file");
in src/Rshapeget.c, where shpnm is the file name given to read.shape
(). So
either the compilation failed, or the file isn't there. Try
getinfo.shape("gshhs_land.shp")
> getinfo.shape("gshhs_land.shp")
Erreur dans getinfo.shape("gshhs_land.shp") :
No such file
Wow, these were R commands I did not know! But I just found the
problem. OS X can "hide" the extension of file names. The names of
the 3 files appeared correct to me, but I fired Terminal and listed
the files that way (which I admit I don't do often), and here what I
got:
LAUchabotd-QC-DFO-MPO-GC-CA:~/Documents/Cartographie/Shapefiles/gshhs
dchabot$ ls -l
total 336928
-rw-r--r-- 1 dchabot dchabot 8823500 Oct 11 13:07
gshhs_land.dbf.txt
-rw-r--r-- 1 dchabot dchabot 162143644 Oct 11 13:07 gshhs_land.shp
-rw-r--r-- 1 dchabot dchabot 1534588 Oct 11 13:07
gshhs_land.shx.txt
2 of the files had a hidden ".txt" which I did not see, but R sure
did. I fixed these 2 file names and read.shape gave me the same
number of shapes you got, and it did not take so long. Sorry for this.
... It is not a shapefile, but in GSHHS format. It is formatted as polygons. The problem with the shapefile as lines is that it cannot be used for filling land or sea, just for drawing coastlines.
So you accessed it using gshhs programs? You mentioned you could work on a way to adapt gshhs programs so they'd be called from R, or bringing the same functionality in R, I'm not sure exactly what you meant. But you can save you the trouble for one program, "gshhs_dp". I think thinPolys, from PBSmapping, implements this already.
OK, so I'll follow this route as well.
Is there any use in writing an R function on top of the gshhs program for extracting polygons from the gshhs databases directly into sp class objects - it doesn't look too difficult, and would give direct access without having to use the shapefile?
I would love to have it. At a minimum, I'd need clear(er) instructions as to how to compile and use the programs provided on GSHHS's site. But until now I valued being able to work with my dataset in R as I am used to, and produce maps for the right subset of data I'm interested in. So I don't plan to use GMT any time soon and would much prefer R functions to access the very nice gshhs databases and take what I need out of them.
I think it can be done.
As an aside, to you know of a way to do map borders in R that look like those in GMT and many paper maps (thick borders, alternating from white to black from one degree of lat or long to the next)? This is a cosmetic issue only but it would be a nice touch to be able to do this in R.
For unprojected it could be done, for projected like GMT is more trouble. If anyone would like to contribute GMTaxis() to work like GMT, please say, contributions welcome!
I wish I could help, but I really do not have the background to do so. I will be happy to test for you though. Denis