Subscribe to R-sig-geo feed
This is an archive for R-sig-geo, not a forum. It is not possible to post through Nabble - you may not start a new thread nor follow up an existing thread. If you wish to post, but are not subscribed to the list through the list homepage, subscribe first (through the list homepage, not via Nabble) and post from your subscribed email address. Until 2015-06-20, subscribers could post through Nabble, but policy has been changed as too many non-subscribers misunderstood the interface.
Updated: 1 hour 21 min ago

Re: drawing a polygon from several lines

Fri, 08/24/2018 - 09:30
It's my pleasure to share Lulla. Thanks for the suggestion, I will adopt it.

Regards,

Antonio

Em sex, 24 de ago de 2018 às 10:58, Vijay Lulla <[hidden email]>
escreveu:

> Very nice, Antonio!!  Thanks for sharing your script.  The only minor
> suggestion I have to offer is that you can make your "cut isobath" portion
> easier to read (at least IMO) by using `with` and `between` (either
> dplyr::between or data.table::between).  Below is how I would've done it.
>
> limmin <- with(list(coords = coordinates(isomin)[[1]][[1]]),
>                subset(coords, between(coords[, 1], lonmin, lonmax) &
>                               between(coords[, 2], latmin, latmax)))
>
> Again, thanks for sharing your problem and solution.  I found it
> instructive!
> Cordially,
> Vijay.
>
>
> On Fri, Aug 24, 2018 at 9:16 AM Antonio Silva <[hidden email]>
> wrote:
>
>> Dear colleagues
>>
>> I finally could select the squares (cells) under the polygon. I extract
>> the
>> points from the selected isobaths that fell within lat / lon limits and
>> with them I created a spatial polygon.
>>
>> Bellow I share the script I wrote. Considering that I don't know much
>> about
>> mapping in R, I would appreciate to hear suggestions to improve the code.
>>
>> Best regards.
>>
>> --
>> Antônio Olinto Ávila da Silva
>> Instituto de Pesca (Fisheries Institute)
>> São Paulo, Brasil
>>
>> ===================
>>
>> library(rgdal)
>> library(rgeos)
>> rm(list = ls())
>>
>> # import shapes
>> # isobath
>> https://www.dropbox.com/sh/us859uyijns36r9/AAC9N4OSwweLuW-fdsBWUy-ra?dl=0
>> isobaths <- readOGR(".","isobath")
>> isobaths <- spTransform(isobaths, CRS("+proj=longlat +datum=WGS84"))
>> proj4string(isobaths)
>> summary(isobaths)
>>
>> # create a spatial grid and polygons
>> grd <-
>>
>> GridTopology(cellcentre.offset=c(-47.75,-25.416667),cellsize=c(10/60,10/60),cells.dim=c(23,12))
>> squares <- as.SpatialPolygons.GridTopology(grd)
>> proj4string(squares) <- CRS("+proj=longlat +datum=WGS84")
>> summary(squares)
>> IDs <- sapply(slot(squares, "polygons"), function(x) slot(x, "ID"))
>>
>> # set limits
>> latmin <- -24.8
>> latmax <- -24.02
>> lonmin <- -47.2
>> lonmax <- -44.8
>> depmin <- -25
>> depmax <- -60
>>
>> # bounding box
>> mat.area <-
>>
>> matrix(c(lonmin,latmax,lonmax,latmax,lonmax,latmin,lonmin,latmin),ncol=2,byrow=T)
>> spp.area <- SpatialPolygons(list(Polygons(list(Polygon(mat.area)) ,
>> ID='1')))
>> proj4string(spp.area) <- CRS("+proj=longlat +datum=WGS84")
>>
>> # plot polygons isolines and bounding box
>> plot(squares, axes=T)
>> #~ text(coordinates(squares),labels=IDs, cex=0.7)
>> plot(isobaths,add=T)
>> plot(spp.area,border="red",add=T)
>>
>> # select isobaths
>> summary(isobaths)
>> isomin <- subset(isobaths, ID %in% depmin)
>> proj4string(isomin) <- CRS("+proj=longlat +datum=WGS84")
>> isomax <- subset(isobaths, ID %in% depmax)
>> proj4string(isomax) <- CRS("+proj=longlat +datum=WGS84")
>>
>> plot(isomin,add=T,col="deepskyblue3",lwd=2)
>> plot(isomax,add=T,col="dodgerblue4",lwd=2)
>>
>> # cut isobaths
>> limmin <- subset(coordinates(isomin)[[1]][[1]],
>>           coordinates(isomin)[[1]][[1]][,2]<=latmax &
>> coordinates(isomin)[[1]][[1]][,2]>=latmin &
>>           coordinates(isomin)[[1]][[1]][,1]<=lonmax &
>> coordinates(isomin)[[1]][[1]][,1]>=lonmin)
>> limmax <- subset(coordinates(isomax)[[1]][[1]],
>>           coordinates(isomax)[[1]][[1]][,2]<=latmax &
>> coordinates(isomax)[[1]][[1]][,2]>=latmin &
>>           coordinates(isomax)[[1]][[1]][,1]<=lonmax &
>> coordinates(isomax)[[1]][[1]][,1]>=lonmin)
>>
>> points(limmin,col="chartreuse3")
>> points(limmax,col="chartreuse4")
>>
>> # create the polygon
>> spp.farea <-
>>
>> SpatialPolygons(list(Polygons(list(Polygon(rbind(limmin,limmax))),ID='1')))
>> proj4string(spp.farea) <- CRS("+proj=longlat +datum=WGS84")
>>
>> plot(spp.farea,col="chocolate3",lwd=2,add=T)
>>
>> # select the squares under the polygon
>> fareasq <- squares[spp.farea,]
>> fareace <- gCentroid(fareasq,byid=TRUE)
>>
>> # final plot
>> plot(fareasq,axes=T)
>> points(fareace,pch=10,col="darkgreen",cex=4)
>> text(coordinates(squares), labels=IDs, cex=0.7)
>> plot(isomin,add=T,col="deepskyblue3",lwd=2)
>> plot(isomax,add=T,col="dodgerblue4",lwd=2)
>> plot(spp.area,border="red",add=T,lty=2)
>>
>>         [[alternative HTML version deleted]]
>>
>> _______________________________________________
>> R-sig-Geo mailing list
>> [hidden email]
>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>
>
>
> --
> Vijay Lulla
>
> Assistant Professor,
> Dept. of Geography, IUPUI
> 425 University Blvd, CA-207C.
> Indianapolis, IN-46202
> [hidden email]
> ORCID: https://orcid.org/0000-0002-0823-2522
> Webpage: http://vijaylulla.com
>
        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: Problems converting rasters from float to integer.

Fri, 08/24/2018 - 09:16
On Fri, 24 Aug 2018, nevil amos wrote:

> Thanks,
>
> I think it did help
> datatype vs dataType was indeed the problem!

The underlying problem is that dataType= is going into ..., and because
nothing needs it, it isn't spotted. Could someone contribute a PR to trap
this, from https://github.com/rspatial/raster?

Roger

