[R-sig-Geo] [FORGED] Calculate each polygon percentage inside a circles

ASANTOS alexandresantosbr at yahoo.com.br
Wed May 25 14:23:59 CEST 2016


Dear Rolf Turner,

            It's much better a clean code with a minimum packages, thank 
you very much for your answer. But "pct" object give me a total polygon 
percentage around each point and I need too an identification (in 
columns) of individual contribution of each polygon. In my simulation, I 
find 50.38001% for the point 1, but this is a total percentage of 
polygons and I'd like to know a percentage contribution for each polygon 
(e.g. ID1 = 0.00000 + ID1 = 10.00000 + ID3 = 0.00000 + ID4 = 40.38001 = 
50.38001 total), this is possible?

Thanks again,

Best wishes,

Alexandre


#==================================================================================
library(spatstat)

bdry <- list(cbind(c(180114, 180553, 181127, 181477, 181294,
                      181007, 180409, 180162, 180114),
                    c(332349, 332057, 332342, 333250, 333558,
                      333676, 332618, 332413, 332349)),

              cbind(rev(c(180042, 180545, 180553, 180314, 179955,
                          179142, 179437, 179524, 179979, 180042)),
                    rev(c(332373, 332026, 331426, 330889, 330683,
                          331133, 331623, 332152, 332357, 332373))),

              cbind(rev(c(179110, 179907, 180433, 180712, 180752,
                          180329, 179875, 179668, 179572, 179269,
                          178879, 178600, 178544, 179046, 179110)),
                    rev(c(331086, 330620, 330494, 330265, 330075,
                          330233, 330336, 330004, 329783, 329665,
                          329720, 329933, 330478, 331062, 331086))),

              cbind(c(180304, 180403,179632,179420,180304),
                    c(332791, 333204, 333635, 333058, 332791)))

W <- owin(poly=bdry)

set.seed(42)
X  <- rpoispp(28/area.owin(W),win=W)
plot(X,cols="blue")

for(i in 1:npoints(X)) plot(disc(radius=600,centre=c(X$x[i],X$y[i])),
      add=TRUE,col=rgb(1,0.5,0,alpha=0.4),border=NA)

# To be fair, calculate the area of the discs using area.owin()
# rather than as pi*600^2.

AD <- area.owin(disc(radius=600))

pct <- numeric(npoints(X))
for(i in 1:npoints(X)) {
     Di <- disc(radius=600,centre=c(X$x[i],X$y[i]))
     pct[i] <- 100*area.owin(intersect.owin(Di,W))/AD
}
#==================================================================================

-- 
======================================================================
Alexandre dos Santos
Proteção Florestal
IFMT - Instituto Federal de Educação, Ciência e Tecnologia de Mato Grosso
Campus Cáceres
Caixa Postal 244
Avenida dos Ramires, s/n
Bairro: Distrito Industrial
Cáceres - MT                      CEP: 78.200-000
Fone: (+55) 65 8132-8112 (TIM)   (+55) 65 9686-6970 (VIVO)

         alexandre.santos at cas.ifmt.edu.br
Lattes: http://lattes.cnpq.br/1360403201088680
OrcID: orcid.org/0000-0001-8232-6722
Researchgate: https://www.researchgate.net/profile/Alexandre_Santos10
LinkedIn: https://br.linkedin.com/in/alexandre-dos-santos-87961635
======================================================================

Em 24/05/2016 23:55, Rolf Turner escreveu:
>
> Dear Alexandre,
>
> Your tasks can all be done very simply in spatstat.  The code follows.
> Note that I had to reverse the point order for sr2 and sr3 because
> spatstat requires that the vertices of (exterior) boundaries of 
> polygonal windows be given in anticlockwise order.
>
> I commented out the plotting of the discs centred at each point of the 
> simulated pattern since I found the plot to be unenlightening and 
> messy looking.
>
> The resulting percentages that you are after are in "pct".
>
> If you want more precision for the disc areas, set npoly equal to a 
> large number than the default 128 (e.g.1024 or 2048) in the calls to 
> disc().
>
> cheers,
>
> Rolf Turner
>


	[[alternative HTML version deleted]]



More information about the R-sig-Geo mailing list