subsetting a spatial polygons

classic Classic list List threaded Threaded
4 messages Options
Reply | Threaded
Open this post in threaded view
|

subsetting a spatial polygons

Antonio Silva
Dear list users

I have a SpatialPolygons with several squares. How to subset it to have
only the squares between given latitudes and longitudes?

In the example

library(sp)
library(raster)
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")
plot(polys,axes=T)

How to select only the squares, let's say, between 24-25°S and 45-46°W?

The farthest I went was:

e <- extent(-45.9,-45.1,-24.9,-24.1) # which is not very elegant
mask <- crop(polys,e)
polys2 <- polys[mask,]
plot(polys2,add=T,col="green")

Thanks a lot. Best regards

--
Antônio 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
Reply | Threaded
Open this post in threaded view
|

Re: subsetting a spatial polygons

Tim Appelhans
Antonio,

I am not sure why you think that your solution is not very elegant.
In case you want to have more visual control over the subsetting, you
could try mapedit:

library(mapedit)
myselection = selectFeatures(polys, mode = "draw")

which will let you draw a e.g. rectangle and only return those features
that intersect it.

Best
Tim


On 09/12/2018 04:18 AM, Antonio Silva wrote:

> Dear list users
>
> I have a SpatialPolygons with several squares. How to subset it to have
> only the squares between given latitudes and longitudes?
>
> In the example
>
> library(sp)
> library(raster)
> 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")
> plot(polys,axes=T)
>
> How to select only the squares, let's say, between 24-25°S and 45-46°W?
>
> The farthest I went was:
>
> e <- extent(-45.9,-45.1,-24.9,-24.1) # which is not very elegant
> mask <- crop(polys,e)
> polys2 <- polys[mask,]
> plot(polys2,add=T,col="green")
>
> Thanks a lot. Best regards
>

_______________________________________________
R-sig-Geo mailing list
[hidden email]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
Reply | Threaded
Open this post in threaded view
|

Re: subsetting a spatial polygons

Michael Sumner-2
crop doesn't select it actually cuts on the boundary, for simple shapes I
use coordinates (for centroids in easy matrix) and select with [ using
tests on x and y.

Cheers, Mike

On Wed, Sep 12, 2018, 15:29 Tim Appelhans <[hidden email]> wrote:

> Antonio,
>
> I am not sure why you think that your solution is not very elegant.
> In case you want to have more visual control over the subsetting, you
> could try mapedit:
>
> library(mapedit)
> myselection = selectFeatures(polys, mode = "draw")
>
> which will let you draw a e.g. rectangle and only return those features
> that intersect it.
>
> Best
> Tim
>
>
> On 09/12/2018 04:18 AM, Antonio Silva wrote:
> > Dear list users
> >
> > I have a SpatialPolygons with several squares. How to subset it to have
> > only the squares between given latitudes and longitudes?
> >
> > In the example
> >
> > library(sp)
> > library(raster)
> > 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")
> > plot(polys,axes=T)
> >
> > How to select only the squares, let's say, between 24-25°S and 45-46°W?
> >
> > The farthest I went was:
> >
> > e <- extent(-45.9,-45.1,-24.9,-24.1) # which is not very elegant
> > mask <- crop(polys,e)
> > polys2 <- polys[mask,]
> > plot(polys2,add=T,col="green")
> >
> > Thanks a lot. Best regards
> >
>
> _______________________________________________
> 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
Reply | Threaded
Open this post in threaded view
|

Re: subsetting a spatial polygons

Antonio Silva
Thanks Michael and Tim, thanks for the attention

The first thing I did was I try the subsets using subset.

subset(polys,
coordinates(polys)[,1]>=-46 & coordinates(polys)[,1]<=-45 &
coordinates(polys)[,2]>=-25 & coordinates(polys)[,2]<=-24)

But it seems to work only with Spatial*DataFrame objects. Then I (unluckily)
gave up the x-y approach.

With Mike's advice I tried again as

plot(
    polys[
        coordinates(polys)[,1]>=-46 & coordinates(polys)[,1]<=-45 &
        coordinates(polys)[,2]>=-25 & coordinates(polys)[,2]<=-24,]
,add=T,col="green")

and I got what I wanted in a more "elegant" fashion.

Once more thanks a lot for the attention.

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

Em qua, 12 de set de 2018 às 05:05, Michael Sumner <[hidden email]>
escreveu:

> crop doesn't select it actually cuts on the boundary, for simple shapes I
> use coordinates (for centroids in easy matrix) and select with [ using
> tests on x and y.
>
> Cheers, Mike
>
> On Wed, Sep 12, 2018, 15:29 Tim Appelhans <[hidden email]> wrote:
>
> > Antonio,
> >
> > I am not sure why you think that your solution is not very elegant.
> > In case you want to have more visual control over the subsetting, you
> > could try mapedit:
> >
> > library(mapedit)
> > myselection = selectFeatures(polys, mode = "draw")
> >
> > which will let you draw a e.g. rectangle and only return those features
> > that intersect it.
> >
> > Best
> > Tim
> >
> >
> > On 09/12/2018 04:18 AM, Antonio Silva wrote:
> > > Dear list users
> > >
> > > I have a SpatialPolygons with several squares. How to subset it to have
> > > only the squares between given latitudes and longitudes?
> > >
> > > In the example
> > >
> > > library(sp)
> > > library(raster)
> > > 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")
> > > plot(polys,axes=T)
> > >
> > > How to select only the squares, let's say, between 24-25°S and 45-46°W?
> > >
> > > The farthest I went was:
> > >
> > > e <- extent(-45.9,-45.1,-24.9,-24.1) # which is not very elegant
> > > mask <- crop(polys,e)
> > > polys2 <- polys[mask,]
> > > plot(polys2,add=T,col="green")
> > >
> > > Thanks a lot. Best regards
> > >
> >
> > _______________________________________________
> > 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