>
> On Fri, 24 Aug 2018 at 22:14, Michael Sumner <[hidden email]> wrote:
>
>> That doesn't help though!
>>
>> Sorry, I thought that was it, will sleep on it.
>>
>> Cheers, Mike
>>
>> On Fri, 24 Aug 2018 at 22:00 Michael Sumner <[hidden email]> wrote:
>>
>>> Sorry, I was wrong - the actual issue is that the writeRaster argument is
>>> "datatype", it gets conflated with the function that is dataType().
>>>
>>> f <- "test_int.tif"
>>> writeRaster(r, f, datatype = "INT2U", overwrite = TRUE)
>>> s <- raster(f)
>>> dataType(s)
>>>
>>> [1] "INT2U"
>>>
>>> Not the first time this has caught me out. I found it by checking
>>> mode(readGDAL(f)[[1]]) and realizing that the TIFF type was Float32, prior
>>> to the next steps.
>>>
>>> HTH
>>>
>>> :)
>>>
>>> Cheers, Mike
>>>
>>> On Fri, 24 Aug 2018 at 18:26 nevil amos <[hidden email]> wrote:
>>>
>>>> Roger,
>>>>
>>>> there was a reproducible example in my post:  it seems to have got a bit
>>>> scrambled in your reply.
>>>>
>>>> from the help in raster I thought that if I included  datatype="INT2S"
>>>> in the writeRaster() command the values would be formatted in the output
>>>> file ( and the subsequent reading using raster() of that file) as integer
>>>> is this not the case?
>>>>
>>>> Here it is the reproducible example
>>>>
>>>> v<-c(rep(1.00000,25),rep(0.00000,50),rep(NA,25))
>>>> m<-matrix(v,10,10)
>>>> r<-raster(m)
>>>> dataType(r)
>>>> writeRaster(r,"test_int.tif",dataType="INT2U",overwrite=T)
>>>> s<-raster("test_int.tif")
>>>> dataType(s)
>>>>
>>>> On Fri, 24 Aug 2018 at 18:05, Roger Bivand <[hidden email]> wrote:
>>>>
>>>>> Thanks, Mike. I agree that the lower level interface in rgdal is
>>>>> flexible enough, but as you say non-trivial. Nevil: could you please
>>>>> provide a small reproducible example to point people in the right direction?
>>>>>
>>>>> Roger
>>>>>
>>>>> Roger Bivand
>>>>> Norwegian School of Economics
>>>>> Bergen, Norway
>>>>>
>>>>>
>>>>>
>>>>> Fra: Michael Sumner
>>>>> Sendt: fredag 24. august, 09.14
>>>>> Emne: Re: [R-sig-Geo] Problems converting rasters from float to integer.
>>>>> Til: nevil amos
>>>>> Kopi: [hidden email]
>>>>>
>>>>>
>>>>> I'm pretty sure readGDAL from rgdal (that raster uses) will keep as
>>>>> integers, so you can build an empty raster and copy the values over. But
>>>>> will need to derive from the rgdal version to catch all the metadata
>>>>> (structure, extent, and crs). I think it's doable and will try later. There
>>>>> are other options but nothing trivial afaik. Cheers, Mike On Fri, 24 Aug
>>>>> 2018, 16:45 nevil amos wrote: > I have a large number of rasters ( tiffs)
>>>>> that contain whole number values > between 0 and 100, and NA values. or
>>>>> 0,1,and NA > they are currently in Float format, I am trying to rewrite
>>>>> them as integer > rasters, firstly to save space, and secondly so that I
>>>>> can later read the > values stack of all the rasters into an integer array.
>>>>> using getValues(). > > To do this I am setting the dataType to INT2U in
>>>>> writeRaster, however when > I read the save file back into R the format is
>>>>> not INT2U but FLT8S > > toy example: > > v m r dataType(r) >
>>>>> writeRaster(r,"test_int.tif",dataType="INT2U",overwrite=T) > s dataType(s)
>>>>>>> the result I get: > > > v > m > r > dataType(r) > [1] "FLT4S" > >
>>>>> writeRaster(r,"test_int.tif",dataType="INT2U",overwrite=T) > > s >
>>>>> dataType(s) > [1] "FLT8S" > > > > > Can you suggest how I ensure the values
>>>>> are stored as integer? > > Many thanks > > [[alternative HTML version
>>>>> deleted]] > > _______________________________________________ > R-sig-Geo
>>>>> mailing list > [hidden email] >
>>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo > -- Dr. Michael
>>>>> Sumner Software and Database Engineer Australian Antarctic Division 203
>>>>> Channel Highway Kingston Tasmania 7050 Australia
>>>>> <https://maps.google.com/?q=203+Channel%0D%0A+Highway+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>
>>>>> [[alternative HTML version deleted]]
>>>>> _______________________________________________ R-sig-Geo mailing list
>>>>> [hidden email] https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>>>>
>>>>> --
>>> Dr. Michael Sumner
>>> Software and Database Engineer
>>> Australian Antarctic Division
>>> 203 Channel Highway
>>> Kingston Tasmania 7050 Australia
>>>
>>> --
>> Dr. Michael Sumner
>> Software and Database Engineer
>> Australian Antarctic Division
>> 203 Channel Highway
>> Kingston Tasmania 7050 Australia
>>
>>
>
--
Roger Bivand
Department of Economics, Norwegian School of Economics,
Helleveien 30, N-5045 Bergen, Norway.
voice: +47 55 95 93 55; e-mail: [hidden email]
http://orcid.org/0000-0003-2392-6140
https://scholar.google.no/citations?user=AWeghB0AAAAJ&hl=en

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Roger Bivand
Department of Economics
Norwegian School of Economics
Helleveien 30
N-5045 Bergen, Norway

Re: drawing a polygon from several lines

Fri, 08/24/2018 - 08:58
Very nice, Antonio!!  Thanks for sharing your script.  The only minor
suggestion I have to offer is that you can make your "cut isobath" portion
easier to read (at least IMO) by using `with` and `between` (either
dplyr::between or data.table::between).  Below is how I would've done it.

limmin <- with(list(coords = coordinates(isomin)[[1]][[1]]),
               subset(coords, between(coords[, 1], lonmin, lonmax) &
                              between(coords[, 2], latmin, latmax)))

Again, thanks for sharing your problem and solution.  I found it
instructive!
Cordially,
Vijay.


On Fri, Aug 24, 2018 at 9:16 AM Antonio Silva <[hidden email]> wrote:

> Dear colleagues
>
> I finally could select the squares (cells) under the polygon. I extract the
> points from the selected isobaths that fell within lat / lon limits and
> with them I created a spatial polygon.
>
> Bellow I share the script I wrote. Considering that I don't know much about
> mapping in R, I would appreciate to hear suggestions to improve the code.
>
> Best regards.
>
> --
> Antônio Olinto Ávila da Silva
> Instituto de Pesca (Fisheries Institute)
> São Paulo, Brasil
>
> ===================
>
> library(rgdal)
> library(rgeos)
> rm(list = ls())
>
> # import shapes
> # isobath
> https://www.dropbox.com/sh/us859uyijns36r9/AAC9N4OSwweLuW-fdsBWUy-ra?dl=0
> isobaths <- readOGR(".","isobath")
> isobaths <- spTransform(isobaths, CRS("+proj=longlat +datum=WGS84"))
> proj4string(isobaths)
> summary(isobaths)
>
> # create a spatial grid and polygons
> grd <-
>
> GridTopology(cellcentre.offset=c(-47.75,-25.416667),cellsize=c(10/60,10/60),cells.dim=c(23,12))
> squares <- as.SpatialPolygons.GridTopology(grd)
> proj4string(squares) <- CRS("+proj=longlat +datum=WGS84")
> summary(squares)
> IDs <- sapply(slot(squares, "polygons"), function(x) slot(x, "ID"))
>
> # set limits
> latmin <- -24.8
> latmax <- -24.02
> lonmin <- -47.2
> lonmax <- -44.8
> depmin <- -25
> depmax <- -60
>
> # bounding box
> mat.area <-
>
> matrix(c(lonmin,latmax,lonmax,latmax,lonmax,latmin,lonmin,latmin),ncol=2,byrow=T)
> spp.area <- SpatialPolygons(list(Polygons(list(Polygon(mat.area)) ,
> ID='1')))
> proj4string(spp.area) <- CRS("+proj=longlat +datum=WGS84")
>
> # plot polygons isolines and bounding box
> plot(squares, axes=T)
> #~ text(coordinates(squares),labels=IDs, cex=0.7)
> plot(isobaths,add=T)
> plot(spp.area,border="red",add=T)
>
> # select isobaths
> summary(isobaths)
> isomin <- subset(isobaths, ID %in% depmin)
> proj4string(isomin) <- CRS("+proj=longlat +datum=WGS84")
> isomax <- subset(isobaths, ID %in% depmax)
> proj4string(isomax) <- CRS("+proj=longlat +datum=WGS84")
>
> plot(isomin,add=T,col="deepskyblue3",lwd=2)
> plot(isomax,add=T,col="dodgerblue4",lwd=2)
>
> # cut isobaths
> limmin <- subset(coordinates(isomin)[[1]][[1]],
>           coordinates(isomin)[[1]][[1]][,2]<=latmax &
> coordinates(isomin)[[1]][[1]][,2]>=latmin &
>           coordinates(isomin)[[1]][[1]][,1]<=lonmax &
> coordinates(isomin)[[1]][[1]][,1]>=lonmin)
> limmax <- subset(coordinates(isomax)[[1]][[1]],
>           coordinates(isomax)[[1]][[1]][,2]<=latmax &
> coordinates(isomax)[[1]][[1]][,2]>=latmin &
>           coordinates(isomax)[[1]][[1]][,1]<=lonmax &
> coordinates(isomax)[[1]][[1]][,1]>=lonmin)
>
> points(limmin,col="chartreuse3")
> points(limmax,col="chartreuse4")
>
> # create the polygon
> spp.farea <-
> SpatialPolygons(list(Polygons(list(Polygon(rbind(limmin,limmax))),ID='1')))
> proj4string(spp.farea) <- CRS("+proj=longlat +datum=WGS84")
>
> plot(spp.farea,col="chocolate3",lwd=2,add=T)
>
> # select the squares under the polygon
> fareasq <- squares[spp.farea,]
> fareace <- gCentroid(fareasq,byid=TRUE)
>
> # final plot
> plot(fareasq,axes=T)
> points(fareace,pch=10,col="darkgreen",cex=4)
> text(coordinates(squares), labels=IDs, cex=0.7)
> plot(isomin,add=T,col="deepskyblue3",lwd=2)
> plot(isomax,add=T,col="dodgerblue4",lwd=2)
> plot(spp.area,border="red",add=T,lty=2)
>
>         [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> [hidden email]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>

--
Vijay Lulla

Assistant Professor,
Dept. of Geography, IUPUI
425 University Blvd, CA-207C.
Indianapolis, IN-46202
[hidden email]
ORCID: https://orcid.org/0000-0002-0823-2522
Webpage: http://vijaylulla.com

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: Problems converting rasters from float to integer.

Fri, 08/24/2018 - 08:18
Thanks,

I think it did help
datatype vs dataType was indeed the problem!

On Fri, 24 Aug 2018 at 22:14, Michael Sumner <[hidden email]> wrote:

> That doesn't help though!
>
> Sorry, I thought that was it, will sleep on it.
>
> Cheers, Mike
>
> On Fri, 24 Aug 2018 at 22:00 Michael Sumner <[hidden email]> wrote:
>
>> Sorry, I was wrong - the actual issue is that the writeRaster argument is
>> "datatype", it gets conflated with the function that is dataType().
>>
>> f <- "test_int.tif"
>> writeRaster(r, f, datatype = "INT2U", overwrite = TRUE)
>> s <- raster(f)
>> dataType(s)
>>
>> [1] "INT2U"
>>
>> Not the first time this has caught me out. I found it by checking
>> mode(readGDAL(f)[[1]]) and realizing that the TIFF type was Float32, prior
>> to the next steps.
>>
>> HTH
>>
>> :)
>>
>> Cheers, Mike
>>
>> On Fri, 24 Aug 2018 at 18:26 nevil amos <[hidden email]> wrote:
>>
>>> Roger,
>>>
>>> there was a reproducible example in my post:  it seems to have got a bit
>>> scrambled in your reply.
>>>
>>> from the help in raster I thought that if I included  datatype="INT2S"
>>> in the writeRaster() command the values would be formatted in the output
>>> file ( and the subsequent reading using raster() of that file) as integer
>>> is this not the case?
>>>
>>> Here it is the reproducible example
>>>
>>> v<-c(rep(1.00000,25),rep(0.00000,50),rep(NA,25))
>>> m<-matrix(v,10,10)
>>> r<-raster(m)
>>> dataType(r)
>>> writeRaster(r,"test_int.tif",dataType="INT2U",overwrite=T)
>>> s<-raster("test_int.tif")
>>> dataType(s)
>>>
>>> On Fri, 24 Aug 2018 at 18:05, Roger Bivand <[hidden email]> wrote:
>>>
>>>> Thanks, Mike. I agree that the lower level interface in rgdal is
>>>> flexible enough, but as you say non-trivial. Nevil: could you please
>>>> provide a small reproducible example to point people in the right direction?
>>>>
>>>> Roger
>>>>
>>>> Roger Bivand
>>>> Norwegian School of Economics
>>>> Bergen, Norway
>>>>
>>>>
>>>>
>>>> Fra: Michael Sumner
>>>> Sendt: fredag 24. august, 09.14
>>>> Emne: Re: [R-sig-Geo] Problems converting rasters from float to integer.
>>>> Til: nevil amos
>>>> Kopi: [hidden email]
>>>>
>>>>
>>>> I'm pretty sure readGDAL from rgdal (that raster uses) will keep as
>>>> integers, so you can build an empty raster and copy the values over. But
>>>> will need to derive from the rgdal version to catch all the metadata
>>>> (structure, extent, and crs). I think it's doable and will try later. There
>>>> are other options but nothing trivial afaik. Cheers, Mike On Fri, 24 Aug
>>>> 2018, 16:45 nevil amos wrote: > I have a large number of rasters ( tiffs)
>>>> that contain whole number values > between 0 and 100, and NA values. or
>>>> 0,1,and NA > they are currently in Float format, I am trying to rewrite
>>>> them as integer > rasters, firstly to save space, and secondly so that I
>>>> can later read the > values stack of all the rasters into an integer array.
>>>> using getValues(). > > To do this I am setting the dataType to INT2U in
>>>> writeRaster, however when > I read the save file back into R the format is
>>>> not INT2U but FLT8S > > toy example: > > v m r dataType(r) >
>>>> writeRaster(r,"test_int.tif",dataType="INT2U",overwrite=T) > s dataType(s)
>>>> > > the result I get: > > > v > m > r > dataType(r) > [1] "FLT4S" > >
>>>> writeRaster(r,"test_int.tif",dataType="INT2U",overwrite=T) > > s >
>>>> dataType(s) > [1] "FLT8S" > > > > > Can you suggest how I ensure the values
>>>> are stored as integer? > > Many thanks > > [[alternative HTML version
>>>> deleted]] > > _______________________________________________ > R-sig-Geo
>>>> mailing list > [hidden email] >
>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo > -- Dr. Michael
>>>> Sumner Software and Database Engineer Australian Antarctic Division 203
>>>> Channel Highway Kingston Tasmania 7050 Australia
>>>> <https://maps.google.com/?q=203+Channel%0D%0A+Highway+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>
>>>> [[alternative HTML version deleted]]
>>>> _______________________________________________ R-sig-Geo mailing list
>>>> [hidden email] https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>>>
>>>> --
>> Dr. Michael Sumner
>> Software and Database Engineer
>> Australian Antarctic Division
>> 203 Channel Highway
>> Kingston Tasmania 7050 Australia
>>
>> --
> Dr. Michael Sumner
> Software and Database Engineer
> Australian Antarctic Division
> 203 Channel Highway
> Kingston Tasmania 7050 Australia
>
>
        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: drawing a polygon from several lines

Fri, 08/24/2018 - 08:15
Dear colleagues

I finally could select the squares (cells) under the polygon. I extract the
points from the selected isobaths that fell within lat / lon limits and
with them I created a spatial polygon.

Bellow I share the script I wrote. Considering that I don't know much about
mapping in R, I would appreciate to hear suggestions to improve the code.

Best regards.

--
Antônio Olinto Ávila da Silva
Instituto de Pesca (Fisheries Institute)
São Paulo, Brasil

===================

library(rgdal)
library(rgeos)
rm(list = ls())

# import shapes
# isobath
https://www.dropbox.com/sh/us859uyijns36r9/AAC9N4OSwweLuW-fdsBWUy-ra?dl=0
isobaths <- readOGR(".","isobath")
isobaths <- spTransform(isobaths, CRS("+proj=longlat +datum=WGS84"))
proj4string(isobaths)
summary(isobaths)

# create a spatial grid and polygons
grd <-
GridTopology(cellcentre.offset=c(-47.75,-25.416667),cellsize=c(10/60,10/60),cells.dim=c(23,12))
squares <- as.SpatialPolygons.GridTopology(grd)
proj4string(squares) <- CRS("+proj=longlat +datum=WGS84")
summary(squares)
IDs <- sapply(slot(squares, "polygons"), function(x) slot(x, "ID"))

# set limits
latmin <- -24.8
latmax <- -24.02
lonmin <- -47.2
lonmax <- -44.8
depmin <- -25
depmax <- -60

# bounding box
mat.area <-
matrix(c(lonmin,latmax,lonmax,latmax,lonmax,latmin,lonmin,latmin),ncol=2,byrow=T)
spp.area <- SpatialPolygons(list(Polygons(list(Polygon(mat.area)) ,
ID='1')))
proj4string(spp.area) <- CRS("+proj=longlat +datum=WGS84")

# plot polygons isolines and bounding box
plot(squares, axes=T)
#~ text(coordinates(squares),labels=IDs, cex=0.7)
plot(isobaths,add=T)
plot(spp.area,border="red",add=T)

# select isobaths
summary(isobaths)
isomin <- subset(isobaths, ID %in% depmin)
proj4string(isomin) <- CRS("+proj=longlat +datum=WGS84")
isomax <- subset(isobaths, ID %in% depmax)
proj4string(isomax) <- CRS("+proj=longlat +datum=WGS84")

plot(isomin,add=T,col="deepskyblue3",lwd=2)
plot(isomax,add=T,col="dodgerblue4",lwd=2)

# cut isobaths
limmin <- subset(coordinates(isomin)[[1]][[1]],
          coordinates(isomin)[[1]][[1]][,2]<=latmax &
coordinates(isomin)[[1]][[1]][,2]>=latmin &
          coordinates(isomin)[[1]][[1]][,1]<=lonmax &
coordinates(isomin)[[1]][[1]][,1]>=lonmin)
limmax <- subset(coordinates(isomax)[[1]][[1]],
          coordinates(isomax)[[1]][[1]][,2]<=latmax &
coordinates(isomax)[[1]][[1]][,2]>=latmin &
          coordinates(isomax)[[1]][[1]][,1]<=lonmax &
coordinates(isomax)[[1]][[1]][,1]>=lonmin)

points(limmin,col="chartreuse3")
points(limmax,col="chartreuse4")

# create the polygon
spp.farea <-
SpatialPolygons(list(Polygons(list(Polygon(rbind(limmin,limmax))),ID='1')))
proj4string(spp.farea) <- CRS("+proj=longlat +datum=WGS84")

plot(spp.farea,col="chocolate3",lwd=2,add=T)

# select the squares under the polygon
fareasq <- squares[spp.farea,]
fareace <- gCentroid(fareasq,byid=TRUE)

# final plot
plot(fareasq,axes=T)
points(fareace,pch=10,col="darkgreen",cex=4)
text(coordinates(squares), labels=IDs, cex=0.7)
plot(isomin,add=T,col="deepskyblue3",lwd=2)
plot(isomax,add=T,col="dodgerblue4",lwd=2)
plot(spp.area,border="red",add=T,lty=2)

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: Problems converting rasters from float to integer.

Fri, 08/24/2018 - 07:13
That doesn't help though!

Sorry, I thought that was it, will sleep on it.

Cheers, Mike

On Fri, 24 Aug 2018 at 22:00 Michael Sumner <[hidden email]> wrote:

> Sorry, I was wrong - the actual issue is that the writeRaster argument is
> "datatype", it gets conflated with the function that is dataType().
>
> f <- "test_int.tif"
> writeRaster(r, f, datatype = "INT2U", overwrite = TRUE)
> s <- raster(f)
> dataType(s)
>
> [1] "INT2U"
>
> Not the first time this has caught me out. I found it by checking
> mode(readGDAL(f)[[1]]) and realizing that the TIFF type was Float32, prior
> to the next steps.
>
> HTH
>
> :)
>
> Cheers, Mike
>
> On Fri, 24 Aug 2018 at 18:26 nevil amos <[hidden email]> wrote:
>
>> Roger,
>>
>> there was a reproducible example in my post:  it seems to have got a bit
>> scrambled in your reply.
>>
>> from the help in raster I thought that if I included  datatype="INT2S" in
>> the writeRaster() command the values would be formatted in the output file
>> ( and the subsequent reading using raster() of that file) as integer  is
>> this not the case?
>>
>> Here it is the reproducible example
>>
>> v<-c(rep(1.00000,25),rep(0.00000,50),rep(NA,25))
>> m<-matrix(v,10,10)
>> r<-raster(m)
>> dataType(r)
>> writeRaster(r,"test_int.tif",dataType="INT2U",overwrite=T)
>> s<-raster("test_int.tif")
>> dataType(s)
>>
>> On Fri, 24 Aug 2018 at 18:05, Roger Bivand <[hidden email]> wrote:
>>
>>> Thanks, Mike. I agree that the lower level interface in rgdal is
>>> flexible enough, but as you say non-trivial. Nevil: could you please
>>> provide a small reproducible example to point people in the right direction?
>>>
>>> Roger
>>>
>>> Roger Bivand
>>> Norwegian School of Economics
>>> Bergen, Norway
>>>
>>>
>>>
>>> Fra: Michael Sumner
>>> Sendt: fredag 24. august, 09.14
>>> Emne: Re: [R-sig-Geo] Problems converting rasters from float to integer.
>>> Til: nevil amos
>>> Kopi: [hidden email]
>>>
>>>
>>> I'm pretty sure readGDAL from rgdal (that raster uses) will keep as
>>> integers, so you can build an empty raster and copy the values over. But
>>> will need to derive from the rgdal version to catch all the metadata
>>> (structure, extent, and crs). I think it's doable and will try later. There
>>> are other options but nothing trivial afaik. Cheers, Mike On Fri, 24 Aug
>>> 2018, 16:45 nevil amos wrote: > I have a large number of rasters ( tiffs)
>>> that contain whole number values > between 0 and 100, and NA values. or
>>> 0,1,and NA > they are currently in Float format, I am trying to rewrite
>>> them as integer > rasters, firstly to save space, and secondly so that I
>>> can later read the > values stack of all the rasters into an integer array.
>>> using getValues(). > > To do this I am setting the dataType to INT2U in
>>> writeRaster, however when > I read the save file back into R the format is
>>> not INT2U but FLT8S > > toy example: > > v m r dataType(r) >
>>> writeRaster(r,"test_int.tif",dataType="INT2U",overwrite=T) > s dataType(s)
>>> > > the result I get: > > > v > m > r > dataType(r) > [1] "FLT4S" > >
>>> writeRaster(r,"test_int.tif",dataType="INT2U",overwrite=T) > > s >
>>> dataType(s) > [1] "FLT8S" > > > > > Can you suggest how I ensure the values
>>> are stored as integer? > > Many thanks > > [[alternative HTML version
>>> deleted]] > > _______________________________________________ > R-sig-Geo
>>> mailing list > [hidden email] >
>>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo > -- Dr. Michael Sumner
>>> Software and Database Engineer Australian Antarctic Division 203
>>> Channel Highway Kingston Tasmania 7050 Australia
>>> <https://maps.google.com/?q=203+Channel%0D%0A+Highway+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>
>>> [[alternative HTML version deleted]]
>>> _______________________________________________ R-sig-Geo mailing list
>>> [hidden email] https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>>
>>> --
> Dr. Michael Sumner
> Software and Database Engineer
> Australian Antarctic Division
> 203 Channel Highway
> Kingston Tasmania 7050 Australia
>
> -- Dr. Michael Sumner
Software and Database Engineer
Australian Antarctic Division
203 Channel Highway
Kingston Tasmania 7050 Australia

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: Problems converting rasters from float to integer.

Fri, 08/24/2018 - 07:00
Sorry, I was wrong - the actual issue is that the writeRaster argument is
"datatype", it gets conflated with the function that is dataType().

f <- "test_int.tif"
writeRaster(r, f, datatype = "INT2U", overwrite = TRUE)
s <- raster(f)
dataType(s)

[1] "INT2U"

Not the first time this has caught me out. I found it by checking
mode(readGDAL(f)[[1]]) and realizing that the TIFF type was Float32, prior
to the next steps.

HTH

:)

Cheers, Mike

On Fri, 24 Aug 2018 at 18:26 nevil amos <[hidden email]> wrote:

> Roger,
>
> there was a reproducible example in my post:  it seems to have got a bit
> scrambled in your reply.
>
> from the help in raster I thought that if I included  datatype="INT2S" in
> the writeRaster() command the values would be formatted in the output file
> ( and the subsequent reading using raster() of that file) as integer  is
> this not the case?
>
> Here it is the reproducible example
>
> v<-c(rep(1.00000,25),rep(0.00000,50),rep(NA,25))
> m<-matrix(v,10,10)
> r<-raster(m)
> dataType(r)
> writeRaster(r,"test_int.tif",dataType="INT2U",overwrite=T)
> s<-raster("test_int.tif")
> dataType(s)
>
> On Fri, 24 Aug 2018 at 18:05, Roger Bivand <[hidden email]> wrote:
>
>> Thanks, Mike. I agree that the lower level interface in rgdal is flexible
>> enough, but as you say non-trivial. Nevil: could you please provide a small
>> reproducible example to point people in the right direction?
>>
>> Roger
>>
>> Roger Bivand
>> Norwegian School of Economics
>> Bergen, Norway
>>
>>
>>
>> Fra: Michael Sumner
>> Sendt: fredag 24. august, 09.14
>> Emne: Re: [R-sig-Geo] Problems converting rasters from float to integer.
>> Til: nevil amos
>> Kopi: [hidden email]
>>
>>
>> I'm pretty sure readGDAL from rgdal (that raster uses) will keep as
>> integers, so you can build an empty raster and copy the values over. But
>> will need to derive from the rgdal version to catch all the metadata
>> (structure, extent, and crs). I think it's doable and will try later. There
>> are other options but nothing trivial afaik. Cheers, Mike On Fri, 24 Aug
>> 2018, 16:45 nevil amos wrote: > I have a large number of rasters ( tiffs)
>> that contain whole number values > between 0 and 100, and NA values. or
>> 0,1,and NA > they are currently in Float format, I am trying to rewrite
>> them as integer > rasters, firstly to save space, and secondly so that I
>> can later read the > values stack of all the rasters into an integer array.
>> using getValues(). > > To do this I am setting the dataType to INT2U in
>> writeRaster, however when > I read the save file back into R the format is
>> not INT2U but FLT8S > > toy example: > > v m r dataType(r) >
>> writeRaster(r,"test_int.tif",dataType="INT2U",overwrite=T) > s dataType(s)
>> > > the result I get: > > > v > m > r > dataType(r) > [1] "FLT4S" > >
>> writeRaster(r,"test_int.tif",dataType="INT2U",overwrite=T) > > s >
>> dataType(s) > [1] "FLT8S" > > > > > Can you suggest how I ensure the values
>> are stored as integer? > > Many thanks > > [[alternative HTML version
>> deleted]] > > _______________________________________________ > R-sig-Geo
>> mailing list > [hidden email] >
>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo > -- Dr. Michael Sumner
>> Software and Database Engineer Australian Antarctic Division 203 Channel
>> Highway Kingston Tasmania 7050 Australia
>> <https://maps.google.com/?q=203+Channel%0D%0A+Highway+Kingston+Tasmania+7050+Australia&entry=gmail&source=g>
>> [[alternative HTML version deleted]]
>> _______________________________________________ R-sig-Geo mailing list
>> [hidden email] https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>
>> -- Dr. Michael Sumner
Software and Database Engineer
Australian Antarctic Division
203 Channel Highway
Kingston Tasmania 7050 Australia

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: Problems converting rasters from float to integer.

Fri, 08/24/2018 - 03:26
Roger,

there was a reproducible example in my post:  it seems to have got a bit
scrambled in your reply.

from the help in raster I thought that if I included  datatype="INT2S" in
the writeRaster() command the values would be formatted in the output file
( and the subsequent reading using raster() of that file) as integer  is
this not the case?

Here it is the reproducible example

v<-c(rep(1.00000,25),rep(0.00000,50),rep(NA,25))
m<-matrix(v,10,10)
r<-raster(m)
dataType(r)
writeRaster(r,"test_int.tif",dataType="INT2U",overwrite=T)
s<-raster("test_int.tif")
dataType(s)

On Fri, 24 Aug 2018 at 18:05, Roger Bivand <[hidden email]> wrote:

> Thanks, Mike. I agree that the lower level interface in rgdal is flexible
> enough, but as you say non-trivial. Nevil: could you please provide a small
> reproducible example to point people in the right direction?
>
> Roger
>
> Roger Bivand
> Norwegian School of Economics
> Bergen, Norway
>
>
>
> Fra: Michael Sumner
> Sendt: fredag 24. august, 09.14
> Emne: Re: [R-sig-Geo] Problems converting rasters from float to integer.
> Til: nevil amos
> Kopi: [hidden email]
>
>
> I'm pretty sure readGDAL from rgdal (that raster uses) will keep as
> integers, so you can build an empty raster and copy the values over. But
> will need to derive from the rgdal version to catch all the metadata
> (structure, extent, and crs). I think it's doable and will try later. There
> are other options but nothing trivial afaik. Cheers, Mike On Fri, 24 Aug
> 2018, 16:45 nevil amos wrote: > I have a large number of rasters ( tiffs)
> that contain whole number values > between 0 and 100, and NA values. or
> 0,1,and NA > they are currently in Float format, I am trying to rewrite
> them as integer > rasters, firstly to save space, and secondly so that I
> can later read the > values stack of all the rasters into an integer array.
> using getValues(). > > To do this I am setting the dataType to INT2U in
> writeRaster, however when > I read the save file back into R the format is
> not INT2U but FLT8S > > toy example: > > v m r dataType(r) >
> writeRaster(r,"test_int.tif",dataType="INT2U",overwrite=T) > s dataType(s)
> > > the result I get: > > > v > m > r > dataType(r) > [1] "FLT4S" > >
> writeRaster(r,"test_int.tif",dataType="INT2U",overwrite=T) > > s >
> dataType(s) > [1] "FLT8S" > > > > > Can you suggest how I ensure the values
> are stored as integer? > > Many thanks > > [[alternative HTML version
> deleted]] > > _______________________________________________ > R-sig-Geo
> mailing list > [hidden email] >
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo > -- Dr. Michael Sumner
> Software and Database Engineer Australian Antarctic Division 203 Channel
> Highway Kingston Tasmania 7050 Australia [[alternative HTML version
> deleted]] _______________________________________________ R-sig-Geo mailing
> list [hidden email]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>
>
        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: Problems converting rasters from float to integer.

Fri, 08/24/2018 - 03:05
Thanks, Mike. I agree that the lower level interface in rgdal is flexible enough, but as you say non-trivial. Nevil: could you please provide a small reproducible example to point people in the right direction?

Roger

Roger Bivand
Norwegian School of Economics
Bergen, Norway



Fra: Michael Sumner
Sendt: fredag 24. august, 09.14
Emne: Re: [R-sig-Geo] Problems converting rasters from float to integer.
Til: nevil amos
Kopi: [hidden email]


I'm pretty sure readGDAL from rgdal (that raster uses) will keep as integers, so you can build an empty raster and copy the values over. But will need to derive from the rgdal version to catch all the metadata (structure, extent, and crs). I think it's doable and will try later. There are other options but nothing trivial afaik. Cheers, Mike On Fri, 24 Aug 2018, 16:45 nevil amos wrote: > I have a large number of rasters ( tiffs) that contain whole number values > between 0 and 100, and NA values. or 0,1,and NA > they are currently in Float format, I am trying to rewrite them as integer > rasters, firstly to save space, and secondly so that I can later read the > values stack of all the rasters into an integer array. using getValues(). > > To do this I am setting the dataType to INT2U in writeRaster, however when > I read the save file back into R the format is not INT2U but FLT8S > > toy example: > > v m r dataType(r) > writeRaster(r,"test_int.tif",dataType="INT2U",overwrite=T) > s d
 ataType(s) > > the result I get: > > > v > m > r > dataType(r) > [1] "FLT4S" > > writeRaster(r,"test_int.tif",dataType="INT2U",overwrite=T) > > s > dataType(s) > [1] "FLT8S" > > > > > Can you suggest how I ensure the values are stored as integer? > > Many thanks > > [[alternative HTML version deleted]] > > _______________________________________________ > R-sig-Geo mailing list > [hidden email] > https://stat.ethz.ch/mailman/listinfo/r-sig-geo > -- Dr. Michael Sumner Software and Database Engineer Australian Antarctic Division 203 Channel Highway Kingston Tasmania 7050 Australia [[alternative HTML version deleted]] _______________________________________________ R-sig-Geo mailing list [hidden email] https://stat.ethz.ch/mailman/listinfo/r-sig-geo


        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Roger Bivand
Department of Economics
Norwegian School of Economics
Helleveien 30
N-5045 Bergen, Norway

Re: Problems converting rasters from float to integer.

Fri, 08/24/2018 - 02:13
I'm pretty sure readGDAL from rgdal (that raster uses) will keep as
integers, so you can build an empty raster and copy the values over. But
will need to derive from the rgdal version to catch all the metadata
(structure, extent, and crs).

I think it's doable and will try later. There are other options but nothing
trivial afaik.

Cheers, Mike

On Fri, 24 Aug 2018, 16:45 nevil amos <[hidden email]> wrote:

> I have a large number of rasters ( tiffs) that contain whole number values
> between 0 and 100, and NA values. or 0,1,and NA
> they are currently in Float format, I am trying to rewrite them as integer
> rasters, firstly to save space, and secondly so that I can later read the
> values stack of all the rasters into an integer array. using getValues().
>
> To do this I am setting the dataType to INT2U in writeRaster, however when
> I read the save file back into R the format is not INT2U but FLT8S
>
> toy example:
>
> v<-c(rep(1.00000,25),rep(0.00000,50),rep(NA,25))
> m<-matrix(v,10,10)
> r<-raster(m)
> dataType(r)
> writeRaster(r,"test_int.tif",dataType="INT2U",overwrite=T)
> s<-raster("test_int.tif")
> dataType(s)
>
> the result I get:
>
> > v<-c(rep(1.00000,25),rep(0.00000,50),rep(NA,25))
> > m<-matrix(v,10,10)
> > r<-raster(m)
> > dataType(r)
> [1] "FLT4S"
> > writeRaster(r,"test_int.tif",dataType="INT2U",overwrite=T)
> > s<-raster("test_int.tif")
> > dataType(s)
> [1] "FLT8S"
> >
>
>
> Can you suggest how I ensure the values are stored as integer?
>
> Many thanks
>
>         [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> [hidden email]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> --
Dr. Michael Sumner
Software and Database Engineer
Australian Antarctic Division
203 Channel Highway
Kingston Tasmania 7050 Australia

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Problems converting rasters from float to integer.

Fri, 08/24/2018 - 01:44
I have a large number of rasters ( tiffs) that contain whole number values
between 0 and 100, and NA values. or 0,1,and NA
they are currently in Float format, I am trying to rewrite them as integer
rasters, firstly to save space, and secondly so that I can later read the
values stack of all the rasters into an integer array. using getValues().

To do this I am setting the dataType to INT2U in writeRaster, however when
I read the save file back into R the format is not INT2U but FLT8S

toy example:

v<-c(rep(1.00000,25),rep(0.00000,50),rep(NA,25))
m<-matrix(v,10,10)
r<-raster(m)
dataType(r)
writeRaster(r,"test_int.tif",dataType="INT2U",overwrite=T)
s<-raster("test_int.tif")
dataType(s)

the result I get:

> v<-c(rep(1.00000,25),rep(0.00000,50),rep(NA,25))
> m<-matrix(v,10,10)
> r<-raster(m)
> dataType(r)
[1] "FLT4S"
> writeRaster(r,"test_int.tif",dataType="INT2U",overwrite=T)
> s<-raster("test_int.tif")
> dataType(s)
[1] "FLT8S"
>


Can you suggest how I ensure the values are stored as integer?

Many thanks

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

R-movement packages survey

Thu, 08/23/2018 - 09:59
Hi everybody,

I'm preparing a review on R packages related to biologging and movement.
If you have used a package related to those subjects please participate
in this survey https://urldefense.proofpoint.com/v2/url?u=https-3A__mablab.org_survey_index.php_759774-3Flang-3Den&d=DwIDaQ&c=pZJPUDQ3SB9JplYbifm4nt2lEVG5pWx2KikqINpWlZM&r=OnBRYbn-fPC7OCR-yjUnqiUnj2gfCrtQn15ueJyGfIM&m=xGMFhejJGmRHQO806DDAuDcHsu0xqrrQOb1f-ZUdOro&s=ApjY8blfgbDLnrphzCM4BJGRd_4-hBqHhwwI9Kzyxbg&e= that
will be very useful for the review.

Thanks in advance.

Best,

--
Rocío Joo

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

drawing a polygon from several lines

Wed, 08/22/2018 - 16:34
Hello all,

I want to select the cells of a grid that in a range of latitudes and
depths.
Latitudes limits are indicated by two lines and the depths by two isobaths
from a shapefile.
The shapefile can be downloaded at
https://www.dropbox.com/sh/us859uyijns36r9/AAC9N4OSwweLuW-fdsBWUy-ra?dl=0

Below are the steps I've already taken. I could not find a way to create a
polygon from the lines of depth and latitude. With this polygon I could
easily select the overlapped cells.

I would really appreciate any help. Best regards

----
library(rgeos)
library(sp)
rm(list = ls())

# import shape
isobaths <- readOGR(".","isobath")
isobaths <- spTransform(bath, CRS("+proj=longlat +datum=WGS84"))
proj4string(isobaths)
summary(isobaths)

# create a spatial grid and polygons
grd <-
GridTopology(cellcentre.offset=c(-47.75,-25.416667),cellsize=c(10/60,10/60),cells.dim=c(23,12))
polys <- as.SpatialPolygons.GridTopology(grd)
proj4string(polys) <- CRS("+proj=longlat +datum=WGS84")
summary(polys)
IDs <- sapply(slot(polys, "polygons"), function(x) slot(x, "ID"))

# create spatial lines
# southern limit
x <- c(-47.5,-45.5)
y <- c(-24.8,-24.8)
ls <- SpatialLines(list(Lines(Line(cbind(x,y)),
ID="b")),proj4string=CRS("+proj=longlat +datum=WGS84"))

# northern limit
x <- c(-46.2,-44.5)
y <- c(-24.02,-24.02)
ln <- SpatialLines(list(Lines(Line(cbind(x,y)),
ID="c")),proj4string=CRS("+proj=longlat +datum=WGS84"))

# plot isolines polygons and lines
plot(polys)
axis(1);axis(4)
text(coordinates(polys), labels=IDs, cex=0.7)
plot(isobaths,add=T)
plot(ls,add=T,col="blue",lwd=2)
plot(ln,add=T,col="orange",lwd=2)

# select isobaths
summary(isobaths)
dmin <- subset(isobaths, ID %in% "-25")
proj4string(dmin) <- CRS("+proj=longlat +datum=WGS84")

dmax <- subset(isobaths, ID %in% "-60")
proj4string(dmax) <- CRS("+proj=longlat +datum=WGS84")

plot(dmin,add=T,col="red",lwd=2)
plot(dmax,add=T,col="red",lwd=2)

# from here I'd like to have a polygon to select the cells.
na.omit(over(polys,dmin))

--
Antonio Olinto Ávila da Silva
Instituto de Pesca (Fisheries Institute)
São Paulo, Brasil

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: parallelize distance matrix

Tue, 08/21/2018 - 17:12
You'll have to do clusterExport(cl, c("ll")) before you call parRapply.
Please see ?`parRapply` and pages 10 and 11 of the parallel::parallel
vignette.
HTH,
Vijay.

On Tue, Aug 21, 2018 at 2:54 PM Roman Luštrik <[hidden email]>
wrote:

> Cross posted on SO:
>
> https://stackoverflow.com/questions/51952776/parallelization-apply-to-parrapply
>
> On Tue, Aug 21, 2018 at 8:15 PM Ariel Fuentesdi <[hidden email]>
> wrote:
>
> > Hi, I want to parallelize the calculation of a distance matrix, I've
> tried
> > but I wasn't succesful, this is what I did with a training data.
> >
> > My data set is:
> >
> > ll <- matrix(c(5, 6, 60, 60), ncol=2)
> >
> > And I use the function *spDistsN1* from the library *"sp"* to obtain a
> > distance matrix with *apply*:
> >
> > apply(ll, 1, function(x) spDistsN1(as.matrix(ll), x, longlat = T))
> >
> > But I want to do it with parallelization, so for that:
> >
> > library(parallel)
> > ncore <- detectCores()
> > cl <- makeCluster(ncore)
> > clusterEvalQ(cl = cl, expr = c(library(sp)))
> > parRapply(cl = cl, x = ll, FUN =  function(x) spDistsN1(as.matrix(ll), x,
> > longlat = T))
> >
> > It shows the following error:
> >
> > *Error in checkForRemoteErrors(val) : 4 nodes produced errors; first
> error:
> > object 'll' not found*
> >
> > How do I fix it?
> > Regards,
> > Ariel Fuentes
> >
> >         [[alternative HTML version deleted]]
> >
> > _______________________________________________
> > R-sig-Geo mailing list
> > [hidden email]
> > https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> >
>
>
> --
> In God we trust, all others bring data.
>
>         [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> [hidden email]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>
        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: parallelize distance matrix

Tue, 08/21/2018 - 13:53
Cross posted on SO:
https://stackoverflow.com/questions/51952776/parallelization-apply-to-parrapply

On Tue, Aug 21, 2018 at 8:15 PM Ariel Fuentesdi <[hidden email]>
wrote:

> Hi, I want to parallelize the calculation of a distance matrix, I've tried
> but I wasn't succesful, this is what I did with a training data.
>
> My data set is:
>
> ll <- matrix(c(5, 6, 60, 60), ncol=2)
>
> And I use the function *spDistsN1* from the library *"sp"* to obtain a
> distance matrix with *apply*:
>
> apply(ll, 1, function(x) spDistsN1(as.matrix(ll), x, longlat = T))
>
> But I want to do it with parallelization, so for that:
>
> library(parallel)
> ncore <- detectCores()
> cl <- makeCluster(ncore)
> clusterEvalQ(cl = cl, expr = c(library(sp)))
> parRapply(cl = cl, x = ll, FUN =  function(x) spDistsN1(as.matrix(ll), x,
> longlat = T))
>
> It shows the following error:
>
> *Error in checkForRemoteErrors(val) : 4 nodes produced errors; first error:
> object 'll' not found*
>
> How do I fix it?
> Regards,
> Ariel Fuentes
>
>         [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-Geo mailing list
> [hidden email]
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>

--
In God we trust, all others bring data.

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

parallelize distance matrix

Tue, 08/21/2018 - 13:15
Hi, I want to parallelize the calculation of a distance matrix, I've tried
but I wasn't succesful, this is what I did with a training data.

My data set is:

ll <- matrix(c(5, 6, 60, 60), ncol=2)

And I use the function *spDistsN1* from the library *"sp"* to obtain a
distance matrix with *apply*:

apply(ll, 1, function(x) spDistsN1(as.matrix(ll), x, longlat = T))

But I want to do it with parallelization, so for that:

library(parallel)
ncore <- detectCores()
cl <- makeCluster(ncore)
clusterEvalQ(cl = cl, expr = c(library(sp)))
parRapply(cl = cl, x = ll, FUN =  function(x) spDistsN1(as.matrix(ll), x,
longlat = T))

It shows the following error:

*Error in checkForRemoteErrors(val) : 4 nodes produced errors; first error:
object 'll' not found*

How do I fix it?
Regards,
Ariel Fuentes

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: www.asdar-book.org site issues

Tue, 08/21/2018 - 07:22
The website seems to be up again.

On 08/11/2018 01:36 PM, Roger Bivand wrote:
> The ASDAR book site is not being shown for viewing, but the content is
> still there. Anyone needing content (several inquiries off-list) may try:
>
> ASDAR_BOOK <- "http://www.asdar-book.org/bundles2ed"
> chapters <- c("hello", "cm", "vis", "die", "cm2", "std", "sppa", "geos",
>  "lat", "dismap")
> for (i in chapters) {fn <- paste(i, "bundle.zip", sep="_");
> download.file(paste(ASDAR_BOOK, fn, sep = "/"), fn)}
>
> to download the chapter bundles. Summer maintenance takes time, and we
> hope to fix the underlying problem in a week or so.
>
> Roger
>
--
Edzer Pebesma
Institute for Geoinformatics
Heisenbergstrasse 2, 48151 Muenster, Germany
Phone: +49 251 8333081

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

[jobs] postdoc in remote sensing with USDA-ARS

Mon, 08/20/2018 - 15:06
The USDA-Agricultural Research Service, Southeast Watershed Research Unit
in Tifton, GA, is seeking a postdoctoral research associate for a
full-time, 2 year appointment. The research associate will develop,
calibrate, and validate a semi-empirical model for estimating cotton
biomass from fields within the southeastern United States. Model
development will involve the use of optical and radar remote sensing
products along with field data for calibration and validation. The research
associate will conduct field work and will collaborate with model
developers in Agriculture and Agri-Food Canada (AAFC) as well as other
partners in ARS and within the Long Term Agroecosystems Research (LTAR)
network.



Qualifications include a recent Ph.D. (within the last 4 years) in
Agricultural Engineering, Agronomy, Geography or a related discipline;
ability to participate as a team member in the collection of field data of
soil moisture, biomass, phenology and leaf area index (LAI); and knowledge
of remote sensing and geospatial analysis for using radar and optical data
from earth orbiting satellites. This job is open to U.S. Citizens and
Permanent Residents seeking US Citizenship.



The vacancy announcement (number RA-18-001-HA) can be accessed by visiting:
https://www.usajobs.gov/GetJob/ViewDetails/507599500 .



Please direct any questions about this position to Dr. Alisa Coffin,
[hidden email], 229-386-3665.



Please distribute this announcement or the attached advertisement to
interested contacts and potential candidates.

        [[alternative HTML version deleted]]

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: rgdal compile confusion [RESOLVED]

Mon, 08/20/2018 - 11:19
On Mon, 20 Aug 2018, Shaun Walbridge wrote:

> It looks like the maintainer of the slackbuild recently downgraded GDAL
>  back to the 2.2.4 version because of the problems with the 2.3.0 release:
>  https://urldefense.proofpoint.com/v2/url?u=https-3A__git.slackbuilds.org_slackbuilds_commit_gis_gdal-3Fh-3D14.2-26id-3Dc97bd17ac886c33fb90b3954945e2c3b07907acc&d=DwIGaQ&c=n6-cguzQvX_tUIrZOS_4Og&r=nxHWzsbFn4X_N0ruafqevZE4LSBsSvWC_0wkjR1NVTY&m=ATemW5xbQSNOMLNPRJ7Qhs2z8VgsMWup3KQpGwj65HU&s=V_lu-cvt5T7I-SaxG-5K3zvdqiCDT8vQ-VDgV0r4_Ak&e=

> I'd recommend you reinstall the newest (2.2.4) slackbuild version of the
> package, and try rebuilding rgdal while using that.

Shaun,

   How interesting. I did not see that notice on the SBo mail list and David
did not point me to it when I wrote him about the problem.

   Yes, gdal-2.2.4 supports the building of rgdal-1.3-4. The compatibilty
problem must be subtle as configure finds C++11 on my system yet does not
build the latest gdal or rgdal.

Thanks very much,

Rich

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Re: rgdal compile confusion

Mon, 08/20/2018 - 10:37
It looks like the maintainer of the slackbuild recently downgraded GDAL back to the 2.2.4 version because of the problems with the 2.3.0 release:
  https://urldefense.proofpoint.com/v2/url?u=https-3A__git.slackbuilds.org_slackbuilds_commit_gis_gdal-3Fh-3D14.2-26id-3Dc97bd17ac886c33fb90b3954945e2c3b07907acc&d=DwIGaQ&c=n6-cguzQvX_tUIrZOS_4Og&r=fCPRb7QX-vd5bnO9gIJHCiX852SVUtyYX--xtCKtpfk&m=vscOV_5PwQCvnbNef9uvIzp5QCsmoJ43pg4dGDIoXfM&s=FGqho_qYJ7kdidW3uV6ru6kZrzlMbMBU4ujJC0cq4NY&e=

With the line "The update to version 2.3.0 broke lots of dependent packages that don't use c++11. Sorry!"

I'd recommend you reinstall the newest (2.2.4) slackbuild version of the package, and try rebuilding rgdal while using that.

On 8/20/18, 11:28 AM, "R-sig-Geo on behalf of Rich Shepard" <[hidden email] on behalf of [hidden email]> wrote:

    On Mon, 20 Aug 2018, Roger Bivand wrote:
   
    > And the output of sessionInfo() - we don't know your platform. How was
    > GDAL itself installed: from source or binary? Was that binary built on the
    > same platform? Are there multiple GDAL installations on your system? If
    > so, the configure step may see one version but install sees another.
   
    Roger,
   
    > sessionInfo()
    R version 3.5.1 (2018-07-02)
    Platform: i586-slackware-linux-gnu (32-bit)
    Running under: Slackware 14.2
   
    Matrix products: default
    BLAS: /usr/lib/R/lib/libRblas.so
    LAPACK: /usr/lib/R/lib/libRlapack.so
   
    locale:
      [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C
      [3] LC_TIME=en_US.UTF-8        LC_COLLATE=C
      [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8
      [7] LC_PAPER=en_US.UTF-8       LC_NAME=C
      [9] LC_ADDRESS=C               LC_TELEPHONE=C
    [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
   
    attached base packages:
    [1] stats     graphics  grDevices utils     datasets  methods   base
   
    loaded via a namespace (and not attached):
    [1] compiler_3.5.1 tools_3.5.1    tcltk_3.5.1
   
       gdal is built from source using the SlackBuiilds.org
    <https://urldefense.proofpoint.com/v2/url?u=https-3A__www.slackbuilds.org_&d=DwIGAA&c=n6-cguzQvX_tUIrZOS_4Og&r=YFaRLkcUCdDkLrpTbNOUV9J1CwYBCTMwgm5tdQkRSm4&m=wxTHWsd_1Dg_NP86zMVoYKF9rSG2FjVW_DFJf3-svAU&s=CGzyyiturMPUgBqnwm8a991UlrEj6v7Greh6EUG9L1U&e=> build script.
   
       gdal-2.3.1 is the only gdal package installed.
   
       With the upgraded (2.3.0 -> 2.3.1) gdal version installed fails the same
    way in response to both update.packages() and install.packages('rgdal'):
   
    > install.packages('rgdal')
    trying URL 'https://urldefense.proofpoint.com/v2/url?u=https-3A__ftp.osuosl.org_pub_cran_src_contrib_rgdal-5F1.3-2D4.tar.gz&d=DwIGAA&c=n6-cguzQvX_tUIrZOS_4Og&r=YFaRLkcUCdDkLrpTbNOUV9J1CwYBCTMwgm5tdQkRSm4&m=wxTHWsd_1Dg_NP86zMVoYKF9rSG2FjVW_DFJf3-svAU&s=FgWPmtd9KO4To56icX1JJSQCZ8DGlUHxUCdkHjzY7X0&e='
    Content type 'application/x-gzip' length 1664774 bytes (1.6 MB)
    ==================================================
    downloaded 1.6 MB
   
    * installing *source* package ¡rgdal¢ ...
    ** package ¡rgdal¢ successfully unpacked and MD5 sums checked
    configure: R_HOME: /usr/lib/R
    configure: CC: gcc
    configure: CXX: g++
    configure: C++11 support available
    configure: rgdal: 1.3-4
    checking for /usr/bin/svnversion... yes
    configure: svn revision: 766
    checking for gdal-config... /usr/bin/gdal-config
    checking gdal-config usability... yes
    configure: GDAL: 2.3.1
    checking C++11 support for GDAL >= 2.3.0... yes
    checking GDAL version >= 1.11.4... yes
    checking gdal: linking with --libs only... no
    checking gdal: linking with --libs and --dep-libs... no
    In file included from /usr/include/gdal.h:45:0,
                      from gdal_test.cc:1:
    /usr/include/cpl_port.h:187:6: error: #error Must have C++11 or newer.
      #    error Must have C++11 or newer.
           ^
    In file included from /usr/include/gdal.h:49:0,
                      from gdal_test.cc:1:
    /usr/include/cpl_minixml.h:202:47: error: expected template-name before '<' token
      class CPLXMLTreeCloser: public std::unique_ptr<CPLXMLNode, CPLXMLTreeCloserDeleter>
                                                    ^
    /usr/include/cpl_minixml.h:202:47: error: expected '{' before '<' token
    /usr/include/cpl_minixml.h:202:47: error: expected unqualified-id before '<' token
    In file included from /usr/include/ogr_api.h:45:0,
                      from /usr/include/gdal.h:50,
                      from gdal_test.cc:1:
    /usr/include/ogr_core.h:79:28: error: expected '}' before end of line
    /usr/include/ogr_core.h:79:28: error: expected declaration before end of line
    In file included from /usr/include/gdal.h:45:0,
                      from gdal_test.cc:1:
    /usr/include/cpl_port.h:187:6: error: #error Must have C++11 or newer.
      #    error Must have C++11 or newer.
           ^
    In file included from /usr/include/gdal.h:49:0,
                      from gdal_test.cc:1:
    /usr/include/cpl_minixml.h:202:47: error: expected template-name before '<' token
      class CPLXMLTreeCloser: public std::unique_ptr<CPLXMLNode, CPLXMLTreeCloserDeleter>
                                                    ^
    /usr/include/cpl_minixml.h:202:47: error: expected '{' before '<' token
    /usr/include/cpl_minixml.h:202:47: error: expected unqualified-id before '<' token
    In file included from /usr/include/ogr_api.h:45:0,
                      from /usr/include/gdal.h:50,
                      from gdal_test.cc:1:
    /usr/include/ogr_core.h:79:28: error: expected '}' before end of line
    /usr/include/ogr_core.h:79:28: error: expected declaration before end of line
    configure: Install failure: compilation and/or linkage problems.
    configure: error: GDALAllRegister not found in libgdal.
    ERROR: configuration failed for package ¡rgdal¢
    * removing ¡/usr/lib/R/library/rgdal¢
    * restoring previous ¡/usr/lib/R/library/rgdal¢
   
    The downloaded source packages are in
      ¡/tmp/RtmpfJQHZG/downloaded_packages¢
    Updating HTML index of packages in '.Library'
    Making 'packages.html' ... done
    Warning message:
    In install.packages("rgdal") :
       installation of package ¡rgdal¢ had non-zero exit status
   
       What other information can I provide?
   
    Regards,
   
    Rich
   
    _______________________________________________
    R-sig-Geo mailing list
    [hidden email]
    https://urldefense.proofpoint.com/v2/url?u=https-3A__stat.ethz.ch_mailman_listinfo_r-2Dsig-2Dgeo&d=DwIGAA&c=n6-cguzQvX_tUIrZOS_4Og&r=YFaRLkcUCdDkLrpTbNOUV9J1CwYBCTMwgm5tdQkRSm4&m=wxTHWsd_1Dg_NP86zMVoYKF9rSG2FjVW_DFJf3-svAU&s=jvOPSfwlAz48alb_wKku-UEk-ob5pOUZ5WHG6BIgn1w&e=
   


_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

Pages