From 2f2ebbb7b8a8a03e8e14186bb91128a267edeb8a Mon Sep 17 00:00:00 2001 From: Philippe Verley <philippe.verley@ird.fr> Date: Thu, 29 Mar 2018 14:03:56 +0000 Subject: [PATCH] In src/*.h cosmetic changes to declare the same way all functions. In fads.R there was a wrong call to ripley_tr_disq_ic instead of coor_tr_disq_ic. Deleted comments en plot.vads.R that were problematic according to devtools:check() function. Updated DESCRIPTION, NAMESPACE and MD5 files --- DESCRIPTION | 8 +- MD5 | 122 +-- NAMESPACE | 4 +- R/fads.R | 2510 ++++++++++++++++++++++++------------------------- R/plot.vads.R | 40 +- R/vads.R | 601 ++++++------ src/Zlibs.h | 248 ++--- src/adssub.h | 51 +- 8 files changed, 1734 insertions(+), 1850 deletions(-) diff --git a/DESCRIPTION b/DESCRIPTION index 5e7e7f0..4c11bcf 100755 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -2,10 +2,11 @@ Package: ads Type: Package Title: Spatial point patterns analysis Version: 1.5-3 -Date: 2018-03-28 -Authors@R: c(person("Raphael", "Pelissier", role=c("aut", "cre"), email="raphael.pelissier@ird.fr"), +Date: 2018-03-29 +Authors@R: c(person("Raphael", "Pelissier", role=c("aut"), email="raphael.pelissier@ird.fr"), person("Francois", "Goreau", role="aut"), - person("Philippe", "Verley", role="ctb")) + person("Philippe", "Verley", role=c("ctb", "cre"), email="philippe.verley@ird.fr")) +Author: Raphael Pelissier [aut], Francois Goreau [aut], Philippe Verley [ctb, cre] Maintainer: Raphael Pelissier <raphael.pelissier@ird.fr> Imports: ade4, spatstat Description: Perform first- and second-order multi-scale analyses derived from Ripley K-function, for univariate, @@ -16,3 +17,4 @@ Packaged: 2015-01-13 12:09:18 UTC; root NeedsCompilation: yes Repository: CRAN Date/Publication: 2015-01-13 14:49:23 +RoxygenNote: 6.0.1 diff --git a/MD5 b/MD5 index a20da8f..83b6579 100644 --- a/MD5 +++ b/MD5 @@ -1,55 +1,67 @@ -2755553481e78a8f5b6ef00d43e983a5 *DESCRIPTION -68a58538a504c7cfa73472d87053be45 *INDEX -7478c6d745db573c9d03ebc59bdc969c *NAMESPACE -f6aac0fb4ef187df6521cbb2d9ffc9cf *R/fads.R -7a71372e86fd8aadfc770b261cd953ca *R/mimetic.R -e87bd01e08a83a7b1c52b2a4a7e5df35 *R/plot.fads.R -95d2bc01bd2779736c826a21c83b62f3 *R/plot.vads.R -6317d29d7e9045d55b1a41906e867462 *R/print.fads.R -d0f79442970e383cf8c3bd1677320b53 *R/print.vads.R -c4b216a6fb57acb020f5e21682f7ed13 *R/spp.R -96dada922fee237e190207f544de8ed3 *R/summary.vads.R -c5d76ac2aa5f6fa68c0cfd08f197989e *R/swin.R -34b15c6ed440f88d868d3ab3b93db3b8 *R/triangulate.R -a322bba8d53aff07f9a8ce9a69c89aef *R/util.R -f2c452832a53eb1ed6d168b59918d52a *R/vads.R -31fa89b542936dac1031133b12ef530c *data/Allogny.rda -5450b84c345240671b3410af7c70bc44 *data/BPoirier.rda -24fea786746fc897f8918300f2f2c544 *data/Couepia.rda -a014cac4bc9abd4f40123a762939c8ca *data/Paracou15.rda -674867657e9a3df07b6eced02aa022a1 *data/demopat.rda -c7cd00087730e79e06e8c0d937d0b1f7 *inst/CITATION -94a8b147b3d2730b0c3869a7e1a2592f *man/Allogny.Rd -1755cf31329228337e0b8039af2b6daf *man/BPoirier.Rd -d912451f743e11a9fbe50c8a6ef13b15 *man/Couepia.Rd -4f272a7a2e7528da43c0d2ef8685921b *man/Paracou15.Rd -b3c89a3bcb5db0d5fc9e93d8305df8c8 *man/area.swin.Rd -d7b568c5332aad81bf4fb9f2bd4efed7 *man/demopat.Rd -9b09cc667ab5d522e4a3a77c8b711159 *man/dval.Rd -cd9ec1a82e0ee4e372e7ebd9c11233e4 *man/inside.swin.Rd -e85df8cba02e1d7ae44b447ef05f0e3f *man/internal.Rd -18f7794bc004e9b9599486a725f7da5b *man/k12fun.Rd -26a82fbf3be668f7c74dd4f1b9fd93e4 *man/k12val.Rd -f4c594ec01933acd9d9cae7e797685fe *man/kdfun.Rd -aa725a3e6b7183290f7dd758d594810c *man/kfun.Rd -ea3edb1e20ecd5aa794f48c9fa5ee0e5 *man/kmfun.Rd -0e606a71b6ec6f730bdd0d1498834dbc *man/kp.fun.Rd -500b63ba993de0b1f8a7d7f92765403b *man/kpqfun.Rd -d388b950333aad22313c76d82f51c749 *man/krfun.Rd -9970f6e4573f2e241de12d62ab201e23 *man/ksfun.Rd -d5d294338860db9a406208f219f05f5e *man/kval.Rd -1037e7421f3a50f15fcb07e88dcc260c *man/mimetic.Rd -5bff5bb53402d742900fe12ed542877b *man/plot.fads.Rd -7e88d95ac70e452257a939a70f9a1a76 *man/plot.spp.Rd -c30babaf1b550ec45fa0f1bd89f967e5 *man/plot.vads.Rd -2c61afcc548ad6440240664252a645a5 *man/spp.Rd -cfe40689a11ec983e9b9bb4b759aaa26 *man/swin.Rd -c04839dd6cd8eb396d512260a23dc7fd *man/triangulate.Rd -22fae7567dd07dddf2ea00cb023c2be5 *src/Zlibs.c -109615f7005a5a6e02b98bcf1a904551 *src/Zlibs.h -3d7a3f0a98ac75ad014cc91fbe2d0579 *src/adssub.c -bd642dbad07c62b2f39fae4c5640f93a *src/adssub.h -8a3dad68f1826270eef7ea08e098dfc5 *src/spatstatsub.f -5aa9f5862ef5b0e8bbcc327021e1489a *src/triangulate.c -e482b9d18794f53adaed3f2c98edd19a *src/triangulate.h -fb07aec2cf6396cab654966e2757ab5c *src/util.c +f6aac0fb4ef187df6521cbb2d9ffc9cf R/fads.R +7a71372e86fd8aadfc770b261cd953ca R/mimetic.R +7603ca1f27dbbaf85b7b2ec19e059d7d R/plot.fads.R +154b26509aa63d97c2e0cc646e6681a8 R/plot.vads.R +b85a2aa0125801e94f5b4ce8b0f22a12 R/print.fads.R +1a420e10243f4c0df000ac22e0bab60e R/print.vads.R +c4b216a6fb57acb020f5e21682f7ed13 R/spp.R +bbb6ddd2e55fabe14c9e9dfcdfc495e7 R/summary.vads.R +c5d76ac2aa5f6fa68c0cfd08f197989e R/swin.R +318472a4c160e2ac070b5fb390add04e R/triangulate.R +a322bba8d53aff07f9a8ce9a69c89aef R/util.R +9f002a8b5ed00aaf26b4d7dc2f914f40 R/vads.R +31fa89b542936dac1031133b12ef530c data/Allogny.rda +5450b84c345240671b3410af7c70bc44 data/BPoirier.rda +24fea786746fc897f8918300f2f2c544 data/Couepia.rda +674867657e9a3df07b6eced02aa022a1 data/demopat.rda +a014cac4bc9abd4f40123a762939c8ca data/Paracou15.rda +94a8b147b3d2730b0c3869a7e1a2592f man/Allogny.Rd +b3c89a3bcb5db0d5fc9e93d8305df8c8 man/area.swin.Rd +1755cf31329228337e0b8039af2b6daf man/BPoirier.Rd +d912451f743e11a9fbe50c8a6ef13b15 man/Couepia.Rd +d7b568c5332aad81bf4fb9f2bd4efed7 man/demopat.Rd +9b09cc667ab5d522e4a3a77c8b711159 man/dval.Rd +cd9ec1a82e0ee4e372e7ebd9c11233e4 man/inside.swin.Rd +e85df8cba02e1d7ae44b447ef05f0e3f man/internal.Rd +18f7794bc004e9b9599486a725f7da5b man/k12fun.Rd +26a82fbf3be668f7c74dd4f1b9fd93e4 man/k12val.Rd +f4c594ec01933acd9d9cae7e797685fe man/kdfun.Rd +aa725a3e6b7183290f7dd758d594810c man/kfun.Rd +ea3edb1e20ecd5aa794f48c9fa5ee0e5 man/kmfun.Rd +0e606a71b6ec6f730bdd0d1498834dbc man/kp.fun.Rd +b45ac5f58b213884c5b1d9c79d5d1599 man/kpfun.Rd +500b63ba993de0b1f8a7d7f92765403b man/kpqfun.Rd +d388b950333aad22313c76d82f51c749 man/krfun.Rd +9970f6e4573f2e241de12d62ab201e23 man/ksfun.Rd +d5d294338860db9a406208f219f05f5e man/kval.Rd +1037e7421f3a50f15fcb07e88dcc260c man/mimetic.Rd +4f272a7a2e7528da43c0d2ef8685921b man/Paracou15.Rd +5bff5bb53402d742900fe12ed542877b man/plot.fads.Rd +7e88d95ac70e452257a939a70f9a1a76 man/plot.spp.Rd +c30babaf1b550ec45fa0f1bd89f967e5 man/plot.vads.Rd +2c61afcc548ad6440240664252a645a5 man/spp.Rd +cfe40689a11ec983e9b9bb4b759aaa26 man/swin.Rd +c04839dd6cd8eb396d512260a23dc7fd man/triangulate.Rd +f6aac0fb4ef187df6521cbb2d9ffc9cf R/fads.R +7a71372e86fd8aadfc770b261cd953ca R/mimetic.R +7603ca1f27dbbaf85b7b2ec19e059d7d R/plot.fads.R +154b26509aa63d97c2e0cc646e6681a8 R/plot.vads.R +b85a2aa0125801e94f5b4ce8b0f22a12 R/print.fads.R +1a420e10243f4c0df000ac22e0bab60e R/print.vads.R +c4b216a6fb57acb020f5e21682f7ed13 R/spp.R +bbb6ddd2e55fabe14c9e9dfcdfc495e7 R/summary.vads.R +c5d76ac2aa5f6fa68c0cfd08f197989e R/swin.R +318472a4c160e2ac070b5fb390add04e R/triangulate.R +a322bba8d53aff07f9a8ce9a69c89aef R/util.R +9f002a8b5ed00aaf26b4d7dc2f914f40 R/vads.R +3d7a3f0a98ac75ad014cc91fbe2d0579 src/adssub.c +bd642dbad07c62b2f39fae4c5640f93a src/adssub.h +8a3dad68f1826270eef7ea08e098dfc5 src/spatstatsub.f +5aa9f5862ef5b0e8bbcc327021e1489a src/triangulate.c +e482b9d18794f53adaed3f2c98edd19a src/triangulate.h +fb07aec2cf6396cab654966e2757ab5c src/util.c +f44e3b71761cdf40a551c450517d20cd src/Zlibs.c +109615f7005a5a6e02b98bcf1a904551 src/Zlibs.h +d9ae8cc82c07551b180991f6d508aedf NAMESPACE +a32d955da46938603fc72e218ce86101 DESCRIPTION +68a58538a504c7cfa73472d87053be45 INDEX diff --git a/NAMESPACE b/NAMESPACE index c225fe8..1c9cfcd 100755 --- a/NAMESPACE +++ b/NAMESPACE @@ -4,8 +4,8 @@ importFrom(spatstat,"border","bounding.box.xy","area.owin") importFrom("graphics", "abline", "barplot", "layout", "lines", "par", "plot", "plot.default", "points", "polygon", "symbols", "text") - importFrom("stats", "as.dist", "var") - importFrom("utils", "read.table", "str") +importFrom("stats", "as.dist", "var") +importFrom("utils", "read.table", "str") # Export all names (should be improved in the future) export( diff --git a/R/fads.R b/R/fads.R index dd947c2..90e6f3e 100755 --- a/R/fads.R +++ b/R/fads.R @@ -1,1255 +1,1255 @@ -kfun<-function(p,upto,by,nsim=0,prec=0.01,alpha=0.01) { - # checking for input parameters - stopifnot(inherits(p,"spp")) - if(p$type!="univariate") - warning(paste(p$type,"point pattern has been considered to be univariate\n")) - stopifnot(is.numeric(upto)) - stopifnot(upto>=1) - stopifnot(is.numeric(by)) - stopifnot(by>0) - r<-seq(by,upto,by) - tmax<-length(r) - stopifnot(is.numeric(nsim)) - stopifnot(nsim>=0) - nsim<-testInteger(nsim) - stopifnot(is.numeric(prec)) - stopifnot(prec>=0) - stopifnot(is.numeric(alpha)) - stopifnot(alpha>=0) - if(nsim>0) testIC(nsim,alpha) - - if("rectangle"%in%p$window$type) { - cas<-1 - xmin<-p$window$xmin - xmax<-p$window$xmax - ymin<-p$window$ymin - ymax<-p$window$ymax - stopifnot(upto<=(0.5*max((xmax-xmin),(ymax-ymin)))) - if ("complex"%in%p$window$type) { - cas<-3 - tri<-p$window$triangles - nbTri<-nrow(tri) - } - } - else if("circle"%in%p$window$type) { - cas<-2 - x0<-p$window$x0 - y0<-p$window$y0 - r0<-p$window$r0 - stopifnot(upto<=r0) - if ("complex"%in%p$window$type) { - cas<-4 - tri<-p$window$triangles - nbTri<-nrow(tri) - } - } - else - stop("invalid window type") - intensity<-p$n/area.swin(p$window) - - if(cas==1) { #rectangle - if(nsim==0) { #without CI - res<-.C("ripley_rect", - as.integer(p$n),as.double(p$x),as.double(p$y), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(tmax),as.double(by), - g=double(tmax),k=double(tmax), - PACKAGE="ads") - } - else { #with CI - res<-.C("ripley_rect_ic", - as.integer(p$n),as.double(p$x),as.double(p$y), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax),as.double(intensity), - as.integer(tmax),as.double(by), - as.integer(nsim),as.double(prec),as.double(alpha), - g=double(tmax),k=double(tmax), - gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), - gval=double(tmax),kval=double(tmax),lval=double(tmax),nval=double(tmax), - PACKAGE="ads") - } - } - else if(cas==2) { #circle - if(nsim==0) { #without CI - res<-.C("ripley_disq", - as.integer(p$n),as.double(p$x),as.double(p$y), - as.double(x0),as.double(y0),as.double(r0), - as.integer(tmax),as.double(by), - g=double(tmax),k=double(tmax), - PACKAGE="ads") - } - else { #with CI - res<-.C("ripley_disq_ic", - as.integer(p$n),as.double(p$x),as.double(p$y), - as.double(x0),as.double(y0),as.double(r0),as.double(intensity), - as.integer(tmax),as.double(by), - as.integer(nsim),as.double(prec),as.double(alpha), - g=double(tmax),k=double(tmax), - gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), - gval=double(tmax),kval=double(tmax),lval=double(tmax),nval=double(tmax), - PACKAGE="ads") - } - } - else if(cas==3) { #complex within rectangle - if(nsim==0) { #without CI - res<-.C("ripley_tr_rect", - as.integer(p$n),as.double(p$x),as.double(p$y), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), - g=double(tmax),k=double(tmax), - PACKAGE="ads") - } - else { #with CI - res<-.C("ripley_tr_rect_ic", - as.integer(p$n),as.double(p$x),as.double(p$y), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax),as.double(intensity), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), - as.integer(nsim),as.double(prec),as.double(alpha), - g=double(tmax),k=double(tmax), - gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), - gval=double(tmax),kval=double(tmax),lval=double(tmax),nval=double(tmax), - PACKAGE="ads") - } - } - else if(cas==4) { #complex within circle - if(nsim==0) { #without CI - res<-.C("ripley_tr_disq", - as.integer(p$n),as.double(p$x),as.double(p$y), - as.double(x0),as.double(y0),as.double(r0), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), - g=double(tmax),k=double(tmax), - PACKAGE="ads") - } - else { #with CI - res<-.C("ripley_tr_disq_ic", - as.integer(p$n),as.double(p$x),as.double(p$y), - as.double(x0),as.double(y0),as.double(r0),as.double(intensity), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), - as.integer(nsim),as.double(prec),as.double(alpha), - g=double(tmax),k=double(tmax), - gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), - gval=double(tmax),kval=double(tmax),lval=double(tmax),nval=double(tmax), - PACKAGE="ads") - } - } - # formatting results - ds<-c(pi,diff(pi*r^2)) - g<-data.frame(obs=res$g/(intensity*ds),theo=rep(1,tmax)) - n<-data.frame(obs=res$k/(pi*r^2),theo=rep(intensity,tmax)) - k<-data.frame(obs=res$k/intensity,theo=pi*r^2) - l<-data.frame(obs=sqrt(res$k/(intensity*pi))-r,theo=rep(0,tmax)) - if(nsim>0) { - g<-cbind(g,sup=res$gic1/(intensity*ds),inf=res$gic2/(intensity*ds),pval=res$gval/(nsim+1)) - n<-cbind(n,sup=res$kic1/(pi*r^2),inf=res$kic2/(pi*r^2),pval=res$nval/(nsim+1)) - k<-cbind(k,sup=res$kic1/intensity,inf=res$kic2/intensity,pval=res$kval/(nsim+1)) - l<-cbind(l,sup=sqrt(res$kic1/(intensity*pi))-r,inf=sqrt(res$kic2/(intensity*pi))-r,pval=res$lval/(nsim+1)) - } - call<-match.call() - res<-list(call=call,r=r,g=g,n=n,k=k,l=l) - class(res)<-c("fads","kfun") - return(res) -} - -k12fun<-function(p,upto,by,nsim=0,H0=c("pitor","pimim","rl"),prec=0.01,nsimax=3000,conv=50,rep=10,alpha=0.01,marks) { - # checking for input parameters - options( CBoundsCheck = TRUE ) - # regle les problemes pour 32-bit - stopifnot(inherits(p,"spp")) - stopifnot(p$type=="multivariate") - stopifnot(is.numeric(upto)) - stopifnot(upto>=1) - stopifnot(is.numeric(by)) - stopifnot(by>0) - r<-seq(by,upto,by) - tmax<-length(r) - stopifnot(is.numeric(nsim)) - stopifnot(nsim>=0) - nsim<-testInteger(nsim) - H0<-H0[1] - stopifnot(H0=="pitor" || H0=="pimim" || H0=="rl") - if(H0=="rl") H0<-1 - else if(H0=="pitor") H0<-2 - else H0<-3 - stopifnot(is.numeric(prec)) - stopifnot(prec>=0) - stopifnot(is.numeric(alpha)) - stopifnot(alpha>=0) - if(nsim>0) testIC(nsim,alpha) - if(missing(marks)) - marks<-c(1,2) - stopifnot(length(marks)==2) - stopifnot(marks[1]!=marks[2]) - mark1<-marks[1] - mark2<-marks[2] - if(is.numeric(mark1)) - mark1<-levels(p$marks)[testInteger(mark1)] - else if(!mark1%in%p$marks) stop(paste("mark \'",mark1,"\' doesn\'t exist",sep="")) - if(is.numeric(mark2)) - mark2<-levels(p$marks)[testInteger(mark2)] - else if(!mark2%in%p$marks) stop(paste("mark \'",mark2,"\' doesn\'t exist",sep="")) - - if("rectangle"%in%p$window$type) { - cas<-1 - xmin<-p$window$xmin - xmax<-p$window$xmax - ymin<-p$window$ymin - ymax<-p$window$ymax - stopifnot(upto<=(0.5*max((xmax-xmin),(ymax-ymin)))) - if ("complex"%in%p$window$type) { - cas<-3 - tri<-p$window$triangles - nbTri<-nrow(tri) - } - } - else if("circle"%in%p$window$type) { - cas<-2 - x0<-p$window$x0 - y0<-p$window$y0 - r0<-p$window$r0 - stopifnot(upto<=r0) - if ("complex"%in%p$window$type) { - cas<-4 - tri<-p$window$triangles - nbTri<-nrow(tri) - } - } - else - stop("invalid window type") - surface<-area.swin(p$window) - x1<-p$x[p$marks==mark1] - y1<-p$y[p$marks==mark1] - x2<-p$x[p$marks==mark2] - y2<-p$y[p$marks==mark2] - nbPts1<-length(x1) - nbPts2<-length(x2) - intensity2<-nbPts2/surface -# intensity<-(nbPts1+nbPts2)/surface - - # computing intertype functions - if(cas==1) { #rectangle - if(nsim==0) { #without CI - res<-.C("intertype_rect", - as.integer(nbPts1),as.double(x1),as.double(y1), - as.integer(nbPts2),as.double(x2),as.double(y2), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(tmax),as.double(by), - g=double(tmax),k=double(tmax), - PACKAGE="ads") - } - else { #with CI - res<-.C("intertype_rect_ic", - as.integer(nbPts1),as.double(x1),as.double(y1), - as.integer(nbPts2),as.double(x2),as.double(y2), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax),as.double(surface), - as.integer(tmax),as.double(by), - as.integer(nsim),as.integer(H0),as.double(prec),as.integer(nsimax),as.integer(conv),as.integer(rep),as.double(alpha), - g=double(tmax),k=double(tmax), - gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), - gval=double(tmax),kval=double(tmax),lval=double(tmax),nval=double(tmax), - PACKAGE="ads") - } - } - else if(cas==2) { #circle - if(nsim==0) { #without CI - res<-.C("intertype_disq", - as.integer(nbPts1),as.double(x1),as.double(y1), - as.integer(nbPts2),as.double(x2),as.double(y2), - as.double(x0),as.double(y0),as.double(r0), - as.integer(tmax),as.double(by), - g=double(tmax),k=double(tmax), - PACKAGE="ads") - } - else { #with CI - res<-.C("intertype_disq_ic", - as.integer(nbPts1),as.double(x1),as.double(y1), - as.integer(nbPts2),as.double(x2),as.double(y2), - as.double(x0),as.double(y0),as.double(r0),as.double(surface), - as.integer(tmax),as.double(by), - as.integer(nsim),as.integer(H0),as.double(prec),as.integer(nsimax),as.integer(conv),as.integer(rep),as.double(alpha), - g=double(tmax),k=double(tmax), - gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), - gval=double(tmax),kval=double(tmax),lval=double(tmax),nval=double(tmax), - PACKAGE="ads") - } - } - else if(cas==3) { #complex within rectangle - if(nsim==0) { #without CI - res<-.C("intertype_tr_rect", - as.integer(nbPts1),as.double(x1),as.double(y1), - as.integer(nbPts2),as.double(x2),as.double(y2), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), - g=double(tmax),k=double(tmax), - PACKAGE="ads") - } - else { #with CI - res<-.C("intertype_tr_rect_ic", - as.integer(nbPts1),as.double(x1),as.double(y1), - as.integer(nbPts2),as.double(x2),as.double(y2), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax),as.double(surface), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), - as.integer(nsim),as.integer(H0),as.double(prec),as.integer(nsimax),as.integer(conv),as.integer(rep),as.double(alpha), - g=double(tmax),k=double(tmax), - gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), - gval=double(tmax),kval=double(tmax),lval=double(tmax),nval=double(tmax), - PACKAGE="ads") - } - } - else if(cas==4) { #complex within circle - if(nsim==0) { #without CI - res<-.C("intertype_tr_disq", - as.integer(nbPts1),as.double(x1),as.double(y1), - as.integer(nbPts2),as.double(x2),as.double(y2), - as.double(x0),as.double(y0),as.double(r0), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), - g=double(tmax),k=double(tmax), - PACKAGE="ads") - } - else { #with CI - res<-.C("intertype_tr_disq_ic", - as.integer(nbPts1),as.double(x1),as.double(y1), - as.integer(nbPts2),as.double(x2),as.double(y2), - as.double(x0),as.double(y0),as.double(r0),as.double(surface), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), - as.integer(nsim),as.integer(H0),as.double(prec),as.integer(nsimax),as.integer(conv),as.integer(rep),as.double(alpha), - g=double(tmax),k=double(tmax), - gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), - gval=double(tmax),kval=double(tmax),lval=double(tmax),nval=double(tmax), - PACKAGE="ads") - } - } - # formatting results - ds<-c(pi,diff(pi*r^2)) - g<-res$g/(intensity2*ds) - n<-res$k/(pi*r^2) - k<-res$k/intensity2 - l<-sqrt(res$k/(intensity2*pi))-r - if(H0==1) { - rip<-kfun(spp(c(x1,x2),c(y1,y2),p$window),upto,by) - theo<-list(g=rip$g$obs,n=intensity2*rip$k$obs/(pi*r^2),k=rip$k$obs,l=rip$l$obs) - } - else if (H0==2||H0==3) - theo<-list(g=rep(1,tmax),n=rep(intensity2,tmax),k=pi*r^2,l=rep(0,tmax)) - g<-data.frame(obs=g,theo=theo$g) - n<-data.frame(obs=n,theo=theo$n) - k<-data.frame(obs=k,theo=theo$k) - l<-data.frame(obs=l,theo=theo$l) - if(nsim>0) { - g<-cbind(g,sup=res$gic1/(intensity2*ds),inf=res$gic2/(intensity2*ds),pval=res$gval/(nsim+1)) - n<-cbind(n,sup=res$kic1/(pi*r^2),inf=res$kic2/(pi*r^2),pval=res$nval/(nsim+1)) - k<-cbind(k,sup=res$kic1/intensity2,inf=res$kic2/intensity2,pval=res$kval/(nsim+1)) - l<-cbind(l,sup=sqrt(res$kic1/(intensity2*pi))-r,inf=sqrt(res$kic2/(intensity2*pi))-r,pval=res$lval/(nsim+1)) - } - call<-match.call() - res<-list(call=call,r=r,g12=g,n12=n,k12=k,l12=l,marks=c(mark1,mark2)) - class(res)<-c("fads","k12fun") - return(res) -} - -kijfun<-kpqfun<-function(p,upto,by) { -# checking for input parameters - stopifnot(inherits(p,"spp")) - stopifnot(p$type=="multivariate") - stopifnot(is.numeric(upto)) - stopifnot(upto>=1) - stopifnot(is.numeric(by)) - stopifnot(by>0) - r<-seq(by,upto,by) - tmax<-length(r) - - if("rectangle"%in%p$window$type) { - cas<-1 - xmin<-p$window$xmin - xmax<-p$window$xmax - ymin<-p$window$ymin - ymax<-p$window$ymax - stopifnot(upto<=(0.5*max((xmax-xmin),(ymax-ymin)))) - if ("complex"%in%p$window$type) { - cas<-3 - tri<-p$window$triangles - nbTri<-nrow(tri) - } - } - else if("circle"%in%p$window$type) { - cas<-2 - x0<-p$window$x0 - y0<-p$window$y0 - r0<-p$window$r0 - stopifnot(upto<=r0) - if ("complex"%in%p$window$type) { - cas<-4 - tri<-p$window$triangles - nbTri<-nrow(tri) - } - } - else - stop("invalid window type") - surface<-area.swin(p$window) - - tabMarks<-levels(p$marks) - nbMarks<-length(tabMarks) - mpt_nb<-summary(p$marks) -# computing RipleyClass - gij<-double(tmax*nbMarks^2) - kij<-double(tmax*nbMarks^2) - lij<-double(tmax*nbMarks^2) - nij<-double(tmax*nbMarks^2) - for(i in 1:nbMarks) { - x1<-p$x[p$marks==tabMarks[i]] - y1<-p$y[p$marks==tabMarks[i]] - if(cas==1) { #rectangle - res<-.C("ripley_rect", - as.integer(mpt_nb[i]),as.double(x1),as.double(y1), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(tmax),as.double(by), - g=double(tmax),k=double(tmax), - PACKAGE="ads") - } - else if(cas==2) { #circle - res<-.C("ripley_disq", - as.integer(mpt_nb[i]),as.double(x1),as.double(y1), - as.double(x0),as.double(y0),as.double(r0), - as.integer(tmax),as.double(by), - g=double(tmax),k=double(tmax), - PACKAGE="ads") - } - else if(cas==3) { #complex within rectangle - res<-.C("ripley_tr_rect", - as.integer(mpt_nb[i]),as.double(x1),as.double(y1), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), - g=double(tmax),k=double(tmax), - PACKAGE="ads") - } - else if(cas==4) { #complex within circle - res<-.C("ripley_tr_disq", - as.integer(mpt_nb[i]),as.double(x1),as.double(y1), - as.double(x0),as.double(y0),as.double(r0), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), - g=double(tmax),k=double(tmax), - PACKAGE="ads") - } - intensity1<-mpt_nb[i]/surface - matcol<-(i-1)*nbMarks+i-1 - j<-(matcol*tmax+1):(matcol*tmax+tmax) - ds<-c(pi,diff(pi*r^2)) - gij[j]<-res$g/(intensity1*ds) - nij[j]<-res$k/(pi*r^2) - kij[j]<-res$k/intensity1 - lij[j]<-sqrt(res$k/(intensity1*pi))-r - if(i<nbMarks) { - for(j in (i+1):nbMarks) { - x2<-p$x[p$marks==tabMarks[j]] - y2<-p$y[p$marks==tabMarks[j]] - if(cas==1) { #rectangle - res<-.C("intertype_rect", - as.integer(mpt_nb[i]),as.double(x1),as.double(y1), - as.integer(mpt_nb[j]),as.double(x2),as.double(y2), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(tmax),as.double(by), - g=double(tmax),k=double(tmax), - PACKAGE="ads") - } - else if(cas==2) { #circle - res<-.C("intertype_disq", - as.integer(mpt_nb[i]),as.double(x1),as.double(y1), - as.integer(mpt_nb[j]),as.double(x2),as.double(y2), - as.double(x0),as.double(y0),as.double(r0), - as.integer(tmax),as.double(by), - g=double(tmax),k=double(tmax), - PACKAGE="ads") - } - else if(cas==3) { #complex within rectangle - res<-.C("intertype_tr_rect", - as.integer(mpt_nb[i]),as.double(x1),as.double(y1), - as.integer(mpt_nb[j]),as.double(x2),as.double(y2), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), - g=double(tmax),k=double(tmax), - PACKAGE="ads") - } - else if(cas==4) { #complex within circle - res<-.C("intertype_tr_disq", - as.integer(mpt_nb[i]),as.double(x1),as.double(y1), - as.integer(mpt_nb[j]),as.double(x2),as.double(y2), - as.double(x0),as.double(y0),as.double(r0), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), - g=double(tmax),k=double(tmax), - PACKAGE="ads") - } - intensity2<-mpt_nb[j]/surface - matcol<-(i-1)*nbMarks+j-1 - k<-(matcol*tmax+1):(matcol*tmax+tmax) - gij[k]<-res$g/(intensity2*ds) - nij[k]<-res$k/(pi*r^2) - kij[k]<-res$k/intensity2 - lij[k]<-sqrt(res$k/(intensity2*pi))-r - if(cas==1) { #rectangle - res<-.C("intertype_rect", - as.integer(mpt_nb[j]),as.double(x2),as.double(y2), - as.integer(mpt_nb[i]),as.double(x1),as.double(y1), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(tmax),as.double(by), - g=double(tmax),k=double(tmax), - PACKAGE="ads") - } - else if(cas==2) { #circle - res<-.C("intertype_disq", - as.integer(mpt_nb[j]),as.double(x2),as.double(y2), - as.integer(mpt_nb[i]),as.double(x1),as.double(y1), - as.double(x0),as.double(y0),as.double(r0), - as.integer(tmax),as.double(by), - g=double(tmax),k=double(tmax), - PACKAGE="ads") - } - else if(cas==3) { #complex within rectangle - res<-.C("intertype_tr_rect", - as.integer(mpt_nb[j]),as.double(x2),as.double(y2), - as.integer(mpt_nb[i]),as.double(x1),as.double(y1), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), - g=double(tmax),k=double(tmax), - PACKAGE="ads") - } - else if(cas==4) { #complex within circle - res<-.C("intertype_tr_disq", - as.integer(mpt_nb[j]),as.double(x2),as.double(y2), - as.integer(mpt_nb[i]),as.double(x1),as.double(y1), - as.double(x0),as.double(y0),as.double(r0), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), - g=double(tmax),k=double(tmax), - PACKAGE="ads") - } - matcol<-(j-1)*nbMarks+i-1 - k<-(matcol*tmax+1):(matcol*tmax+tmax) - gij[k]<-res$g/(intensity1*ds) - nij[k]<-res$k/(pi*r^2) - kij[k]<-res$k/intensity1 - lij[k]<-sqrt(res$k/(intensity1*pi))-r - } - } - } - labij<-paste(rep(tabMarks,each=nbMarks),rep(tabMarks,nbMarks),sep="-") - gij<-matrix(gij,nrow=tmax,ncol=nbMarks^2) - kij<-matrix(kij,nrow=tmax,ncol=nbMarks^2) - nij<-matrix(nij,nrow=tmax,ncol=nbMarks^2) - lij<-matrix(lij,nrow=tmax,ncol=nbMarks^2) - call<-match.call() - res<-list(call=call,r=r,labpq=labij,gij=gij,kpq=kij,lpq=lij,npq=nij,intensity=summary(p)$intensity) - class(res)<-c("fads","kpqfun") - return(res) -} - -ki.fun<-kp.fun<-function(p,upto,by) { - # checking for input parameters - stopifnot(inherits(p,"spp")) - stopifnot(p$type=="multivariate") - stopifnot(is.numeric(upto)) - stopifnot(upto>=1) - stopifnot(is.numeric(by)) - stopifnot(by>0) - r<-seq(by,upto,by) - tmax<-length(r) - - if("rectangle"%in%p$window$type) { - cas<-1 - xmin<-p$window$xmin - xmax<-p$window$xmax - ymin<-p$window$ymin - ymax<-p$window$ymax - stopifnot(upto<=(0.5*max((xmax-xmin),(ymax-ymin)))) - if ("complex"%in%p$window$type) { - cas<-3 - tri<-p$window$triangles - nbTri<-nrow(tri) - } - } - else if("circle"%in%p$window$type) { - cas<-2 - x0<-p$window$x0 - y0<-p$window$y0 - r0<-p$window$r0 - stopifnot(upto<=r0) - if ("complex"%in%p$window$type) { - cas<-4 - tri<-p$window$triangles - nbTri<-nrow(tri) - } - } - else - stop("invalid window type") - surface<-area.swin(p$window) - - tabMarks<-levels(p$marks) - nbMarks<-length(tabMarks) - mpt_nb<-summary(p$marks) - #computing RipleyAll - gis<-double(tmax*nbMarks) - kis<-double(tmax*nbMarks) - lis<-double(tmax*nbMarks) - nis<-double(tmax*nbMarks) - for(i in 1:nbMarks) { - x1<-p$x[p$marks==tabMarks[i]] - y1<-p$y[p$marks==tabMarks[i]] - x2<-p$x[p$marks!=tabMarks[i]] - y2<-p$y[p$marks!=tabMarks[i]] - nbPts1<-mpt_nb[i] - nbPts2<-sum(mpt_nb)-nbPts1 - if(cas==1) { #rectangle - res<-.C("intertype_rect", - as.integer(nbPts1),as.double(x1),as.double(y1), - as.integer(nbPts2),as.double(x2),as.double(y2), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(tmax),as.double(by), - g=double(tmax),k=double(tmax), - PACKAGE="ads") - } - else if(cas==2) { #circle - res<-.C("intertype_disq", - as.integer(nbPts1),as.double(x1),as.double(y1), - as.integer(nbPts2),as.double(x2),as.double(y2), - as.double(x0),as.double(y0),as.double(r0), - as.integer(tmax),as.double(by), - g=double(tmax),k=double(tmax), - PACKAGE="ads") - } - else if(cas==3) { #complex within rectangle - res<-.C("intertype_tr_rect", - as.integer(nbPts1),as.double(x1),as.double(y1), - as.integer(nbPts2),as.double(x2),as.double(y2), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), - g=double(tmax),k=double(tmax), - PACKAGE="ads") - } - else if(cas==4) { #complex within circle - res<-.C("intertype_tr_disq", - as.integer(nbPts1),as.double(x1),as.double(y1), - as.integer(nbPts2),as.double(x2),as.double(y2), - as.double(x0),as.double(y0),as.double(r0), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), - g=double(tmax),k=double(tmax), - PACKAGE="ads") - } - intensity2<-nbPts2/surface - j<-((i-1)*tmax+1):((i-1)*tmax+tmax) - ds<-c(pi,diff(pi*r^2)) - gis[j]<-res$g/(intensity2*ds) - nis[j]<-res$k/(pi*r^2) - kis[j]<-res$k/intensity2 - lis[j]<-sqrt(res$k/(intensity2*pi))-r - } - # formatting results - labi<-tabMarks - gis<-matrix(gis,nrow=tmax,ncol=nbMarks) - kis<-matrix(kis,nrow=tmax,ncol=nbMarks) - nis<-matrix(nis,nrow=tmax,ncol=nbMarks) - lis<-matrix(lis,nrow=tmax,ncol=nbMarks) - call<-match.call() - res<-list(call=call,r=r,labp=labi,gp.=gis,kp.=kis,lp.=lis,np.=nis,intensity=summary(p)$intensity) - class(res)<-c("fads","kp.fun") - return(res) -} - -kmfun<-function(p,upto,by,nsim=0,alpha=0.01) { - # checking for input parameters - stopifnot(inherits(p,"spp")) - stopifnot(p$type=="marked") - stopifnot(is.numeric(upto)) - stopifnot(upto>=1) - stopifnot(is.numeric(by)) - stopifnot(by>0) - r<-seq(by,upto,by) - tmax<-length(r) - stopifnot(is.numeric(nsim)) - stopifnot(nsim>=0) - nsim<-testInteger(nsim) - #cmoy<-mean(p$marks) - cvar<-var(p$marks) - stopifnot(is.numeric(alpha)) - stopifnot(alpha>=0) - if(nsim>0) testIC(nsim,alpha) - - if("rectangle"%in%p$window$type) { - cas<-1 - xmin<-p$window$xmin - xmax<-p$window$xmax - ymin<-p$window$ymin - ymax<-p$window$ymax - stopifnot(upto<=(0.5*max((xmax-xmin),(ymax-ymin)))) - if ("complex"%in%p$window$type) { - cas<-3 - tri<-p$window$triangles - nbTri<-nrow(tri) - } - } - else if("circle"%in%p$window$type) { - cas<-2 - x0<-p$window$x0 - y0<-p$window$y0 - r0<-p$window$r0 - stopifnot(upto<=r0) - if ("complex"%in%p$window$type) { - cas<-4 - tri<-p$window$triangles - nbTri<-nrow(tri) - } - } - else - stop("invalid window type") - intensity<-p$n/area.swin(p$window) - - if(cas==1) { #rectangle - if(nsim==0) { #without CI - res<-.C("corr_rect", - as.integer(p$n),as.double(p$x),as.double(p$y),as.double(p$marks), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(tmax),as.double(by), - gm=double(tmax),km=double(tmax), - PACKAGE="ads") - } - else { #with CI - res<-.C("corr_rect_ic", - as.integer(p$n),as.double(p$x),as.double(p$y),as.double(p$marks), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(tmax),as.double(by), - as.integer(nsim),as.double(alpha), - gm=double(tmax),km=double(tmax), - gmic1=double(tmax),gmic2=double(tmax),kmic1=double(tmax),kmic2=double(tmax), - gmval=double(tmax),kmval=double(tmax), - PACKAGE="ads") - } - } - else if(cas==2) { #circle - if(nsim==0) { #without CI - res<-.C("corr_disq", - as.integer(p$n),as.double(p$x),as.double(p$y),as.double(p$marks), - as.double(x0),as.double(y0),as.double(r0), - as.integer(tmax),as.double(by), - gm=double(tmax),km=double(tmax), - PACKAGE="ads") - } - else { #with CI - res<-.C("corr_disq_ic", - as.integer(p$n),as.double(p$x),as.double(p$y),as.double(p$marks), - as.double(x0),as.double(y0),as.double(r0), - as.integer(tmax),as.double(by), - as.integer(nsim),as.double(alpha), - gm=double(tmax),km=double(tmax), - gmic1=double(tmax),gmic2=double(tmax),kmic1=double(tmax),kmic2=double(tmax), - gmval=double(tmax),kmval=double(tmax), - PACKAGE="ads") - } - } - else if(cas==3) { #complex within rectangle - if(nsim==0) { #without CI - res<-.C("corr_tr_rect", - as.integer(p$n),as.double(p$x),as.double(p$y),as.double(p$marks), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), - gm=double(tmax),km=double(tmax), - PACKAGE="ads") - } - else { #with CI - res<-.C("corr_tr_rect_ic", - as.integer(p$n),as.double(p$x),as.double(p$y),as.double(p$marks), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), - as.integer(nsim),as.double(alpha), - gm=double(tmax),km=double(tmax), - gmic1=double(tmax),gmic2=double(tmax),kmic1=double(tmax),kmic2=double(tmax), - gmval=double(tmax),kmval=double(tmax), - PACKAGE="ads") - } - } - else if(cas==4) { #complex within circle - if(nsim==0) { #without CI - res<-.C("corr_tr_disq", - as.integer(p$n),as.double(p$x),as.double(p$y),as.double(p$marks), - as.double(x0),as.double(y0),as.double(r0), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), - gm=double(tmax),km=double(tmax), - PACKAGE="ads") - } - else { #with CI - res<-.C("ripley_tr_disq_ic", - as.integer(p$n),as.double(p$x),as.double(p$y),as.double(p$marks), - as.double(x0),as.double(y0),as.double(r0), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), - as.integer(nsim),as.double(alpha), - gm=double(tmax),km=double(tmax), - gmic1=double(tmax),gmic2=double(tmax),kmic1=double(tmax),kmic2=double(tmax), - gmval=double(tmax),kmval=double(tmax), - PACKAGE="ads") - } - } - # formatting results - gm<-data.frame(obs=res$gm,theo=rep(0,tmax)) - km<-data.frame(obs=res$km,theo=rep(0,tmax)) - if(nsim>0) { - gm<-cbind(gm,sup=res$gmic1,inf=res$gmic2,pval=res$gmval/(nsim+1)) - km<-cbind(km,sup=res$kmic1,inf=res$kmic2,pval=res$kmval/(nsim+1)) - } - call<-match.call() - res<-list(call=call,r=r,gm=gm,km=km) - class(res)<-c("fads","kmfun") - return(res) -} - -ksfun<-function(p,upto,by,nsim=0,alpha=0.01) { -# checking for input parameters - #options( CBoundsCheck = TRUE ) - stopifnot(inherits(p,"spp")) - stopifnot(p$type=="multivariate") - stopifnot(is.numeric(upto)) - stopifnot(upto>=1) - stopifnot(is.numeric(by)) - stopifnot(by>0) - r<-seq(by,upto,by) - tmax<-length(r) - stopifnot(is.numeric(nsim)) - stopifnot(nsim>=0) - nsim<-testInteger(nsim) - stopifnot(is.numeric(alpha)) - stopifnot(alpha>=0) - if(nsim>0) testIC(nsim,alpha) - -###faire test sur les marks - - if("rectangle"%in%p$window$type) { - cas<-1 - xmin<-p$window$xmin - xmax<-p$window$xmax - ymin<-p$window$ymin - ymax<-p$window$ymax - stopifnot(upto<=(0.5*max((xmax-xmin),(ymax-ymin)))) - if ("complex"%in%p$window$type) { - cas<-3 - tri<-p$window$triangles - nbTri<-nrow(tri) - } - } - else if("circle"%in%p$window$type) { - cas<-2 - x0<-p$window$x0 - y0<-p$window$y0 - r0<-p$window$r0 - stopifnot(upto<=r0) - if ("complex"%in%p$window$type) { - cas<-4 - tri<-p$window$triangles - nbTri<-nrow(tri) - } - } - else - stop("invalid window type") - surface<-area.swin(p$window) - intensity<-p$n/surface - tabMarks<-levels(p$marks) - nbMarks<-nlevels(p$marks) -#nbMarks<-length(tabMarks) - marks<-as.numeric(p$marks) - freq<-as.vector(table(p$marks)) - D<-1-sum(freq*(freq-1))/(p$n*(p$n-1)) - -# computing Shimatani - if(cas==1) { #rectangle - if(nsim==0) { #without CI - res<-.C("shimatani_rect", - as.integer(p$n),as.double(p$x),as.double(p$y), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(tmax),as.double(by), - as.integer(nbMarks),as.integer(marks),as.double(surface),gg=double(tmax),kk=double(tmax),erreur=integer(tmax), - PACKAGE="ads") - } - else { #with CI - res<-.C("shimatani_rect_ic", - as.integer(p$n),as.double(p$x),as.double(p$y),as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(tmax),as.double(by), as.integer(nsim), as.double(alpha), - as.integer(nbMarks),as.integer(marks),as.double(surface),as.double(D), - gg=double(tmax),kk=double(tmax),gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), - gval=double(tmax),kval=double(tmax),erreur=integer(tmax), - PACKAGE="ads") - } - } - else if(cas==2) { #circle - if(nsim==0) { #without CI - res<-.C("shimatani_disq", - as.integer(p$n),as.double(p$x),as.double(p$y), - as.double(x0),as.double(y0),as.double(r0), - as.integer(tmax),as.double(by), - as.integer(nbMarks),as.integer(marks),as.double(surface),gg=double(tmax),kk=double(tmax),erreur=integer(tmax), - PACKAGE="ads") - } - else { #with CI - res<-.C("shimatani_disq_ic", - as.integer(p$n),as.double(p$x),as.double(p$y),as.double(x0),as.double(y0),as.double(r0), - as.integer(tmax),as.double(by), as.integer(nsim), as.double(alpha), - as.integer(nbMarks),as.integer(marks),as.double(surface),as.double(D), - gg=double(tmax),kk=double(tmax),gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), - gval=double(tmax),kval=double(tmax),erreur=integer(tmax), - PACKAGE="ads") - } - } - else if(cas==3) { #complex within rectangle - if(nsim==0) { #without CI - res<-.C("shimatani_tr_rect", - as.integer(p$n),as.double(p$x),as.double(p$y), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), - as.integer(nbMarks),as.integer(marks),as.double(surface),gg=double(tmax),kk=double(tmax),erreur=integer(tmax), - PACKAGE="ads") - } - else { #with CI - res<-.C("shimatani_tr_rect_ic", - as.integer(p$n),as.double(p$x),as.double(p$y),as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), as.integer(nsim), as.double(alpha), - as.integer(nbMarks),as.integer(marks),as.double(surface),as.double(D), - gg=double(tmax),kk=double(tmax),gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), - gval=double(tmax),kval=double(tmax),erreur=integer(tmax), - PACKAGE="ads") - } - } - else if(cas==4) { #complex within circle - if(nsim==0) { #without CI - res<-.C("shimatani_tr_disq", - as.integer(p$n),as.double(p$x),as.double(p$y), - as.double(x0),as.double(y0),as.double(r0), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), - as.integer(nbMarks),as.integer(marks),as.double(surface),gg=double(tmax),kk=double(tmax),erreur=integer(tmax), - PACKAGE="ads") - } - else { #with CI - res<-.C("shimatani_tr_disq_ic", - as.integer(p$n),as.double(p$x),as.double(p$y),as.double(x0),as.double(y0),as.double(r0), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), as.integer(nsim), as.double(alpha), - as.integer(nbMarks),as.integer(marks),as.double(surface),as.double(D), - gg=double(tmax),kk=double(tmax),gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), - gval=double(tmax),kval=double(tmax),erreur=integer(tmax), - PACKAGE="ads") - } - } - if(sum(res$erreur>0)){ - message("Error in ", appendLF=F) - print(match.call()) - message("No neigbors within distance intervals:") - print(paste(by*(res$erreur[res$erreur>0]-1),"-",by*res$erreur[res$erreur>0])) - message("Increase argument 'by'") - return(res=NULL) - } - gs<-data.frame(obs=res$gg/D,theo=rep(1,tmax)) - ks<-data.frame(obs=res$kk/D,theo=rep(1,tmax)) - if(nsim>0) { - gs<-cbind(gs,sup=res$gic1/D,inf=res$gic2/D,pval=res$gval/(nsim+1)) - ks<-cbind(ks,sup=res$kic1/D,inf=res$kic2/D,pval=res$kval/(nsim+1)) - } - call<-match.call() - res<-list(call=call,r=r,gs=gs,ks=ks) - class(res)<-c("fads","ksfun") - return(res) -} - - -################# -#V2 that calls K12fun -############## -krfun<-function(p,upto,by,nsim=0,dis=NULL,H0=c("rl","se"),alpha=0.01) { -# checking for input parameters - stopifnot(inherits(p,"spp")) - stopifnot(p$type=="multivariate") - stopifnot(is.numeric(upto)) - stopifnot(upto>=1) - stopifnot(is.numeric(by)) - stopifnot(by>0) - r<-seq(by,upto,by) - tmax<-length(r) - H0<-H0[1] - stopifnot(H0=="se" || H0=="rl") - ifelse(H0=="se",H0<-2,H0<-1) - if(is.null(dis)) { - stopifnot(H0==1) - dis<-as.dist(matrix(1,nlevels(p$marks),nlevels(p$marks))) - attr(dis,"Labels")<-levels(p$marks) - } - stopifnot(inherits(dis,"dist")) - stopifnot(attr(dis,"Diag")==FALSE) - stopifnot(attr(dis,"Upper")==FALSE) - stopifnot(suppressWarnings(is.euclid(dis))) -###revoir tests sur dis - if(length(levels(p$marks))!=length(labels(dis))) { - stopifnot(all(levels(p$marks)%in%labels(dis))) -#dis<-subsetdist(dis,which(labels(dis)%in%levels(p$marks))) - dis<-subsetdist(dis,levels(p$marks)) - warning("matrix 'dis' have been subsetted to match with levels(p$marks)") - } -#else if(any(labels(dis)!=levels(p$marks))) { -# attr(dis,"Labels")<-levels(p$marks) -# warning("levels(p$marks) have been assigned to attr(dis, ''Labels'')") -# } - stopifnot(is.numeric(nsim)) - stopifnot(nsim>=0) - nsim<-testInteger(nsim) - stopifnot(is.numeric(alpha)) - stopifnot(alpha>=0) - if(nsim>0) testIC(nsim,alpha) - -###faire test sur les marks - - if("rectangle"%in%p$window$type) { - cas<-1 - xmin<-p$window$xmin - xmax<-p$window$xmax - ymin<-p$window$ymin - ymax<-p$window$ymax - stopifnot(upto<=(0.5*max((xmax-xmin),(ymax-ymin)))) - if ("complex"%in%p$window$type) { - cas<-3 - tri<-p$window$triangles - nbTri<-nrow(tri) - } - } - else if("circle"%in%p$window$type) { - cas<-2 - x0<-p$window$x0 - y0<-p$window$y0 - r0<-p$window$r0 - stopifnot(upto<=r0) - if ("complex"%in%p$window$type) { - cas<-4 - tri<-p$window$triangles - nbTri<-nrow(tri) - } - } - else - stop("invalid window type") - surface<-area.swin(p$window) - intensity<-p$n/surface - nbMarks<-nlevels(p$marks) - marks<-as.numeric(p$marks) # => position du label dans levels(p$marks) - dis<-as.dist(sortmat(dis,levels(p$marks))) - HD<-suppressWarnings(divc(as.data.frame(unclass(table(p$marks))),sqrt(2*dis),scale=F)[1,1]) - HD<-HD*p$n/(p$n-1) - dis<-as.vector(dis) - -# computing Rao - if(cas==1) { #rectangle - if(nsim==0) { #without CI - res<-.C("rao_rect", - as.integer(p$n),as.double(p$x),as.double(p$y), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(tmax),as.double(by),as.integer(H0), - as.integer(nbMarks),as.integer(marks),as.double(dis),as.double(surface),as.double(HD), - gg=double(tmax),kk=double(tmax),gs=double(tmax),ks=double(tmax),erreur=integer(tmax), - PACKAGE="ads") - } - else { #with CI - res<-.C("rao_rect_ic", - as.integer(p$n),as.double(p$x),as.double(p$y),as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(tmax),as.double(by), as.integer(nsim),as.integer(H0),as.double(alpha), - as.integer(nbMarks),as.integer(marks),as.double(dis),as.double(surface),as.double(HD), - gg=double(tmax),kk=double(tmax),gs=double(tmax),ks=double(tmax),gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), - gval=double(tmax),kval=double(tmax),erreur=integer(tmax), - PACKAGE="ads") - } - } - else if(cas==2) { #circle - if(nsim==0) { #without CI - res<-.C("rao_disq", - as.integer(p$n),as.double(p$x),as.double(p$y), - as.double(x0),as.double(y0),as.double(r0), - as.integer(tmax),as.double(by),as.integer(H0), - as.integer(nbMarks),as.integer(marks),as.double(dis),as.double(surface),as.double(HD), - gg=double(tmax),kk=double(tmax),gs=double(tmax),ks=double(tmax),erreur=integer(tmax), - PACKAGE="ads") - } - else { #with CI - res<-.C("rao_disq_ic", - as.integer(p$n),as.double(p$x),as.double(p$y),as.double(x0),as.double(y0),as.double(r0), - as.integer(tmax),as.double(by), as.integer(nsim),as.integer(H0),as.double(alpha), - as.integer(nbMarks),as.integer(marks),as.double(dis),as.double(surface),as.double(HD), - gg=double(tmax),kk=double(tmax),gs=double(tmax),ks=double(tmax),gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), - gval=double(tmax),kval=double(tmax),erreur=integer(tmax), - PACKAGE="ads") - } - } - else if(cas==3) { #complex within rectangle - if(nsim==0) { #without CI - res<-.C("rao_tr_rect", - as.integer(p$n),as.double(p$x),as.double(p$y), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by),as.integer(H0), - as.integer(nbMarks),as.integer(marks),as.double(dis),as.double(surface),as.double(HD), - gg=double(tmax),kk=double(tmax),gs=double(tmax),ks=double(tmax),erreur=integer(tmax), - PACKAGE="ads") - } - else { #with CI - res<-.C("rao_tr_rect_ic", - as.integer(p$n),as.double(p$x),as.double(p$y), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by),as.integer(nsim), as.integer(H0),as.double(alpha), - as.integer(nbMarks),as.integer(marks),as.double(dis),as.double(surface),as.double(HD), - gg=double(tmax),kk=double(tmax),gs=double(tmax),ks=double(tmax),gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), - gval=double(tmax),kval=double(tmax),erreur=integer(tmax), - PACKAGE="ads") - } - } - else if(cas==4) { #complex within circle - if(nsim==0) { #without CI - res<-.C("rao_tr_disq", - as.integer(p$n),as.double(p$x),as.double(p$y), - as.double(x0),as.double(y0),as.double(r0), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by),as.integer(H0), - as.integer(nbMarks),as.integer(marks),as.double(dis),as.double(surface),as.double(HD), - gg=double(tmax),kk=double(tmax),gs=double(tmax),ks=double(tmax),erreur=integer(tmax), - PACKAGE="ads") - } - else { #with CI (not based on K12) - res<-.C("rao_tr_disq_ic", - as.integer(p$n),as.double(p$x),as.double(p$y), - as.double(x0),as.double(y0),as.double(r0), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by),as.integer(nsim), as.integer(H0),as.double(alpha), - as.integer(nbMarks),as.integer(marks),as.double(dis),as.double(surface),as.double(HD), - gg=double(tmax),kk=double(tmax),gs=double(tmax),ks=double(tmax),gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), - gval=double(tmax),kval=double(tmax),serreur=integer(tmax), - PACKAGE="ads") - } - } - if(sum(res$erreur>0)){ - message("Error in ", appendLF=F) - print(match.call()) - message("No neigbors within distance intervals:") - print(paste(by*(res$erreur[res$erreur>0]-1),"-",by*res$erreur[res$erreur>0])) - message("Increase argument 'by'") - return(res=NULL) - } - if(H0==1) { - theog<-rep(1,tmax) - theok<-rep(1,tmax) - } - if(H0==2) { - theog<-res$gs - theok<-res$ks - } - gr<-data.frame(obs=res$gg,theo=theog) - kr<-data.frame(obs=res$kk,theo=theok) - if(nsim>0) { - gr<-cbind(gr,sup=res$gic1,inf=res$gic2,pval=res$gval/(nsim+1)) - kr<-cbind(kr,sup=res$kic1,inf=res$kic2,pval=res$kval/(nsim+1)) - } - call<-match.call() - res<-list(call=call,r=r,gr=gr,kr=kr) - class(res)<-c("fads","krfun") - return(res) -} - -kdfun<-function(p,upto,by,dis,nsim=0,alpha=0.01) { -# checking for input parameters - stopifnot(inherits(p,"spp")) - stopifnot(p$type=="multivariate") - if(min(p$x)<0) - p$x<-p$x+abs(min(p$x)) - if(min(p$y)<0) - p$y<-p$y+abs(min(p$y)) - stopifnot(is.numeric(upto)) - stopifnot(upto>=1) - stopifnot(is.numeric(by)) - stopifnot(by>0) - r<-seq(by,upto,by) - tmax<-length(r) - stopifnot(inherits(dis,"dist")) - stopifnot(attr(dis,"Diag")==FALSE) - stopifnot(attr(dis,"Upper")==FALSE) - stopifnot(suppressWarnings(is.euclid(dis))) -###revoir tests sur dis - if(length(levels(p$marks))!=length(labels(dis))) { - stopifnot(all(levels(p$marks)%in%labels(dis))) -#dis<-subsetdist(dis,which(labels(dis)%in%levels(p$marks))) - dis<-subsetdist(dis,levels(p$marks)) - warning("matrix 'dis' have been subsetted to match with levels(p$marks)") - } -#else if(any(labels(dis)!=levels(p$marks))) { -# attr(dis,"Labels")<-levels(p$marks) -# warning("levels(p$marks) have been assigned to attr(dis, ''Labels'')") -# } - stopifnot(is.numeric(nsim)) - stopifnot(nsim>=0) - nsim<-testInteger(nsim) - stopifnot(is.numeric(alpha)) - stopifnot(alpha>=0) - if(nsim>0) testIC(nsim,alpha) - - surface<-area.swin(p$window) - intensity<-p$n/surface - nbMarks<-nlevels(p$marks) - marks<-as.numeric(p$marks) # => position du label dans levels(p$marks) - dis<-as.dist(sortmat(dis,levels(p$marks))) - HD<-suppressWarnings(divc(as.data.frame(unclass(table(p$marks))),sqrt(2*dis),scale=F)[1,1]) - HD<-HD*p$n/(p$n-1) - dis<-as.vector(dis) - -###faire test sur les marks - - if(nsim==0) { #without CI - res<-.C("shen", - as.integer(p$n),as.double(p$x),as.double(p$y), - as.integer(tmax),as.double(by), - as.integer(nbMarks),as.integer(marks),as.double(dis),as.double(surface),as.double(HD), - gd=double(tmax),kd=double(tmax),erreur=integer(tmax), - PACKAGE="ads") - } - else { #with CI - res<-.C("shen_ic", - as.integer(p$n),as.double(p$x),as.double(p$y), - as.integer(tmax),as.double(by), as.integer(nsim),as.double(alpha), - as.integer(nbMarks),as.integer(marks),as.double(dis),as.double(surface),as.double(HD), - gd=double(tmax),kd=double(tmax),gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), - gval=double(tmax),kval=double(tmax),erreur=integer(tmax), - PACKAGE="ads") - } - - if(sum(res$erreur>0)){ - message("Error in ", appendLF=F) - print(match.call()) - message("No neigbors within distance intervals:") - print(paste(by*(res$erreur[res$erreur>0]-1),"-",by*res$erreur[res$erreur>0])) - message("Increase argument 'by'") - return(res=NULL) - } - gd<-data.frame(obs=res$gd,theo=rep(1,tmax)) - kd<-data.frame(obs=res$kd,theo=rep(1,tmax)) - if(nsim>0) { - gd<-cbind(gd,sup=res$gic1,inf=res$gic2,pval=res$gval/(nsim+1)) - kd<-cbind(kd,sup=res$kic1,inf=res$kic2,pval=res$kval/(nsim+1)) - } - call<-match.call() - res<-list(call=call,r=r,gd=gd,kd=kd) - class(res)<-c("fads","kdfun") - return(res) -} +kfun<-function(p,upto,by,nsim=0,prec=0.01,alpha=0.01) { + # checking for input parameters + stopifnot(inherits(p,"spp")) + if(p$type!="univariate") + warning(paste(p$type,"point pattern has been considered to be univariate\n")) + stopifnot(is.numeric(upto)) + stopifnot(upto>=1) + stopifnot(is.numeric(by)) + stopifnot(by>0) + r<-seq(by,upto,by) + tmax<-length(r) + stopifnot(is.numeric(nsim)) + stopifnot(nsim>=0) + nsim<-testInteger(nsim) + stopifnot(is.numeric(prec)) + stopifnot(prec>=0) + stopifnot(is.numeric(alpha)) + stopifnot(alpha>=0) + if(nsim>0) testIC(nsim,alpha) + + if("rectangle"%in%p$window$type) { + cas<-1 + xmin<-p$window$xmin + xmax<-p$window$xmax + ymin<-p$window$ymin + ymax<-p$window$ymax + stopifnot(upto<=(0.5*max((xmax-xmin),(ymax-ymin)))) + if ("complex"%in%p$window$type) { + cas<-3 + tri<-p$window$triangles + nbTri<-nrow(tri) + } + } + else if("circle"%in%p$window$type) { + cas<-2 + x0<-p$window$x0 + y0<-p$window$y0 + r0<-p$window$r0 + stopifnot(upto<=r0) + if ("complex"%in%p$window$type) { + cas<-4 + tri<-p$window$triangles + nbTri<-nrow(tri) + } + } + else + stop("invalid window type") + intensity<-p$n/area.swin(p$window) + + if(cas==1) { #rectangle + if(nsim==0) { #without CI + res<-.C("ripley_rect", + as.integer(p$n),as.double(p$x),as.double(p$y), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(tmax),as.double(by), + g=double(tmax),k=double(tmax), + PACKAGE="ads") + } + else { #with CI + res<-.C("ripley_rect_ic", + as.integer(p$n),as.double(p$x),as.double(p$y), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax),as.double(intensity), + as.integer(tmax),as.double(by), + as.integer(nsim),as.double(prec),as.double(alpha), + g=double(tmax),k=double(tmax), + gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), + gval=double(tmax),kval=double(tmax),lval=double(tmax),nval=double(tmax), + PACKAGE="ads") + } + } + else if(cas==2) { #circle + if(nsim==0) { #without CI + res<-.C("ripley_disq", + as.integer(p$n),as.double(p$x),as.double(p$y), + as.double(x0),as.double(y0),as.double(r0), + as.integer(tmax),as.double(by), + g=double(tmax),k=double(tmax), + PACKAGE="ads") + } + else { #with CI + res<-.C("ripley_disq_ic", + as.integer(p$n),as.double(p$x),as.double(p$y), + as.double(x0),as.double(y0),as.double(r0),as.double(intensity), + as.integer(tmax),as.double(by), + as.integer(nsim),as.double(prec),as.double(alpha), + g=double(tmax),k=double(tmax), + gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), + gval=double(tmax),kval=double(tmax),lval=double(tmax),nval=double(tmax), + PACKAGE="ads") + } + } + else if(cas==3) { #complex within rectangle + if(nsim==0) { #without CI + res<-.C("ripley_tr_rect", + as.integer(p$n),as.double(p$x),as.double(p$y), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), + g=double(tmax),k=double(tmax), + PACKAGE="ads") + } + else { #with CI + res<-.C("ripley_tr_rect_ic", + as.integer(p$n),as.double(p$x),as.double(p$y), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax),as.double(intensity), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), + as.integer(nsim),as.double(prec),as.double(alpha), + g=double(tmax),k=double(tmax), + gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), + gval=double(tmax),kval=double(tmax),lval=double(tmax),nval=double(tmax), + PACKAGE="ads") + } + } + else if(cas==4) { #complex within circle + if(nsim==0) { #without CI + res<-.C("ripley_tr_disq", + as.integer(p$n),as.double(p$x),as.double(p$y), + as.double(x0),as.double(y0),as.double(r0), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), + g=double(tmax),k=double(tmax), + PACKAGE="ads") + } + else { #with CI + res<-.C("ripley_tr_disq_ic", + as.integer(p$n),as.double(p$x),as.double(p$y), + as.double(x0),as.double(y0),as.double(r0),as.double(intensity), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), + as.integer(nsim),as.double(prec),as.double(alpha), + g=double(tmax),k=double(tmax), + gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), + gval=double(tmax),kval=double(tmax),lval=double(tmax),nval=double(tmax), + PACKAGE="ads") + } + } + # formatting results + ds<-c(pi,diff(pi*r^2)) + g<-data.frame(obs=res$g/(intensity*ds),theo=rep(1,tmax)) + n<-data.frame(obs=res$k/(pi*r^2),theo=rep(intensity,tmax)) + k<-data.frame(obs=res$k/intensity,theo=pi*r^2) + l<-data.frame(obs=sqrt(res$k/(intensity*pi))-r,theo=rep(0,tmax)) + if(nsim>0) { + g<-cbind(g,sup=res$gic1/(intensity*ds),inf=res$gic2/(intensity*ds),pval=res$gval/(nsim+1)) + n<-cbind(n,sup=res$kic1/(pi*r^2),inf=res$kic2/(pi*r^2),pval=res$nval/(nsim+1)) + k<-cbind(k,sup=res$kic1/intensity,inf=res$kic2/intensity,pval=res$kval/(nsim+1)) + l<-cbind(l,sup=sqrt(res$kic1/(intensity*pi))-r,inf=sqrt(res$kic2/(intensity*pi))-r,pval=res$lval/(nsim+1)) + } + call<-match.call() + res<-list(call=call,r=r,g=g,n=n,k=k,l=l) + class(res)<-c("fads","kfun") + return(res) +} + +k12fun<-function(p,upto,by,nsim=0,H0=c("pitor","pimim","rl"),prec=0.01,nsimax=3000,conv=50,rep=10,alpha=0.01,marks) { + # checking for input parameters + options( CBoundsCheck = TRUE ) + # regle les problemes pour 32-bit + stopifnot(inherits(p,"spp")) + stopifnot(p$type=="multivariate") + stopifnot(is.numeric(upto)) + stopifnot(upto>=1) + stopifnot(is.numeric(by)) + stopifnot(by>0) + r<-seq(by,upto,by) + tmax<-length(r) + stopifnot(is.numeric(nsim)) + stopifnot(nsim>=0) + nsim<-testInteger(nsim) + H0<-H0[1] + stopifnot(H0=="pitor" || H0=="pimim" || H0=="rl") + if(H0=="rl") H0<-1 + else if(H0=="pitor") H0<-2 + else H0<-3 + stopifnot(is.numeric(prec)) + stopifnot(prec>=0) + stopifnot(is.numeric(alpha)) + stopifnot(alpha>=0) + if(nsim>0) testIC(nsim,alpha) + if(missing(marks)) + marks<-c(1,2) + stopifnot(length(marks)==2) + stopifnot(marks[1]!=marks[2]) + mark1<-marks[1] + mark2<-marks[2] + if(is.numeric(mark1)) + mark1<-levels(p$marks)[testInteger(mark1)] + else if(!mark1%in%p$marks) stop(paste("mark \'",mark1,"\' doesn\'t exist",sep="")) + if(is.numeric(mark2)) + mark2<-levels(p$marks)[testInteger(mark2)] + else if(!mark2%in%p$marks) stop(paste("mark \'",mark2,"\' doesn\'t exist",sep="")) + + if("rectangle"%in%p$window$type) { + cas<-1 + xmin<-p$window$xmin + xmax<-p$window$xmax + ymin<-p$window$ymin + ymax<-p$window$ymax + stopifnot(upto<=(0.5*max((xmax-xmin),(ymax-ymin)))) + if ("complex"%in%p$window$type) { + cas<-3 + tri<-p$window$triangles + nbTri<-nrow(tri) + } + } + else if("circle"%in%p$window$type) { + cas<-2 + x0<-p$window$x0 + y0<-p$window$y0 + r0<-p$window$r0 + stopifnot(upto<=r0) + if ("complex"%in%p$window$type) { + cas<-4 + tri<-p$window$triangles + nbTri<-nrow(tri) + } + } + else + stop("invalid window type") + surface<-area.swin(p$window) + x1<-p$x[p$marks==mark1] + y1<-p$y[p$marks==mark1] + x2<-p$x[p$marks==mark2] + y2<-p$y[p$marks==mark2] + nbPts1<-length(x1) + nbPts2<-length(x2) + intensity2<-nbPts2/surface +# intensity<-(nbPts1+nbPts2)/surface + + # computing intertype functions + if(cas==1) { #rectangle + if(nsim==0) { #without CI + res<-.C("intertype_rect", + as.integer(nbPts1),as.double(x1),as.double(y1), + as.integer(nbPts2),as.double(x2),as.double(y2), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(tmax),as.double(by), + g=double(tmax),k=double(tmax), + PACKAGE="ads") + } + else { #with CI + res<-.C("intertype_rect_ic", + as.integer(nbPts1),as.double(x1),as.double(y1), + as.integer(nbPts2),as.double(x2),as.double(y2), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax),as.double(surface), + as.integer(tmax),as.double(by), + as.integer(nsim),as.integer(H0),as.double(prec),as.integer(nsimax),as.integer(conv),as.integer(rep),as.double(alpha), + g=double(tmax),k=double(tmax), + gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), + gval=double(tmax),kval=double(tmax),lval=double(tmax),nval=double(tmax), + PACKAGE="ads") + } + } + else if(cas==2) { #circle + if(nsim==0) { #without CI + res<-.C("intertype_disq", + as.integer(nbPts1),as.double(x1),as.double(y1), + as.integer(nbPts2),as.double(x2),as.double(y2), + as.double(x0),as.double(y0),as.double(r0), + as.integer(tmax),as.double(by), + g=double(tmax),k=double(tmax), + PACKAGE="ads") + } + else { #with CI + res<-.C("intertype_disq_ic", + as.integer(nbPts1),as.double(x1),as.double(y1), + as.integer(nbPts2),as.double(x2),as.double(y2), + as.double(x0),as.double(y0),as.double(r0),as.double(surface), + as.integer(tmax),as.double(by), + as.integer(nsim),as.integer(H0),as.double(prec),as.integer(nsimax),as.integer(conv),as.integer(rep),as.double(alpha), + g=double(tmax),k=double(tmax), + gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), + gval=double(tmax),kval=double(tmax),lval=double(tmax),nval=double(tmax), + PACKAGE="ads") + } + } + else if(cas==3) { #complex within rectangle + if(nsim==0) { #without CI + res<-.C("intertype_tr_rect", + as.integer(nbPts1),as.double(x1),as.double(y1), + as.integer(nbPts2),as.double(x2),as.double(y2), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), + g=double(tmax),k=double(tmax), + PACKAGE="ads") + } + else { #with CI + res<-.C("intertype_tr_rect_ic", + as.integer(nbPts1),as.double(x1),as.double(y1), + as.integer(nbPts2),as.double(x2),as.double(y2), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax),as.double(surface), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), + as.integer(nsim),as.integer(H0),as.double(prec),as.integer(nsimax),as.integer(conv),as.integer(rep),as.double(alpha), + g=double(tmax),k=double(tmax), + gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), + gval=double(tmax),kval=double(tmax),lval=double(tmax),nval=double(tmax), + PACKAGE="ads") + } + } + else if(cas==4) { #complex within circle + if(nsim==0) { #without CI + res<-.C("intertype_tr_disq", + as.integer(nbPts1),as.double(x1),as.double(y1), + as.integer(nbPts2),as.double(x2),as.double(y2), + as.double(x0),as.double(y0),as.double(r0), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), + g=double(tmax),k=double(tmax), + PACKAGE="ads") + } + else { #with CI + res<-.C("intertype_tr_disq_ic", + as.integer(nbPts1),as.double(x1),as.double(y1), + as.integer(nbPts2),as.double(x2),as.double(y2), + as.double(x0),as.double(y0),as.double(r0),as.double(surface), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), + as.integer(nsim),as.integer(H0),as.double(prec),as.integer(nsimax),as.integer(conv),as.integer(rep),as.double(alpha), + g=double(tmax),k=double(tmax), + gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), + gval=double(tmax),kval=double(tmax),lval=double(tmax),nval=double(tmax), + PACKAGE="ads") + } + } + # formatting results + ds<-c(pi,diff(pi*r^2)) + g<-res$g/(intensity2*ds) + n<-res$k/(pi*r^2) + k<-res$k/intensity2 + l<-sqrt(res$k/(intensity2*pi))-r + if(H0==1) { + rip<-kfun(spp(c(x1,x2),c(y1,y2),p$window),upto,by) + theo<-list(g=rip$g$obs,n=intensity2*rip$k$obs/(pi*r^2),k=rip$k$obs,l=rip$l$obs) + } + else if (H0==2||H0==3) + theo<-list(g=rep(1,tmax),n=rep(intensity2,tmax),k=pi*r^2,l=rep(0,tmax)) + g<-data.frame(obs=g,theo=theo$g) + n<-data.frame(obs=n,theo=theo$n) + k<-data.frame(obs=k,theo=theo$k) + l<-data.frame(obs=l,theo=theo$l) + if(nsim>0) { + g<-cbind(g,sup=res$gic1/(intensity2*ds),inf=res$gic2/(intensity2*ds),pval=res$gval/(nsim+1)) + n<-cbind(n,sup=res$kic1/(pi*r^2),inf=res$kic2/(pi*r^2),pval=res$nval/(nsim+1)) + k<-cbind(k,sup=res$kic1/intensity2,inf=res$kic2/intensity2,pval=res$kval/(nsim+1)) + l<-cbind(l,sup=sqrt(res$kic1/(intensity2*pi))-r,inf=sqrt(res$kic2/(intensity2*pi))-r,pval=res$lval/(nsim+1)) + } + call<-match.call() + res<-list(call=call,r=r,g12=g,n12=n,k12=k,l12=l,marks=c(mark1,mark2)) + class(res)<-c("fads","k12fun") + return(res) +} + +kijfun<-kpqfun<-function(p,upto,by) { +# checking for input parameters + stopifnot(inherits(p,"spp")) + stopifnot(p$type=="multivariate") + stopifnot(is.numeric(upto)) + stopifnot(upto>=1) + stopifnot(is.numeric(by)) + stopifnot(by>0) + r<-seq(by,upto,by) + tmax<-length(r) + + if("rectangle"%in%p$window$type) { + cas<-1 + xmin<-p$window$xmin + xmax<-p$window$xmax + ymin<-p$window$ymin + ymax<-p$window$ymax + stopifnot(upto<=(0.5*max((xmax-xmin),(ymax-ymin)))) + if ("complex"%in%p$window$type) { + cas<-3 + tri<-p$window$triangles + nbTri<-nrow(tri) + } + } + else if("circle"%in%p$window$type) { + cas<-2 + x0<-p$window$x0 + y0<-p$window$y0 + r0<-p$window$r0 + stopifnot(upto<=r0) + if ("complex"%in%p$window$type) { + cas<-4 + tri<-p$window$triangles + nbTri<-nrow(tri) + } + } + else + stop("invalid window type") + surface<-area.swin(p$window) + + tabMarks<-levels(p$marks) + nbMarks<-length(tabMarks) + mpt_nb<-summary(p$marks) +# computing RipleyClass + gij<-double(tmax*nbMarks^2) + kij<-double(tmax*nbMarks^2) + lij<-double(tmax*nbMarks^2) + nij<-double(tmax*nbMarks^2) + for(i in 1:nbMarks) { + x1<-p$x[p$marks==tabMarks[i]] + y1<-p$y[p$marks==tabMarks[i]] + if(cas==1) { #rectangle + res<-.C("ripley_rect", + as.integer(mpt_nb[i]),as.double(x1),as.double(y1), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(tmax),as.double(by), + g=double(tmax),k=double(tmax), + PACKAGE="ads") + } + else if(cas==2) { #circle + res<-.C("ripley_disq", + as.integer(mpt_nb[i]),as.double(x1),as.double(y1), + as.double(x0),as.double(y0),as.double(r0), + as.integer(tmax),as.double(by), + g=double(tmax),k=double(tmax), + PACKAGE="ads") + } + else if(cas==3) { #complex within rectangle + res<-.C("ripley_tr_rect", + as.integer(mpt_nb[i]),as.double(x1),as.double(y1), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), + g=double(tmax),k=double(tmax), + PACKAGE="ads") + } + else if(cas==4) { #complex within circle + res<-.C("ripley_tr_disq", + as.integer(mpt_nb[i]),as.double(x1),as.double(y1), + as.double(x0),as.double(y0),as.double(r0), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), + g=double(tmax),k=double(tmax), + PACKAGE="ads") + } + intensity1<-mpt_nb[i]/surface + matcol<-(i-1)*nbMarks+i-1 + j<-(matcol*tmax+1):(matcol*tmax+tmax) + ds<-c(pi,diff(pi*r^2)) + gij[j]<-res$g/(intensity1*ds) + nij[j]<-res$k/(pi*r^2) + kij[j]<-res$k/intensity1 + lij[j]<-sqrt(res$k/(intensity1*pi))-r + if(i<nbMarks) { + for(j in (i+1):nbMarks) { + x2<-p$x[p$marks==tabMarks[j]] + y2<-p$y[p$marks==tabMarks[j]] + if(cas==1) { #rectangle + res<-.C("intertype_rect", + as.integer(mpt_nb[i]),as.double(x1),as.double(y1), + as.integer(mpt_nb[j]),as.double(x2),as.double(y2), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(tmax),as.double(by), + g=double(tmax),k=double(tmax), + PACKAGE="ads") + } + else if(cas==2) { #circle + res<-.C("intertype_disq", + as.integer(mpt_nb[i]),as.double(x1),as.double(y1), + as.integer(mpt_nb[j]),as.double(x2),as.double(y2), + as.double(x0),as.double(y0),as.double(r0), + as.integer(tmax),as.double(by), + g=double(tmax),k=double(tmax), + PACKAGE="ads") + } + else if(cas==3) { #complex within rectangle + res<-.C("intertype_tr_rect", + as.integer(mpt_nb[i]),as.double(x1),as.double(y1), + as.integer(mpt_nb[j]),as.double(x2),as.double(y2), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), + g=double(tmax),k=double(tmax), + PACKAGE="ads") + } + else if(cas==4) { #complex within circle + res<-.C("intertype_tr_disq", + as.integer(mpt_nb[i]),as.double(x1),as.double(y1), + as.integer(mpt_nb[j]),as.double(x2),as.double(y2), + as.double(x0),as.double(y0),as.double(r0), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), + g=double(tmax),k=double(tmax), + PACKAGE="ads") + } + intensity2<-mpt_nb[j]/surface + matcol<-(i-1)*nbMarks+j-1 + k<-(matcol*tmax+1):(matcol*tmax+tmax) + gij[k]<-res$g/(intensity2*ds) + nij[k]<-res$k/(pi*r^2) + kij[k]<-res$k/intensity2 + lij[k]<-sqrt(res$k/(intensity2*pi))-r + if(cas==1) { #rectangle + res<-.C("intertype_rect", + as.integer(mpt_nb[j]),as.double(x2),as.double(y2), + as.integer(mpt_nb[i]),as.double(x1),as.double(y1), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(tmax),as.double(by), + g=double(tmax),k=double(tmax), + PACKAGE="ads") + } + else if(cas==2) { #circle + res<-.C("intertype_disq", + as.integer(mpt_nb[j]),as.double(x2),as.double(y2), + as.integer(mpt_nb[i]),as.double(x1),as.double(y1), + as.double(x0),as.double(y0),as.double(r0), + as.integer(tmax),as.double(by), + g=double(tmax),k=double(tmax), + PACKAGE="ads") + } + else if(cas==3) { #complex within rectangle + res<-.C("intertype_tr_rect", + as.integer(mpt_nb[j]),as.double(x2),as.double(y2), + as.integer(mpt_nb[i]),as.double(x1),as.double(y1), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), + g=double(tmax),k=double(tmax), + PACKAGE="ads") + } + else if(cas==4) { #complex within circle + res<-.C("intertype_tr_disq", + as.integer(mpt_nb[j]),as.double(x2),as.double(y2), + as.integer(mpt_nb[i]),as.double(x1),as.double(y1), + as.double(x0),as.double(y0),as.double(r0), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), + g=double(tmax),k=double(tmax), + PACKAGE="ads") + } + matcol<-(j-1)*nbMarks+i-1 + k<-(matcol*tmax+1):(matcol*tmax+tmax) + gij[k]<-res$g/(intensity1*ds) + nij[k]<-res$k/(pi*r^2) + kij[k]<-res$k/intensity1 + lij[k]<-sqrt(res$k/(intensity1*pi))-r + } + } + } + labij<-paste(rep(tabMarks,each=nbMarks),rep(tabMarks,nbMarks),sep="-") + gij<-matrix(gij,nrow=tmax,ncol=nbMarks^2) + kij<-matrix(kij,nrow=tmax,ncol=nbMarks^2) + nij<-matrix(nij,nrow=tmax,ncol=nbMarks^2) + lij<-matrix(lij,nrow=tmax,ncol=nbMarks^2) + call<-match.call() + res<-list(call=call,r=r,labpq=labij,gij=gij,kpq=kij,lpq=lij,npq=nij,intensity=summary(p)$intensity) + class(res)<-c("fads","kpqfun") + return(res) +} + +ki.fun<-kp.fun<-function(p,upto,by) { + # checking for input parameters + stopifnot(inherits(p,"spp")) + stopifnot(p$type=="multivariate") + stopifnot(is.numeric(upto)) + stopifnot(upto>=1) + stopifnot(is.numeric(by)) + stopifnot(by>0) + r<-seq(by,upto,by) + tmax<-length(r) + + if("rectangle"%in%p$window$type) { + cas<-1 + xmin<-p$window$xmin + xmax<-p$window$xmax + ymin<-p$window$ymin + ymax<-p$window$ymax + stopifnot(upto<=(0.5*max((xmax-xmin),(ymax-ymin)))) + if ("complex"%in%p$window$type) { + cas<-3 + tri<-p$window$triangles + nbTri<-nrow(tri) + } + } + else if("circle"%in%p$window$type) { + cas<-2 + x0<-p$window$x0 + y0<-p$window$y0 + r0<-p$window$r0 + stopifnot(upto<=r0) + if ("complex"%in%p$window$type) { + cas<-4 + tri<-p$window$triangles + nbTri<-nrow(tri) + } + } + else + stop("invalid window type") + surface<-area.swin(p$window) + + tabMarks<-levels(p$marks) + nbMarks<-length(tabMarks) + mpt_nb<-summary(p$marks) + #computing RipleyAll + gis<-double(tmax*nbMarks) + kis<-double(tmax*nbMarks) + lis<-double(tmax*nbMarks) + nis<-double(tmax*nbMarks) + for(i in 1:nbMarks) { + x1<-p$x[p$marks==tabMarks[i]] + y1<-p$y[p$marks==tabMarks[i]] + x2<-p$x[p$marks!=tabMarks[i]] + y2<-p$y[p$marks!=tabMarks[i]] + nbPts1<-mpt_nb[i] + nbPts2<-sum(mpt_nb)-nbPts1 + if(cas==1) { #rectangle + res<-.C("intertype_rect", + as.integer(nbPts1),as.double(x1),as.double(y1), + as.integer(nbPts2),as.double(x2),as.double(y2), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(tmax),as.double(by), + g=double(tmax),k=double(tmax), + PACKAGE="ads") + } + else if(cas==2) { #circle + res<-.C("intertype_disq", + as.integer(nbPts1),as.double(x1),as.double(y1), + as.integer(nbPts2),as.double(x2),as.double(y2), + as.double(x0),as.double(y0),as.double(r0), + as.integer(tmax),as.double(by), + g=double(tmax),k=double(tmax), + PACKAGE="ads") + } + else if(cas==3) { #complex within rectangle + res<-.C("intertype_tr_rect", + as.integer(nbPts1),as.double(x1),as.double(y1), + as.integer(nbPts2),as.double(x2),as.double(y2), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), + g=double(tmax),k=double(tmax), + PACKAGE="ads") + } + else if(cas==4) { #complex within circle + res<-.C("intertype_tr_disq", + as.integer(nbPts1),as.double(x1),as.double(y1), + as.integer(nbPts2),as.double(x2),as.double(y2), + as.double(x0),as.double(y0),as.double(r0), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), + g=double(tmax),k=double(tmax), + PACKAGE="ads") + } + intensity2<-nbPts2/surface + j<-((i-1)*tmax+1):((i-1)*tmax+tmax) + ds<-c(pi,diff(pi*r^2)) + gis[j]<-res$g/(intensity2*ds) + nis[j]<-res$k/(pi*r^2) + kis[j]<-res$k/intensity2 + lis[j]<-sqrt(res$k/(intensity2*pi))-r + } + # formatting results + labi<-tabMarks + gis<-matrix(gis,nrow=tmax,ncol=nbMarks) + kis<-matrix(kis,nrow=tmax,ncol=nbMarks) + nis<-matrix(nis,nrow=tmax,ncol=nbMarks) + lis<-matrix(lis,nrow=tmax,ncol=nbMarks) + call<-match.call() + res<-list(call=call,r=r,labp=labi,gp.=gis,kp.=kis,lp.=lis,np.=nis,intensity=summary(p)$intensity) + class(res)<-c("fads","kp.fun") + return(res) +} + +kmfun<-function(p,upto,by,nsim=0,alpha=0.01) { + # checking for input parameters + stopifnot(inherits(p,"spp")) + stopifnot(p$type=="marked") + stopifnot(is.numeric(upto)) + stopifnot(upto>=1) + stopifnot(is.numeric(by)) + stopifnot(by>0) + r<-seq(by,upto,by) + tmax<-length(r) + stopifnot(is.numeric(nsim)) + stopifnot(nsim>=0) + nsim<-testInteger(nsim) + #cmoy<-mean(p$marks) + cvar<-var(p$marks) + stopifnot(is.numeric(alpha)) + stopifnot(alpha>=0) + if(nsim>0) testIC(nsim,alpha) + + if("rectangle"%in%p$window$type) { + cas<-1 + xmin<-p$window$xmin + xmax<-p$window$xmax + ymin<-p$window$ymin + ymax<-p$window$ymax + stopifnot(upto<=(0.5*max((xmax-xmin),(ymax-ymin)))) + if ("complex"%in%p$window$type) { + cas<-3 + tri<-p$window$triangles + nbTri<-nrow(tri) + } + } + else if("circle"%in%p$window$type) { + cas<-2 + x0<-p$window$x0 + y0<-p$window$y0 + r0<-p$window$r0 + stopifnot(upto<=r0) + if ("complex"%in%p$window$type) { + cas<-4 + tri<-p$window$triangles + nbTri<-nrow(tri) + } + } + else + stop("invalid window type") + intensity<-p$n/area.swin(p$window) + + if(cas==1) { #rectangle + if(nsim==0) { #without CI + res<-.C("corr_rect", + as.integer(p$n),as.double(p$x),as.double(p$y),as.double(p$marks), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(tmax),as.double(by), + gm=double(tmax),km=double(tmax), + PACKAGE="ads") + } + else { #with CI + res<-.C("corr_rect_ic", + as.integer(p$n),as.double(p$x),as.double(p$y),as.double(p$marks), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(tmax),as.double(by), + as.integer(nsim),as.double(alpha), + gm=double(tmax),km=double(tmax), + gmic1=double(tmax),gmic2=double(tmax),kmic1=double(tmax),kmic2=double(tmax), + gmval=double(tmax),kmval=double(tmax), + PACKAGE="ads") + } + } + else if(cas==2) { #circle + if(nsim==0) { #without CI + res<-.C("corr_disq", + as.integer(p$n),as.double(p$x),as.double(p$y),as.double(p$marks), + as.double(x0),as.double(y0),as.double(r0), + as.integer(tmax),as.double(by), + gm=double(tmax),km=double(tmax), + PACKAGE="ads") + } + else { #with CI + res<-.C("corr_disq_ic", + as.integer(p$n),as.double(p$x),as.double(p$y),as.double(p$marks), + as.double(x0),as.double(y0),as.double(r0), + as.integer(tmax),as.double(by), + as.integer(nsim),as.double(alpha), + gm=double(tmax),km=double(tmax), + gmic1=double(tmax),gmic2=double(tmax),kmic1=double(tmax),kmic2=double(tmax), + gmval=double(tmax),kmval=double(tmax), + PACKAGE="ads") + } + } + else if(cas==3) { #complex within rectangle + if(nsim==0) { #without CI + res<-.C("corr_tr_rect", + as.integer(p$n),as.double(p$x),as.double(p$y),as.double(p$marks), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), + gm=double(tmax),km=double(tmax), + PACKAGE="ads") + } + else { #with CI + res<-.C("corr_tr_rect_ic", + as.integer(p$n),as.double(p$x),as.double(p$y),as.double(p$marks), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), + as.integer(nsim),as.double(alpha), + gm=double(tmax),km=double(tmax), + gmic1=double(tmax),gmic2=double(tmax),kmic1=double(tmax),kmic2=double(tmax), + gmval=double(tmax),kmval=double(tmax), + PACKAGE="ads") + } + } + else if(cas==4) { #complex within circle + if(nsim==0) { #without CI + res<-.C("corr_tr_disq", + as.integer(p$n),as.double(p$x),as.double(p$y),as.double(p$marks), + as.double(x0),as.double(y0),as.double(r0), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), + gm=double(tmax),km=double(tmax), + PACKAGE="ads") + } + else { #with CI + res<-.C("corr_tr_disq_ic", + as.integer(p$n),as.double(p$x),as.double(p$y),as.double(p$marks), + as.double(x0),as.double(y0),as.double(r0), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), + as.integer(nsim),as.double(alpha), + gm=double(tmax),km=double(tmax), + gmic1=double(tmax),gmic2=double(tmax),kmic1=double(tmax),kmic2=double(tmax), + gmval=double(tmax),kmval=double(tmax), + PACKAGE="ads") + } + } + # formatting results + gm<-data.frame(obs=res$gm,theo=rep(0,tmax)) + km<-data.frame(obs=res$km,theo=rep(0,tmax)) + if(nsim>0) { + gm<-cbind(gm,sup=res$gmic1,inf=res$gmic2,pval=res$gmval/(nsim+1)) + km<-cbind(km,sup=res$kmic1,inf=res$kmic2,pval=res$kmval/(nsim+1)) + } + call<-match.call() + res<-list(call=call,r=r,gm=gm,km=km) + class(res)<-c("fads","kmfun") + return(res) +} + +ksfun<-function(p,upto,by,nsim=0,alpha=0.01) { +# checking for input parameters + #options( CBoundsCheck = TRUE ) + stopifnot(inherits(p,"spp")) + stopifnot(p$type=="multivariate") + stopifnot(is.numeric(upto)) + stopifnot(upto>=1) + stopifnot(is.numeric(by)) + stopifnot(by>0) + r<-seq(by,upto,by) + tmax<-length(r) + stopifnot(is.numeric(nsim)) + stopifnot(nsim>=0) + nsim<-testInteger(nsim) + stopifnot(is.numeric(alpha)) + stopifnot(alpha>=0) + if(nsim>0) testIC(nsim,alpha) + +###faire test sur les marks + + if("rectangle"%in%p$window$type) { + cas<-1 + xmin<-p$window$xmin + xmax<-p$window$xmax + ymin<-p$window$ymin + ymax<-p$window$ymax + stopifnot(upto<=(0.5*max((xmax-xmin),(ymax-ymin)))) + if ("complex"%in%p$window$type) { + cas<-3 + tri<-p$window$triangles + nbTri<-nrow(tri) + } + } + else if("circle"%in%p$window$type) { + cas<-2 + x0<-p$window$x0 + y0<-p$window$y0 + r0<-p$window$r0 + stopifnot(upto<=r0) + if ("complex"%in%p$window$type) { + cas<-4 + tri<-p$window$triangles + nbTri<-nrow(tri) + } + } + else + stop("invalid window type") + surface<-area.swin(p$window) + intensity<-p$n/surface + tabMarks<-levels(p$marks) + nbMarks<-nlevels(p$marks) +#nbMarks<-length(tabMarks) + marks<-as.numeric(p$marks) + freq<-as.vector(table(p$marks)) + D<-1-sum(freq*(freq-1))/(p$n*(p$n-1)) + +# computing Shimatani + if(cas==1) { #rectangle + if(nsim==0) { #without CI + res<-.C("shimatani_rect", + as.integer(p$n),as.double(p$x),as.double(p$y), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(tmax),as.double(by), + as.integer(nbMarks),as.integer(marks),as.double(surface),gg=double(tmax),kk=double(tmax),erreur=integer(tmax), + PACKAGE="ads") + } + else { #with CI + res<-.C("shimatani_rect_ic", + as.integer(p$n),as.double(p$x),as.double(p$y),as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(tmax),as.double(by), as.integer(nsim), as.double(alpha), + as.integer(nbMarks),as.integer(marks),as.double(surface),as.double(D), + gg=double(tmax),kk=double(tmax),gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), + gval=double(tmax),kval=double(tmax),erreur=integer(tmax), + PACKAGE="ads") + } + } + else if(cas==2) { #circle + if(nsim==0) { #without CI + res<-.C("shimatani_disq", + as.integer(p$n),as.double(p$x),as.double(p$y), + as.double(x0),as.double(y0),as.double(r0), + as.integer(tmax),as.double(by), + as.integer(nbMarks),as.integer(marks),as.double(surface),gg=double(tmax),kk=double(tmax),erreur=integer(tmax), + PACKAGE="ads") + } + else { #with CI + res<-.C("shimatani_disq_ic", + as.integer(p$n),as.double(p$x),as.double(p$y),as.double(x0),as.double(y0),as.double(r0), + as.integer(tmax),as.double(by), as.integer(nsim), as.double(alpha), + as.integer(nbMarks),as.integer(marks),as.double(surface),as.double(D), + gg=double(tmax),kk=double(tmax),gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), + gval=double(tmax),kval=double(tmax),erreur=integer(tmax), + PACKAGE="ads") + } + } + else if(cas==3) { #complex within rectangle + if(nsim==0) { #without CI + res<-.C("shimatani_tr_rect", + as.integer(p$n),as.double(p$x),as.double(p$y), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), + as.integer(nbMarks),as.integer(marks),as.double(surface),gg=double(tmax),kk=double(tmax),erreur=integer(tmax), + PACKAGE="ads") + } + else { #with CI + res<-.C("shimatani_tr_rect_ic", + as.integer(p$n),as.double(p$x),as.double(p$y),as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), as.integer(nsim), as.double(alpha), + as.integer(nbMarks),as.integer(marks),as.double(surface),as.double(D), + gg=double(tmax),kk=double(tmax),gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), + gval=double(tmax),kval=double(tmax),erreur=integer(tmax), + PACKAGE="ads") + } + } + else if(cas==4) { #complex within circle + if(nsim==0) { #without CI + res<-.C("shimatani_tr_disq", + as.integer(p$n),as.double(p$x),as.double(p$y), + as.double(x0),as.double(y0),as.double(r0), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), + as.integer(nbMarks),as.integer(marks),as.double(surface),gg=double(tmax),kk=double(tmax),erreur=integer(tmax), + PACKAGE="ads") + } + else { #with CI + res<-.C("shimatani_tr_disq_ic", + as.integer(p$n),as.double(p$x),as.double(p$y),as.double(x0),as.double(y0),as.double(r0), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), as.integer(nsim), as.double(alpha), + as.integer(nbMarks),as.integer(marks),as.double(surface),as.double(D), + gg=double(tmax),kk=double(tmax),gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), + gval=double(tmax),kval=double(tmax),erreur=integer(tmax), + PACKAGE="ads") + } + } + if(sum(res$erreur>0)){ + message("Error in ", appendLF=F) + print(match.call()) + message("No neigbors within distance intervals:") + print(paste(by*(res$erreur[res$erreur>0]-1),"-",by*res$erreur[res$erreur>0])) + message("Increase argument 'by'") + return(res=NULL) + } + gs<-data.frame(obs=res$gg/D,theo=rep(1,tmax)) + ks<-data.frame(obs=res$kk/D,theo=rep(1,tmax)) + if(nsim>0) { + gs<-cbind(gs,sup=res$gic1/D,inf=res$gic2/D,pval=res$gval/(nsim+1)) + ks<-cbind(ks,sup=res$kic1/D,inf=res$kic2/D,pval=res$kval/(nsim+1)) + } + call<-match.call() + res<-list(call=call,r=r,gs=gs,ks=ks) + class(res)<-c("fads","ksfun") + return(res) +} + + +################# +#V2 that calls K12fun +############## +krfun<-function(p,upto,by,nsim=0,dis=NULL,H0=c("rl","se"),alpha=0.01) { +# checking for input parameters + stopifnot(inherits(p,"spp")) + stopifnot(p$type=="multivariate") + stopifnot(is.numeric(upto)) + stopifnot(upto>=1) + stopifnot(is.numeric(by)) + stopifnot(by>0) + r<-seq(by,upto,by) + tmax<-length(r) + H0<-H0[1] + stopifnot(H0=="se" || H0=="rl") + ifelse(H0=="se",H0<-2,H0<-1) + if(is.null(dis)) { + stopifnot(H0==1) + dis<-as.dist(matrix(1,nlevels(p$marks),nlevels(p$marks))) + attr(dis,"Labels")<-levels(p$marks) + } + stopifnot(inherits(dis,"dist")) + stopifnot(attr(dis,"Diag")==FALSE) + stopifnot(attr(dis,"Upper")==FALSE) + stopifnot(suppressWarnings(is.euclid(dis))) +###revoir tests sur dis + if(length(levels(p$marks))!=length(labels(dis))) { + stopifnot(all(levels(p$marks)%in%labels(dis))) +#dis<-subsetdist(dis,which(labels(dis)%in%levels(p$marks))) + dis<-subsetdist(dis,levels(p$marks)) + warning("matrix 'dis' have been subsetted to match with levels(p$marks)") + } +#else if(any(labels(dis)!=levels(p$marks))) { +# attr(dis,"Labels")<-levels(p$marks) +# warning("levels(p$marks) have been assigned to attr(dis, ''Labels'')") +# } + stopifnot(is.numeric(nsim)) + stopifnot(nsim>=0) + nsim<-testInteger(nsim) + stopifnot(is.numeric(alpha)) + stopifnot(alpha>=0) + if(nsim>0) testIC(nsim,alpha) + +###faire test sur les marks + + if("rectangle"%in%p$window$type) { + cas<-1 + xmin<-p$window$xmin + xmax<-p$window$xmax + ymin<-p$window$ymin + ymax<-p$window$ymax + stopifnot(upto<=(0.5*max((xmax-xmin),(ymax-ymin)))) + if ("complex"%in%p$window$type) { + cas<-3 + tri<-p$window$triangles + nbTri<-nrow(tri) + } + } + else if("circle"%in%p$window$type) { + cas<-2 + x0<-p$window$x0 + y0<-p$window$y0 + r0<-p$window$r0 + stopifnot(upto<=r0) + if ("complex"%in%p$window$type) { + cas<-4 + tri<-p$window$triangles + nbTri<-nrow(tri) + } + } + else + stop("invalid window type") + surface<-area.swin(p$window) + intensity<-p$n/surface + nbMarks<-nlevels(p$marks) + marks<-as.numeric(p$marks) # => position du label dans levels(p$marks) + dis<-as.dist(sortmat(dis,levels(p$marks))) + HD<-suppressWarnings(divc(as.data.frame(unclass(table(p$marks))),sqrt(2*dis),scale=F)[1,1]) + HD<-HD*p$n/(p$n-1) + dis<-as.vector(dis) + +# computing Rao + if(cas==1) { #rectangle + if(nsim==0) { #without CI + res<-.C("rao_rect", + as.integer(p$n),as.double(p$x),as.double(p$y), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(tmax),as.double(by),as.integer(H0), + as.integer(nbMarks),as.integer(marks),as.double(dis),as.double(surface),as.double(HD), + gg=double(tmax),kk=double(tmax),gs=double(tmax),ks=double(tmax),erreur=integer(tmax), + PACKAGE="ads") + } + else { #with CI + res<-.C("rao_rect_ic", + as.integer(p$n),as.double(p$x),as.double(p$y),as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(tmax),as.double(by), as.integer(nsim),as.integer(H0),as.double(alpha), + as.integer(nbMarks),as.integer(marks),as.double(dis),as.double(surface),as.double(HD), + gg=double(tmax),kk=double(tmax),gs=double(tmax),ks=double(tmax),gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), + gval=double(tmax),kval=double(tmax),erreur=integer(tmax), + PACKAGE="ads") + } + } + else if(cas==2) { #circle + if(nsim==0) { #without CI + res<-.C("rao_disq", + as.integer(p$n),as.double(p$x),as.double(p$y), + as.double(x0),as.double(y0),as.double(r0), + as.integer(tmax),as.double(by),as.integer(H0), + as.integer(nbMarks),as.integer(marks),as.double(dis),as.double(surface),as.double(HD), + gg=double(tmax),kk=double(tmax),gs=double(tmax),ks=double(tmax),erreur=integer(tmax), + PACKAGE="ads") + } + else { #with CI + res<-.C("rao_disq_ic", + as.integer(p$n),as.double(p$x),as.double(p$y),as.double(x0),as.double(y0),as.double(r0), + as.integer(tmax),as.double(by), as.integer(nsim),as.integer(H0),as.double(alpha), + as.integer(nbMarks),as.integer(marks),as.double(dis),as.double(surface),as.double(HD), + gg=double(tmax),kk=double(tmax),gs=double(tmax),ks=double(tmax),gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), + gval=double(tmax),kval=double(tmax),erreur=integer(tmax), + PACKAGE="ads") + } + } + else if(cas==3) { #complex within rectangle + if(nsim==0) { #without CI + res<-.C("rao_tr_rect", + as.integer(p$n),as.double(p$x),as.double(p$y), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by),as.integer(H0), + as.integer(nbMarks),as.integer(marks),as.double(dis),as.double(surface),as.double(HD), + gg=double(tmax),kk=double(tmax),gs=double(tmax),ks=double(tmax),erreur=integer(tmax), + PACKAGE="ads") + } + else { #with CI + res<-.C("rao_tr_rect_ic", + as.integer(p$n),as.double(p$x),as.double(p$y), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by),as.integer(nsim), as.integer(H0),as.double(alpha), + as.integer(nbMarks),as.integer(marks),as.double(dis),as.double(surface),as.double(HD), + gg=double(tmax),kk=double(tmax),gs=double(tmax),ks=double(tmax),gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), + gval=double(tmax),kval=double(tmax),erreur=integer(tmax), + PACKAGE="ads") + } + } + else if(cas==4) { #complex within circle + if(nsim==0) { #without CI + res<-.C("rao_tr_disq", + as.integer(p$n),as.double(p$x),as.double(p$y), + as.double(x0),as.double(y0),as.double(r0), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by),as.integer(H0), + as.integer(nbMarks),as.integer(marks),as.double(dis),as.double(surface),as.double(HD), + gg=double(tmax),kk=double(tmax),gs=double(tmax),ks=double(tmax),erreur=integer(tmax), + PACKAGE="ads") + } + else { #with CI (not based on K12) + res<-.C("rao_tr_disq_ic", + as.integer(p$n),as.double(p$x),as.double(p$y), + as.double(x0),as.double(y0),as.double(r0), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by),as.integer(nsim), as.integer(H0),as.double(alpha), + as.integer(nbMarks),as.integer(marks),as.double(dis),as.double(surface),as.double(HD), + gg=double(tmax),kk=double(tmax),gs=double(tmax),ks=double(tmax),gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), + gval=double(tmax),kval=double(tmax),serreur=integer(tmax), + PACKAGE="ads") + } + } + if(sum(res$erreur>0)){ + message("Error in ", appendLF=F) + print(match.call()) + message("No neigbors within distance intervals:") + print(paste(by*(res$erreur[res$erreur>0]-1),"-",by*res$erreur[res$erreur>0])) + message("Increase argument 'by'") + return(res=NULL) + } + if(H0==1) { + theog<-rep(1,tmax) + theok<-rep(1,tmax) + } + if(H0==2) { + theog<-res$gs + theok<-res$ks + } + gr<-data.frame(obs=res$gg,theo=theog) + kr<-data.frame(obs=res$kk,theo=theok) + if(nsim>0) { + gr<-cbind(gr,sup=res$gic1,inf=res$gic2,pval=res$gval/(nsim+1)) + kr<-cbind(kr,sup=res$kic1,inf=res$kic2,pval=res$kval/(nsim+1)) + } + call<-match.call() + res<-list(call=call,r=r,gr=gr,kr=kr) + class(res)<-c("fads","krfun") + return(res) +} + +kdfun<-function(p,upto,by,dis,nsim=0,alpha=0.01) { +# checking for input parameters + stopifnot(inherits(p,"spp")) + stopifnot(p$type=="multivariate") + if(min(p$x)<0) + p$x<-p$x+abs(min(p$x)) + if(min(p$y)<0) + p$y<-p$y+abs(min(p$y)) + stopifnot(is.numeric(upto)) + stopifnot(upto>=1) + stopifnot(is.numeric(by)) + stopifnot(by>0) + r<-seq(by,upto,by) + tmax<-length(r) + stopifnot(inherits(dis,"dist")) + stopifnot(attr(dis,"Diag")==FALSE) + stopifnot(attr(dis,"Upper")==FALSE) + stopifnot(suppressWarnings(is.euclid(dis))) +###revoir tests sur dis + if(length(levels(p$marks))!=length(labels(dis))) { + stopifnot(all(levels(p$marks)%in%labels(dis))) +#dis<-subsetdist(dis,which(labels(dis)%in%levels(p$marks))) + dis<-subsetdist(dis,levels(p$marks)) + warning("matrix 'dis' have been subsetted to match with levels(p$marks)") + } +#else if(any(labels(dis)!=levels(p$marks))) { +# attr(dis,"Labels")<-levels(p$marks) +# warning("levels(p$marks) have been assigned to attr(dis, ''Labels'')") +# } + stopifnot(is.numeric(nsim)) + stopifnot(nsim>=0) + nsim<-testInteger(nsim) + stopifnot(is.numeric(alpha)) + stopifnot(alpha>=0) + if(nsim>0) testIC(nsim,alpha) + + surface<-area.swin(p$window) + intensity<-p$n/surface + nbMarks<-nlevels(p$marks) + marks<-as.numeric(p$marks) # => position du label dans levels(p$marks) + dis<-as.dist(sortmat(dis,levels(p$marks))) + HD<-suppressWarnings(divc(as.data.frame(unclass(table(p$marks))),sqrt(2*dis),scale=F)[1,1]) + HD<-HD*p$n/(p$n-1) + dis<-as.vector(dis) + +###faire test sur les marks + + if(nsim==0) { #without CI + res<-.C("shen", + as.integer(p$n),as.double(p$x),as.double(p$y), + as.integer(tmax),as.double(by), + as.integer(nbMarks),as.integer(marks),as.double(dis),as.double(surface),as.double(HD), + gd=double(tmax),kd=double(tmax),erreur=integer(tmax), + PACKAGE="ads") + } + else { #with CI + res<-.C("shen_ic", + as.integer(p$n),as.double(p$x),as.double(p$y), + as.integer(tmax),as.double(by), as.integer(nsim),as.double(alpha), + as.integer(nbMarks),as.integer(marks),as.double(dis),as.double(surface),as.double(HD), + gd=double(tmax),kd=double(tmax),gic1=double(tmax),gic2=double(tmax),kic1=double(tmax),kic2=double(tmax), + gval=double(tmax),kval=double(tmax),erreur=integer(tmax), + PACKAGE="ads") + } + + if(sum(res$erreur>0)){ + message("Error in ", appendLF=F) + print(match.call()) + message("No neigbors within distance intervals:") + print(paste(by*(res$erreur[res$erreur>0]-1),"-",by*res$erreur[res$erreur>0])) + message("Increase argument 'by'") + return(res=NULL) + } + gd<-data.frame(obs=res$gd,theo=rep(1,tmax)) + kd<-data.frame(obs=res$kd,theo=rep(1,tmax)) + if(nsim>0) { + gd<-cbind(gd,sup=res$gic1,inf=res$gic2,pval=res$gval/(nsim+1)) + kd<-cbind(kd,sup=res$kic1,inf=res$kic2,pval=res$kval/(nsim+1)) + } + call<-match.call() + res<-list(call=call,r=r,gd=gd,kd=kd) + class(res)<-c("fads","kdfun") + return(res) +} diff --git a/R/plot.vads.R b/R/plot.vads.R index 9f80ffb..7955715 100755 --- a/R/plot.vads.R +++ b/R/plot.vads.R @@ -28,8 +28,7 @@ plot.vads.dval<-function (x,main,opt=c("dval","cval"),select,chars=c("circles"," val<-data.frame(adjust.marks.size(val,x$window,if(!missing(maxsize)) maxsize)) def.par <- par(no.readonly = TRUE) on.exit(par(def.par)) - #if(options()$device=="windows") - # csize<-0.75*csize + if (missing(main)) main <- deparse(substitute(x)) mylayout<-layout(matrix(c(rep(1,nf),seq(2,((nf*nf)+1),1)),(nf+1),nf,byrow=TRUE)) @@ -54,14 +53,7 @@ plot.vads.dval<-function (x,main,opt=c("dval","cval"),select,chars=c("circles"," text(c(xl[1]+lms[1],xl[2]+lms[2],xl[3]+lms[3]),yl,labels=signif(lm,2),pos=4,cex=1.5) } } - #if(!is.null(main)) { - # mylayout<-layout(matrix(c(rep(1,nf),seq(2,((nf*nf)+1),1)),(nf+1),nf,byrow=TRUE)) - # plot.default(x$xy$x,x$xy$y,type="n",axes=FALSE) - # text(mean(range(x$xy$x)),mean(range(x$xy$y)),pos=3,cex=2,labels=main) - #ajouter une lŽgende ??? - #} - #else - # mylayout<-layout(matrix(seq(1,(nf*nf),1),nf,nf,byrow=TRUE)) + ifelse(missing(cols),cols<-1,cols<-cols[1]) if(!missing(char0)||!missing(col0)) { ifelse(missing(col0),col0<-cols,col0<-col0[1]) @@ -82,7 +74,6 @@ plot.vads.dval<-function (x,main,opt=c("dval","cval"),select,chars=c("circles"," fg=cols,bg=cols,inches=FALSE,add=TRUE,...) } } - ## mŽthode en courbes de niveaux ? } plot.vads.kval<-function (x,main,opt=c("lval","kval","nval","gval"),select,chars=c("circles","squares"),cols,maxsize,char0,col0,legend=TRUE,csize=1,...) { @@ -117,14 +108,12 @@ plot.vads.kval<-function (x,main,opt=c("lval","kval","nval","gval"),select,chars else stopifnot(opt%in%c("lval","kval","nval","gval")) v<-val - #val<-data.frame(adjust.marks.size(val,x$window,if(!missing(maxsize)) maxsize)) val<-data.frame(adjust.marks.size(val,x$window)) if(!missing(maxsize)) val<-val*maxsize def.par <- par(no.readonly = TRUE) on.exit(par(def.par)) - #if(options()$device=="windows") - # csize<-0.75*csize + if (missing(main)) main <- deparse(substitute(x)) mylayout<-layout(matrix(c(rep(1,nf),seq(2,((nf*nf)+1),1)),(nf+1),nf,byrow=TRUE)) @@ -155,14 +144,7 @@ plot.vads.kval<-function (x,main,opt=c("lval","kval","nval","gval"),select,chars } } - #if(!is.null(main)) { - # mylayout<-layout(matrix(c(rep(1,nf),seq(2,((nf*nf)+1),1)),(nf+1),nf,byrow=TRUE)) - # plot.default(x$xy$x,x$xy$y,type="n",axes=FALSE) - # text(mean(range(x$xy$x)),mean(range(x$xy$y)),pos=3,cex=2,labels=main) - #ajouter une lŽgende ??? - #} - #else - # mylayout<-layout(matrix(seq(1,(nf*nf),1),nf,nf,byrow=TRUE)) + ifelse(missing(cols),cols<-1,cols<-cols[1]) if(!missing(char0)||!missing(col0)) { ifelse(missing(col0),col0<-cols,col0<-col0[1]) @@ -194,7 +176,6 @@ plot.vads.kval<-function (x,main,opt=c("lval","kval","nval","gval"),select,chars } } } - ## mŽthode en courbes de niveaux ? } plot.vads.k12val<-function (x,main,opt=c("lval","kval","nval","gval"),select,chars=c("circles","squares"),cols,maxsize,char0,col0,legend=TRUE,csize=1,...) { @@ -235,8 +216,7 @@ plot.vads.k12val<-function (x,main,opt=c("lval","kval","nval","gval"),select,cha val<-val*maxsize def.par <- par(no.readonly = TRUE) on.exit(par(def.par)) - #if(options()$device=="windows") - # csize<-0.75*csize + if (missing(main)) main <- deparse(substitute(x)) mylayout<-layout(matrix(c(rep(1,nf),seq(2,((nf*nf)+1),1)),(nf+1),nf,byrow=TRUE)) @@ -265,14 +245,7 @@ plot.vads.k12val<-function (x,main,opt=c("lval","kval","nval","gval"),select,cha text(c(xl[1]+lms[1],xl[2]+lms[2],xl[3]+lms[3]),yl*0.5,labels=signif(-lm,2),pos=4,cex=1) } } - #if(!is.null(main)) { - # mylayout<-layout(matrix(c(rep(1,nf),seq(2,((nf*nf)+1),1)),(nf+1),nf,byrow=TRUE)) - # plot.default(x$xy$x,x$xy$y,type="n",axes=FALSE) - # text(mean(range(x$xy$x)),mean(range(x$xy$y)),pos=3,cex=2,labels=main) - #ajouter une lŽgende ??? - #} - #else - # mylayout<-layout(matrix(seq(1,(nf*nf),1),nf,nf,byrow=TRUE)) + ifelse(missing(cols),cols<-1,cols<-cols[1]) if(!missing(char0)||!missing(col0)) { ifelse(missing(col0),col0<-cols,col0<-col0[1]) @@ -304,5 +277,4 @@ plot.vads.k12val<-function (x,main,opt=c("lval","kval","nval","gval"),select,cha } } } - ## mŽthode en courbes de niveaux ? } diff --git a/R/vads.R b/R/vads.R index a26c740..0b8cdca 100755 --- a/R/vads.R +++ b/R/vads.R @@ -1,305 +1,296 @@ -dval<-function(p,upto,by,nx,ny) { -#si multivariŽ, choix du type de points ??? - stopifnot(inherits(p,"spp")) - stopifnot(is.numeric(upto)) - stopifnot(is.numeric(by)) - stopifnot(by>0) - r<-seq(by,upto,by) - tmax<-length(r) - stopifnot(is.numeric(nx)) - stopifnot(nx>=1) - nx<-testInteger(nx) - stopifnot(is.numeric(ny)) - stopifnot(ny>=1) - ny<-testInteger(ny) - if("rectangle"%in%p$window$type) { - cas<-1 - xmin<-p$window$xmin - xmax<-p$window$xmax - ymin<-p$window$ymin - ymax<-p$window$ymax - stopifnot(upto<=(0.5*max((xmax-xmin),(ymax-ymin)))) - xsample<-rep(xmin+(seq(1,nx)-0.5)*(xmax-xmin)/nx,each=ny) - ysample<-rep(ymin+(seq(1,ny)-0.5)*(ymax-ymin)/ny,nx) - if ("complex"%in%p$window$type) { - cas<-3 - tri<-p$window$triangles - nbTri<-nrow(tri) - } - } - else if("circle"%in%p$window$type) { - cas<-2 - x0<-p$window$x0 - y0<-p$window$y0 - r0<-p$window$r0 - stopifnot(upto<=r0) - xsample<-rep(x0-r0+(seq(1,nx)-0.5)*2*r0/nx,each=ny) - ysample<-rep(y0-r0+(seq(1,ny)-0.5)*2*r0/ny,nx) - if ("complex"%in%p$window$type) { - cas<-4 - tri<-p$window$triangles - nbTri<-nrow(tri) - } - } - else - stop("invalid window type") - - ok <- inside.swin(xsample, ysample, p$window) - xsample<-xsample[ok] - ysample<-ysample[ok] - stopifnot(length(xsample)==length(ysample)) - nbSample<-length(xsample) - - if(cas==1) { #rectangle - count<-.C("density_rect", - as.integer(p$n),as.double(p$x),as.double(p$y), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(tmax),as.double(by), - as.double(xsample),as.double(ysample),as.integer(nbSample), - count=double(tmax*nbSample), - PACKAGE="ads")$count - } - else if(cas==2) { #circle - count<-.C("density_disq", - as.integer(p$n),as.double(p$x),as.double(p$y), - as.double(x0),as.double(y0),as.double(r0), - as.integer(tmax),as.double(by), - as.double(xsample),as.double(ysample),as.integer(nbSample), - count=double(tmax*nbSample), - PACKAGE="ads")$count - } - else if(cas==3) { #complex within rectangle - count<-.C("density_tr_rect", - as.integer(p$n),as.double(p$x),as.double(p$y), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), - as.double(xsample),as.double(ysample),as.integer(nbSample), - count=double(tmax*nbSample), - PACKAGE="ads")$count - } - else if(cas==4) { #complex within circle - count<-.C("density_tr_disq", - as.integer(p$n),as.double(p$x),as.double(p$y), - as.double(x0),as.double(y0),as.double(r0), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), - as.double(xsample),as.double(ysample),as.integer(nbSample), - count=double(tmax*nbSample), - PACKAGE="ads")$count - } - ## rajouter un indice lorsque les disques ne sont pas indŽpendants - # formatting results - dens<-count/(pi*r^2) - #grid<-matrix(c(xsample,ysample),nrow=nbSample,ncol=2) - count<-matrix(count,nrow=nbSample,ncol=tmax,byrow=TRUE) - dens<-matrix(dens,nrow=nbSample,ncol=tmax,byrow=TRUE) - call<-match.call() - res<-list(call=call,window=p$window,r=r,xy=data.frame(x=xsample,y=ysample),cval=count,dval=dens) - class(res)<-c("vads","dval") - return(res) -} - -kval<-function(p,upto,by) { - # checking for input parameters - stopifnot(inherits(p,"spp")) - if(p$type!="univariate") - warning(paste(p$type,"point pattern has been considered to be univariate\n")) - stopifnot(is.numeric(upto)) - stopifnot(is.numeric(by)) - stopifnot(by>0) - r<-seq(by,upto,by) - tmax<-length(r) - - if("rectangle"%in%p$window$type) { - cas<-1 - xmin<-p$window$xmin - xmax<-p$window$xmax - ymin<-p$window$ymin - ymax<-p$window$ymax - stopifnot(upto<=(0.5*max((xmax-xmin),(ymax-ymin)))) - if ("complex"%in%p$window$type) { - cas<-3 - tri<-p$window$triangles - nbTri<-nrow(tri) - } - } - else if("circle"%in%p$window$type) { - cas<-2 - x0<-p$window$x0 - y0<-p$window$y0 - r0<-p$window$r0 - stopifnot(upto<=r0) - if ("complex"%in%p$window$type) { - cas<-4 - tri<-p$window$triangles - nbTri<-nrow(tri) - } - } - else - stop("invalid window type") - intensity<-p$n/area.swin(p$window) - - #computing ripley local functions - if(cas==1) { #rectangle - res<-.C("ripleylocal_rect", - as.integer(p$n),as.double(p$x),as.double(p$y), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(tmax),as.double(by), - gi=double(p$n*tmax),ki=double(p$n*tmax), - PACKAGE="ads") - } - else if(cas==2) { #circle - res<-.C("ripleylocal_disq", - as.integer(p$n),as.double(p$x),as.double(p$y), - as.double(x0),as.double(y0),as.double(r0), - as.integer(tmax),as.double(by), - gi=double(p$n*tmax),ki=double(p$n*tmax), - PACKAGE="ads") - } - else if(cas==3) { #complex within rectangle - res<-.C("ripleylocal_tr_rect", - as.integer(p$n),as.double(p$x),as.double(p$y), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), - gi=double(p$n*tmax),ki=double(p$n*tmax), - PACKAGE="ads") - } - else if(cas==4) { #complex within circle - res<-.C("ripleylocal_tr_disq", - as.integer(p$n),as.double(p$x),as.double(p$y), - as.double(x0),as.double(y0),as.double(r0), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), - gi=double(p$n*tmax),ki=double(p$n*tmax), - PACKAGE="ads") - } - # formatting results - #coord<-matrix(c(X$x,X$y),nrow=nbPts,ncol=2) - #coord<-data.frame(x=p$x,y=p$y) - #r<-seq(dr,dr*tmax,dr) - #ds<-pi*r^2-pi*seq(0,dr*tmax-dr,dr)^2 - ds<-c(pi,diff(pi*r^2)) - gi<-matrix(res$gi/(intensity*ds),nrow=p$n,ncol=tmax,byrow=TRUE) - ni<-matrix(res$ki/(pi*r^2),nrow=p$n,ncol=tmax,byrow=TRUE) - ki<-matrix(res$ki/intensity,nrow=p$n,ncol=tmax,byrow=TRUE) - li<-matrix(sqrt(res$ki/(intensity*pi))-r,nrow=p$n,ncol=tmax,byrow=TRUE) - call<-match.call() - res<-list(call=call,window=p$window,r=r,xy=data.frame(x=p$x,y=p$y),gval=gi,kval=ki,nval=ni,lval=li) - class(res)<-c("vads","kval") - return(res) -} - -k12val<-function(p,upto,by,marks) { - # checking for input parameters - stopifnot(inherits(p,"spp")) - stopifnot(p$type=="multivariate") - stopifnot(is.numeric(upto)) - stopifnot(is.numeric(by)) - stopifnot(by>0) - r<-seq(by,upto,by) - tmax<-length(r) - if(missing(marks)) - marks<-c(1,2) - stopifnot(length(marks)==2) - stopifnot(marks[1]!=marks[2]) - mark1<-marks[1] - mark2<-marks[2] - if(is.numeric(mark1)) - mark1<-levels(p$marks)[testInteger(mark1)] - else if(!mark1%in%p$marks) stop(paste("mark \'",mark1,"\' doesn\'t exist",sep="")) - if(is.numeric(mark2)) - mark2<-levels(p$marks)[testInteger(mark2)] - else if(!mark2%in%p$marks) stop(paste("mark \'",mark2,"\' doesn\'t exist",sep="")) - # initializing variables - if("rectangle"%in%p$window$type) { - cas<-1 - xmin<-p$window$xmin - xmax<-p$window$xmax - ymin<-p$window$ymin - ymax<-p$window$ymax - stopifnot(upto<=(0.5*max((xmax-xmin),(ymax-ymin)))) - if ("complex"%in%p$window$type) { - cas<-3 - tri<-p$window$triangles - nbTri<-nrow(tri) - } - } - else if("circle"%in%p$window$type) { - cas<-2 - x0<-p$window$x0 - y0<-p$window$y0 - r0<-p$window$r0 - stopifnot(upto<=r0) - if ("complex"%in%p$window$type) { - cas<-4 - tri<-p$window$triangles - nbTri<-nrow(tri) - } - } - else - stop("invalid window type") - surface<-area.swin(p$window) - x1<-p$x[p$marks==mark1] - y1<-p$y[p$marks==mark1] - x2<-p$x[p$marks==mark2] - y2<-p$y[p$marks==mark2] - nbPts1<-length(x1) - nbPts2<-length(x2) - intensity2<-nbPts2/surface - #computing intertype local functions - if(cas==1) { #rectangle - res<-.C("intertypelocal_rect", - as.integer(nbPts1),as.double(x1),as.double(y1), - as.integer(nbPts2),as.double(x2),as.double(y2), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(tmax),as.double(by), - gi=double(nbPts1*tmax),ki=double(nbPts1*tmax), - PACKAGE="ads") - } - else if(cas==2) { #circle - res<-.C("intertypelocal_disq", - as.integer(nbPts1),as.double(x1),as.double(y1), - as.integer(nbPts2),as.double(x2),as.double(y2), - as.double(x0),as.double(y0),as.double(r0), - as.integer(tmax),as.double(by), - gi=double(nbPts1*tmax),ki=double(nbPts1*tmax), - PACKAGE="ads") - } - else if(cas==3) { #complex within rectangle - res<-.C("intertypelocal_tr_rect", - as.integer(nbPts1),as.double(x1),as.double(y1), - as.integer(nbPts2),as.double(x2),as.double(y2), - as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), - gi=double(nbPts1*tmax),ki=double(nbPts1*tmax), - PACKAGE="ads") - } - else if(cas==4) { #complex within circle - res<-.C("intertypelocal_tr_disq", - as.integer(nbPts1),as.double(x1),as.double(y1), - as.integer(nbPts2),as.double(x2),as.double(y2), - as.double(x0),as.double(y0),as.double(r0), - as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), - as.integer(tmax),as.double(by), - gi=double(nbPts1*tmax),ki=double(nbPts1*tmax), - PACKAGE="ads") - } - # formatting results - #coord<-matrix(c(x1,y1),nrow=nbPts1,ncol=2) - #coord<-data.frame(x1=x1,y1=y1) - #r<-seq(dr,dr*tmax,dr) - #ds<-pi*r^2-pi*seq(0,dr*tmax-dr,dr)^2 - ds<-c(pi,diff(pi*r^2)) - gi<-matrix(res$gi/(intensity2*ds),nrow=nbPts1,ncol=tmax,byrow=TRUE) - ni<-matrix(res$ki/(pi*r^2),nrow=nbPts1,ncol=tmax,byrow=TRUE) - ki<-matrix(res$ki/intensity2,nrow=nbPts1,ncol=tmax,byrow=TRUE) - li<-matrix(sqrt(res$ki/(intensity2*pi))-r,nrow=nbPts1,ncol=tmax,byrow=TRUE) - call<-match.call() - res<-list(call=call,window=p$window,r=r,xy=data.frame(x=x1,y=y1),g12val=gi,k12val=ki,n12val=ni,l12val=li,marks=c(mark1,mark2)) - class(res)<-c("vads","k12val") - return(res) -} - +dval<-function(p,upto,by,nx,ny) { + stopifnot(inherits(p,"spp")) + stopifnot(is.numeric(upto)) + stopifnot(is.numeric(by)) + stopifnot(by>0) + r<-seq(by,upto,by) + tmax<-length(r) + stopifnot(is.numeric(nx)) + stopifnot(nx>=1) + nx<-testInteger(nx) + stopifnot(is.numeric(ny)) + stopifnot(ny>=1) + ny<-testInteger(ny) + if("rectangle"%in%p$window$type) { + cas<-1 + xmin<-p$window$xmin + xmax<-p$window$xmax + ymin<-p$window$ymin + ymax<-p$window$ymax + stopifnot(upto<=(0.5*max((xmax-xmin),(ymax-ymin)))) + xsample<-rep(xmin+(seq(1,nx)-0.5)*(xmax-xmin)/nx,each=ny) + ysample<-rep(ymin+(seq(1,ny)-0.5)*(ymax-ymin)/ny,nx) + if ("complex"%in%p$window$type) { + cas<-3 + tri<-p$window$triangles + nbTri<-nrow(tri) + } + } + else if("circle"%in%p$window$type) { + cas<-2 + x0<-p$window$x0 + y0<-p$window$y0 + r0<-p$window$r0 + stopifnot(upto<=r0) + xsample<-rep(x0-r0+(seq(1,nx)-0.5)*2*r0/nx,each=ny) + ysample<-rep(y0-r0+(seq(1,ny)-0.5)*2*r0/ny,nx) + if ("complex"%in%p$window$type) { + cas<-4 + tri<-p$window$triangles + nbTri<-nrow(tri) + } + } + else + stop("invalid window type") + + ok <- inside.swin(xsample, ysample, p$window) + xsample<-xsample[ok] + ysample<-ysample[ok] + stopifnot(length(xsample)==length(ysample)) + nbSample<-length(xsample) + + if(cas==1) { #rectangle + count<-.C("density_rect", + as.integer(p$n),as.double(p$x),as.double(p$y), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(tmax),as.double(by), + as.double(xsample),as.double(ysample),as.integer(nbSample), + count=double(tmax*nbSample), + PACKAGE="ads")$count + } + else if(cas==2) { #circle + count<-.C("density_disq", + as.integer(p$n),as.double(p$x),as.double(p$y), + as.double(x0),as.double(y0),as.double(r0), + as.integer(tmax),as.double(by), + as.double(xsample),as.double(ysample),as.integer(nbSample), + count=double(tmax*nbSample), + PACKAGE="ads")$count + } + else if(cas==3) { #complex within rectangle + count<-.C("density_tr_rect", + as.integer(p$n),as.double(p$x),as.double(p$y), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), + as.double(xsample),as.double(ysample),as.integer(nbSample), + count=double(tmax*nbSample), + PACKAGE="ads")$count + } + else if(cas==4) { #complex within circle + count<-.C("density_tr_disq", + as.integer(p$n),as.double(p$x),as.double(p$y), + as.double(x0),as.double(y0),as.double(r0), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), + as.double(xsample),as.double(ysample),as.integer(nbSample), + count=double(tmax*nbSample), + PACKAGE="ads")$count + } + ## rajouter un indice lorsque les disques ne sont pas independants + # formatting results + dens<-count/(pi*r^2) + #grid<-matrix(c(xsample,ysample),nrow=nbSample,ncol=2) + count<-matrix(count,nrow=nbSample,ncol=tmax,byrow=TRUE) + dens<-matrix(dens,nrow=nbSample,ncol=tmax,byrow=TRUE) + call<-match.call() + res<-list(call=call,window=p$window,r=r,xy=data.frame(x=xsample,y=ysample),cval=count,dval=dens) + class(res)<-c("vads","dval") + return(res) +} + +kval<-function(p,upto,by) { + # checking for input parameters + stopifnot(inherits(p,"spp")) + if(p$type!="univariate") + warning(paste(p$type,"point pattern has been considered to be univariate\n")) + stopifnot(is.numeric(upto)) + stopifnot(is.numeric(by)) + stopifnot(by>0) + r<-seq(by,upto,by) + tmax<-length(r) + + if("rectangle"%in%p$window$type) { + cas<-1 + xmin<-p$window$xmin + xmax<-p$window$xmax + ymin<-p$window$ymin + ymax<-p$window$ymax + stopifnot(upto<=(0.5*max((xmax-xmin),(ymax-ymin)))) + if ("complex"%in%p$window$type) { + cas<-3 + tri<-p$window$triangles + nbTri<-nrow(tri) + } + } + else if("circle"%in%p$window$type) { + cas<-2 + x0<-p$window$x0 + y0<-p$window$y0 + r0<-p$window$r0 + stopifnot(upto<=r0) + if ("complex"%in%p$window$type) { + cas<-4 + tri<-p$window$triangles + nbTri<-nrow(tri) + } + } + else + stop("invalid window type") + intensity<-p$n/area.swin(p$window) + + #computing ripley local functions + if(cas==1) { #rectangle + res<-.C("ripleylocal_rect", + as.integer(p$n),as.double(p$x),as.double(p$y), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(tmax),as.double(by), + gi=double(p$n*tmax),ki=double(p$n*tmax), + PACKAGE="ads") + } + else if(cas==2) { #circle + res<-.C("ripleylocal_disq", + as.integer(p$n),as.double(p$x),as.double(p$y), + as.double(x0),as.double(y0),as.double(r0), + as.integer(tmax),as.double(by), + gi=double(p$n*tmax),ki=double(p$n*tmax), + PACKAGE="ads") + } + else if(cas==3) { #complex within rectangle + res<-.C("ripleylocal_tr_rect", + as.integer(p$n),as.double(p$x),as.double(p$y), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), + gi=double(p$n*tmax),ki=double(p$n*tmax), + PACKAGE="ads") + } + else if(cas==4) { #complex within circle + res<-.C("ripleylocal_tr_disq", + as.integer(p$n),as.double(p$x),as.double(p$y), + as.double(x0),as.double(y0),as.double(r0), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), + gi=double(p$n*tmax),ki=double(p$n*tmax), + PACKAGE="ads") + } + # formatting results + ds<-c(pi,diff(pi*r^2)) + gi<-matrix(res$gi/(intensity*ds),nrow=p$n,ncol=tmax,byrow=TRUE) + ni<-matrix(res$ki/(pi*r^2),nrow=p$n,ncol=tmax,byrow=TRUE) + ki<-matrix(res$ki/intensity,nrow=p$n,ncol=tmax,byrow=TRUE) + li<-matrix(sqrt(res$ki/(intensity*pi))-r,nrow=p$n,ncol=tmax,byrow=TRUE) + call<-match.call() + res<-list(call=call,window=p$window,r=r,xy=data.frame(x=p$x,y=p$y),gval=gi,kval=ki,nval=ni,lval=li) + class(res)<-c("vads","kval") + return(res) +} + +k12val<-function(p,upto,by,marks) { + # checking for input parameters + stopifnot(inherits(p,"spp")) + stopifnot(p$type=="multivariate") + stopifnot(is.numeric(upto)) + stopifnot(is.numeric(by)) + stopifnot(by>0) + r<-seq(by,upto,by) + tmax<-length(r) + if(missing(marks)) + marks<-c(1,2) + stopifnot(length(marks)==2) + stopifnot(marks[1]!=marks[2]) + mark1<-marks[1] + mark2<-marks[2] + if(is.numeric(mark1)) + mark1<-levels(p$marks)[testInteger(mark1)] + else if(!mark1%in%p$marks) stop(paste("mark \'",mark1,"\' doesn\'t exist",sep="")) + if(is.numeric(mark2)) + mark2<-levels(p$marks)[testInteger(mark2)] + else if(!mark2%in%p$marks) stop(paste("mark \'",mark2,"\' doesn\'t exist",sep="")) + # initializing variables + if("rectangle"%in%p$window$type) { + cas<-1 + xmin<-p$window$xmin + xmax<-p$window$xmax + ymin<-p$window$ymin + ymax<-p$window$ymax + stopifnot(upto<=(0.5*max((xmax-xmin),(ymax-ymin)))) + if ("complex"%in%p$window$type) { + cas<-3 + tri<-p$window$triangles + nbTri<-nrow(tri) + } + } + else if("circle"%in%p$window$type) { + cas<-2 + x0<-p$window$x0 + y0<-p$window$y0 + r0<-p$window$r0 + stopifnot(upto<=r0) + if ("complex"%in%p$window$type) { + cas<-4 + tri<-p$window$triangles + nbTri<-nrow(tri) + } + } + else + stop("invalid window type") + surface<-area.swin(p$window) + x1<-p$x[p$marks==mark1] + y1<-p$y[p$marks==mark1] + x2<-p$x[p$marks==mark2] + y2<-p$y[p$marks==mark2] + nbPts1<-length(x1) + nbPts2<-length(x2) + intensity2<-nbPts2/surface + #computing intertype local functions + if(cas==1) { #rectangle + res<-.C("intertypelocal_rect", + as.integer(nbPts1),as.double(x1),as.double(y1), + as.integer(nbPts2),as.double(x2),as.double(y2), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(tmax),as.double(by), + gi=double(nbPts1*tmax),ki=double(nbPts1*tmax), + PACKAGE="ads") + } + else if(cas==2) { #circle + res<-.C("intertypelocal_disq", + as.integer(nbPts1),as.double(x1),as.double(y1), + as.integer(nbPts2),as.double(x2),as.double(y2), + as.double(x0),as.double(y0),as.double(r0), + as.integer(tmax),as.double(by), + gi=double(nbPts1*tmax),ki=double(nbPts1*tmax), + PACKAGE="ads") + } + else if(cas==3) { #complex within rectangle + res<-.C("intertypelocal_tr_rect", + as.integer(nbPts1),as.double(x1),as.double(y1), + as.integer(nbPts2),as.double(x2),as.double(y2), + as.double(xmin),as.double(xmax),as.double(ymin),as.double(ymax), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), + gi=double(nbPts1*tmax),ki=double(nbPts1*tmax), + PACKAGE="ads") + } + else if(cas==4) { #complex within circle + res<-.C("intertypelocal_tr_disq", + as.integer(nbPts1),as.double(x1),as.double(y1), + as.integer(nbPts2),as.double(x2),as.double(y2), + as.double(x0),as.double(y0),as.double(r0), + as.integer(nbTri),as.double(tri$ax),as.double(tri$ay),as.double(tri$bx),as.double(tri$by),as.double(tri$cx),as.double(tri$cy), + as.integer(tmax),as.double(by), + gi=double(nbPts1*tmax),ki=double(nbPts1*tmax), + PACKAGE="ads") + } + # formatting results + ds<-c(pi,diff(pi*r^2)) + gi<-matrix(res$gi/(intensity2*ds),nrow=nbPts1,ncol=tmax,byrow=TRUE) + ni<-matrix(res$ki/(pi*r^2),nrow=nbPts1,ncol=tmax,byrow=TRUE) + ki<-matrix(res$ki/intensity2,nrow=nbPts1,ncol=tmax,byrow=TRUE) + li<-matrix(sqrt(res$ki/(intensity2*pi))-r,nrow=nbPts1,ncol=tmax,byrow=TRUE) + call<-match.call() + res<-list(call=call,window=p$window,r=r,xy=data.frame(x=x1,y=y1),g12val=gi,k12val=ki,n12val=ni,l12val=li,marks=c(mark1,mark2)) + class(res)<-c("vads","k12val") + return(res) +} + diff --git a/src/Zlibs.h b/src/Zlibs.h index 9130b52..920b269 100755 --- a/src/Zlibs.h +++ b/src/Zlibs.h @@ -1,175 +1,87 @@ -double un_point(double,double,double,double,double,double,double,double,double); +int corr_disq_ic(int*,double*,double*,double*,double*,double*,double*,int*,double*,int*,double*,double*,double*,double*,double*,double*,double*,double*,double*); +int corr_disq(int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*); +int corr_rect_ic(int*,double*,double*,double*,double*,double*,double*,double*,int*,double*,int*,double*,double*,double*,double*,double*,double*,double*,double*,double*); +int corr_rect(int*,double*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*); +int corr_tr_disq_ic(int*,double*x,double*,double*,double*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,int*,double*,double*,double*,double*,double*,double*,double*,double*,double*); +int corr_tr_disq(int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*); +int corr_tr_rect_ic(int*,double*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,int*,double*,double*,double*,double*,double*,double*,double*,double*,double*); +int corr_tr_rect(int*,double*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*); +int density_disq(int*,double*,double*,double*,double*,double*,int*,double*,double*,double*,int*,double*); +int density_rect(int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*,int*,double*); +int density_tr_disq(int*,double*,double*,double*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*,int*,double*); +int density_tr_rect(int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*,int*,double*); double deux_point(double,double,double,double,double,double,double,double,double); -double ununun_point(double,double,double,double,double,double,double,double,double); -double trois_point(double,double,double,double,double,double,double,double,double); -double deuxun_point(double,double,double,double,double,double,double,double,double); double deuxbord_point(double,double,double,double,double,double,double,double,double); - +double deuxun_point(double,double,double,double,double,double,double,double,double); +double echange_point_disq(int,double*,double*,double,double,double,double,double,double,double*,int*,double*,double*,double*); +double echange_point_rect(int,double*,double*,double,double,double,double,double,double,double,double*,int*,double*,double*,double*); +double echange_point_tr_disq(int,double*,double*,double,double,double,int*,double*,double*,double*,double*,double*,double*,double,double,double,double*,int*,double*,double*,double*); +double echange_point_tr_rect(int,double*,double*,double,double,double,double,int*,double*,double*,double*,double*,double*,double*,double,double,double,double*,int*,double*,double*,double*); +void ic(int,int,double**,double**,double*,double*,int); int in_droite(double,double,double,double,double,double,double,double,int); int in_triangle(double,double,double,double,double,double,double,double,int); - -void ic(int,int,double **,double **,double *,double *,int); - -double perim_in_rect(double, double, double, double, double, double, double); +int intertype_disq_ic(int*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,int*,int*,double*,int*,int*,int*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*); +int intertype_disq(int*,double*,double*,int*,double*,double*,double*,double*,double*,int*,double*,double*,double*); +int intertype_rect_ic(int*,double*,double*,int*,double*,double*,double*,double*,double*,double*,double*,int*,double*,int*,int*,double*,int*,int*,int*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*); +int intertype_rect(int*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*); +int intertype_tr_disq_ic(int*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,int*,int*,double*,int*,int*,int*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*); +int intertype_tr_disq(int*,double*,double*,int*,double*,double*,double*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*); +int intertype_tr_rect_ic(int*,double*,double*,int*,double*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,int*,int*,double*,int*,int*,int*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*); +int intertype_tr_rect(int*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*); +int intertype(int*,double*,double*,int*,double*,double*,int*,double*,double*,double*); +int intertypelocal_disq(int*,double*,double*,int*,double*,double*,double*,double*,double*,int*,double*,double*,double*); +int intertypelocal_rect(int*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*); +int intertypelocal_tr_disq(int*,double*,double*,int*,double*,double*,double*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*); +int intertypelocal_tr_rect(int*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*); +int mimetic_disq(int*,double*,double*,double*,double*,double*,double*,double*,int*,double*,double*,int*,int*,double*,double*,double*,double*,double*,int*); +int mimetic_rect(int*,double*,double*,double*,double*,double*,double*,double*,double*,int*,double*,double*,int*,int*,double*,double*,double*,double*,double*,int*); +int mimetic_tr_disq(int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*,double*,double*,double*cy,double*,int*,double*,double*,int*,int*,double*,double*,double*,double*,double*,int*); +int mimetic_tr_rect(int*,double*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*,double*,double*,double*,double*,int*,double*,double*,int*,int*,double*,double*,double*,double*,double*,int*); double perim_in_disq(double,double,double,double,double,double); -double perim_triangle(double,double,double,int,double *,double *,double *,double *,double *,double *); - +double perim_in_rect(double,double,double,double,double,double,double); +double perim_triangle(double,double,double,int,double*,double*,double*,double*,double*,double*); +int randlabelling(double*,double*,int,double*,double*,int,double*,double*,int*); +void randmark(int ,double*,double*); +int randomdist(int*,int,double*,double*); +int randomlab(double*,double*,int,int*,int,double**,int*,double**); +int randshifting_disq(int*,double*,double*,int,double*,double*,double,double,double,double); +int randshifting_rect(int*,double*,double*,int,double*,double*,double,double,double,double,double); +int randshifting_tr_disq(int*,double*,double*,int,double*,double*,double,double,double,int,double*,double*,double*,double*,double*,double*,double); +int randshifting_tr_rect(int*,double*,double*,int,double*,double*,double,double,double,double,int,double*,double*,double*,double*,double*,double*,double); +int rao_disq_ic(int*,double*,double*,double*,double*,double*,int*,double*,int*,int*,double*,int*,int*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,int*); +int rao_disq(int*,double*,double*,double*,double*,double*,int*,double*,int*,int*,int*,double*,double*,double*,double*,double*,double*,double*,int*); +int rao_rect_ic(int*,double*,double*,double*,double*,double*,double*,int*,double*,int*,int*,double*,int*,int*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,int*); +int rao_rect(int*,double*,double*,double*,double*,double*,double*,int*,double*,int*,int*,int*,double*,double*,double*,double*,double*,double*,double*,int*); +int rao_tr_disq_ic(int*,double*,double*,double*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,int*,int*,double*,int*,int*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,int*); +int rao_tr_disq(int*,double*,double*,double*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,int*,int*,int*,double*,double*,double*,double*,double*,double*,double*,int*); +int rao_tr_rect_ic(int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,int*,int*,double*,int*,int*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,int*); +int rao_tr_rect(int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,int*,int*,int*,double*,double*,double*,double*,double*,double*,double*,int*); +int ripley_disq_ic(int*,double*,double*,double*,double*,double*,double*,int*,double*,int*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*); +int ripley_disq(int*,double*,double*,double*,double*,double*,int*,double*,double*,double*); +int ripley_rect_ic(int*,double*,double*,double*,double*,double*,double*,double*,int*,double*,int*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*); int ripley_rect(int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*); -int ripley_disq(int *,double *,double *,double *,double *,double *,int *,double *,double *,double *); -int ripley_tr_rect(int *,double *,double *,double *,double *,double *,double *,int *,double *,double *, - double *,double *,double *,double *,int *,double *,double *,double *); -int ripley_tr_disq(int *,double *,double *,double *,double *,double *,int *,double *,double *,double *,double *, - double *,double *,int *,double *,double *,double *); - -int ripley_rect_ic(int *,double *,double *,double *,double *,double *,double *,double *,int *,double *,int *, - double *,double *,double *,double *,double *,double *,double *,double *,double *,double *,double *,double *); -int ripley_disq_ic(int *,double *,double *,double *,double *,double *,double *,int *,double *,int *,double *, - double *,double *,double *,double *,double *,double *,double *,double *,double *,double *,double *); -int ripley_tr_rect_ic(int *,double *,double *,double *,double *,double *,double *,double *,int *, double *, - double *,double *,double *,double *,double *,int *,double *,int *,double *,double *,double *,double *, - double *,double *,double *,double *,double *,double *,double *,double *); -int ripley_tr_disq_ic(int *,double *,double *,double *,double *,double *,double *,int *, double *, double *, - double *,double *,double *,double *,int *,double *,int *,double *,double *,double *,double *,double *, - double *,double *,double *,double *,double *,double *,double *); - -int ripleylocal_rect(int*,double *,double *,double*,double*,double*,double*,int*,double*,double *,double *); -int ripleylocal_disq(int *,double *,double *,double *,double *,double *,int *,double *,double *,double *); -int ripleylocal_tr_rect(int *,double *,double *,double *,double *,double *,double *,int *, - double *,double *,double *,double *,double *,double *,int *,double *,double *,double *); -int ripleylocal_tr_disq(int *,double *,double *,double *,double *,double *,int *, - double *,double *,double *,double *,double *,double *,int *,double *,double *,double *); - -int density_rect(int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*,int*,double*); -int density_disq(int*,double *,double *,double*,double*,double*,int*,double*,double *,double *,int*,double *); -int density_tr_rect(int*,double *,double *,double*,double*,double*,double*,int*,double *,double *,double *, - double *,double *,double *,int*,double*,double *,double *,int*,double *); -int density_tr_disq(int*,double *,double *,double*,double*,double*,int*,double *,double *,double *,double *, - double *,double *,int*,double*,double *,double *,int*,double *); - -int intertype_rect(int *,double *,double *,int *,double *,double *,double *,double *,double *,double *, - int *,double *,double *,double *); -int intertype_disq(int *,double *,double *,int *,double *,double *,double *,double *,double *,int *, - double *,double *,double *); -int intertype_tr_rect(int *,double *,double *,int *,double *,double *,double *,double *,double *,double *,int *, - double *,double *,double *,double *,double *,double *,int *,double *,double *,double *); -int intertype_tr_disq(int *,double *,double *,int *,double *,double *,double *,double *,double *,int *, - double *,double *,double *,double *,double *,double *,int *,double *,double *,double *); - -int intertype_rect_ic(int *,double *,double *,int *,double *,double *,double *,double *,double *,double *, - double *,int *,double *,int *,int *,double *,int *, int *, int *, double *,double *,double *,double *,double *,double *,double *, - double *,double *,double *,double *); -int intertype_disq_ic(int *,double *,double *,int *,double *,double *,double *,double *,double *,double *,int *, - double *,int *,int *,double *,int *,int *,int *,double *,double *,double *,double *,double *,double *,double *,double *, - double *,double *,double *); -int intertype_tr_rect_ic(int *,double *,double *,int *,double *,double *,double *,double *,double *,double *, - double *,int *,double *,double *,double *,double *,double *,double *,int *,double *,int *,int *,double *,int *,int *,int *,double *, - double *,double *,double *,double *,double *,double *,double *,double *,double *,double *); -int intertype_tr_disq_ic(int *,double *,double *,int *, double *, double *,double *,double *,double *,double *, - int *,double *,double *,double *,double *,double *,double *,int *,double *,int *,int *,double *,int *,int *,int*,double *,double *, - double *,double *,double *,double *,double *,double *,double *,double *,double *); - -int intertypelocal_rect(int*,double *,double *,int*,double *,double *,double*,double*,double*,double*, - int*,double*,double *,double *); -int intertypelocal_disq(int *,double *,double *,int *,double *,double *,double *,double *, double *,int *, - double *,double *,double *); -int intertypelocal_tr_rect(int *,double *,double *,int *,double *,double *,double *,double *,double *,double *, - int *,double *,double *,double *,double *,double *,double *,int *,double *,double *,double *); -int intertypelocal_tr_disq(int *,double *,double *,int *,double *,double *,double *,double *,double *,int *, - double *,double *,double *,double *,double *,double *,int *,double *,double *,double *); - +int ripley_tr_disq_ic(int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,int*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*); +int ripley_tr_disq(int*,double*,double*,double*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*); +int ripley_tr_rect_ic(int*,double*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,int*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*); +int ripley_tr_rect(int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*); +int ripleylocal_disq(int*,double*,double*,double*,double*,double*,int*,double*,double*,double*); +int ripleylocal_rect(int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*); +int ripleylocal_tr_disq(int*,double*,double*,double*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*); +int ripleylocal_tr_rect(int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*); +void s_alea_disq(int,double*,double*,double,double,double,double); void s_alea_rect(int,double[],double[],double,double,double,double,double); -void s_alea_disq(int,double *,double *,double,double,double,double); -void s_alea_tr_rect(int,double *,double *,double,double,double,double,int,double *,double *,double *,double *, - double *,double *,double); -void s_alea_tr_disq(int ,double *,double *,double,double,double,int,double *,double *,double *,double *, - double *,double *,double); - -int randlabelling(double *, double *, int, double *, double *,int, double *, double *,int *); -int randshifting_rect(int *,double *,double *,int,double *,double *,double,double,double,double,double); -int randshifting_disq(int *,double *,double *,int,double *,double *,double,double,double,double); -int randshifting_tr_rect(int *,double *,double *,int,double *,double *,double,double,double,double,int, - double *,double *,double *,double *,double *,double *,double); -int randshifting_tr_disq(int *,double *,double *,int,double *,double *,double,double,double,int, - double *,double *,double *,double *,double *,double *,double); -int randomlab(double *,double *,int,int *,int,double **,int *,double **); -void randmark(int ,double *,double *); -int randomdist(int *,int,double *,double *); - -int corr_rect(int *,double *,double *,double *, double *,double *,double *,double *, -int *,double *,double *,double *); -int corr_disq(int *,double *,double *,double *, double *,double *,double *, -int *,double *,double *,double *); -int corr_tr_rect(int *,double *,double *,double *, double *,double *,double *,double *, -int *, double *, double *, double *, double *, double *, double *, -int *,double *,double *,double *); -int corr_tr_disq(int *,double *,double *,double *, double *,double *,double *, -int *,double *,double *,double *,double *,double *,double *, -int *,double *,double *,double *); - -int corr_rect_ic(int *,double *,double *,double *, double *,double *,double *,double *, -int *,double *,int *, double *,double *,double *, -double *,double *, double *,double *, double *, double *); -int corr_disq_ic(int *,double *,double *,double *, double *,double *,double *, -int *,double *,int *, double *,double *,double *, -double *,double *, double *,double *, double *, double *); -int corr_tr_rect_ic(int *,double *,double *,double *, double *,double *,double *,double *, -int *, double *, double *, double *, double *, double *, double *, -int *,double *,int *, double *,double *,double *, -double *,double *, double *,double *, double *, double *); -int corr_tr_disq_ic(int *,double *x,double *,double *, double *,double *,double *, -int *, double *, double *, double *, double *, double *, double *, -int *,double *,int *, double *,double *,double *, -double *,double *, double *,double *, double *, double *); - -int shimatani_rect(int *,double *,double *,double *,double *,double *,double *,int *,double *,int *,int *,double *,double *, double *,int *); -int shimatani_disq(int *,double *,double *, double *,double *,double *,int *,double *,int *,int *,double *,double *, double *,int *); -int shimatani_tr_rect(int *,double *,double *, double *,double *,double *,double *,int *, double *, double *, double *, double *, double *, double *, - int *,double *,int *,int *,double *,double *, double *,int *); -int shimatani_tr_disq(int *,double *,double *, double *,double *,double *,int *, double *, double *, double *, double *, double *, double *, - int *,double *,int *,int *,double *,double *, double *,int *); -int shimatani_rect_ic(int *,double *,double *, double *,double *,double *,double *,int *,double *,int *,double *,int *,int *,double *,double *, - double *, double *,double *,double *,double *,double *,double *,double *,int *); -int shimatani_disq_ic(int *,double *,double *, double *,double *,double *,int *,double *,int *,double *,int *,int *,double *,double *, - double *, double *,double *,double *,double *,double *,double *,double *,int *); -int shimatani_tr_rect_ic(int *,double *,double *, double *,double *,double *,double *,int *, double *, double *, double *, double *, double *, double *, - int *,double *,int *,double *,int *,int *,double *,double *,double *, double *,double *,double *,double *,double *,double *,double *,int *); -int shimatani_tr_disq_ic(int *,double *,double *, double *,double *,double *,int *, double *, double *, double *, double *, double *, double *, - int *,double *,int *,double *,int *,int *,double *,double *,double *, double *,double *,double *,double *,double *,double *,double *,int *); - -int rao_rect(int *,double *,double *,double *,double *,double *,double *,int *,double *,int *,int *,int *,double *,double *,double *,double *, - double *,double *,double *,int *); -int rao_rect_ic(int *,double *,double *,double *,double *,double *,double *,int *,double *,int *,int *,double *,int *,int *,double *,double *, - double *,double *,double *,double *,double *,double *,double *,double *,double *,double *,double *,int *); -int rao_disq(int *,double *,double *,double *,double *,double *,int *,double *,int *,int *,int *,double *,double *,double *,double *, - double *,double *,double *,int *); -int rao_disq_ic(int *,double *,double *,double *,double *,double *,int *,double *,int *,int *,double *,int *,int *,double *,double *, - double *,double *,double *,double *,double *,double *,double *,double *,double *,double *,double *,int *); -int rao_tr_rect(int *,double *,double *,double *,double *,double *,double *,int *,double *,double *,double *,double *,double *,double *, - int *,double *,int *,int *,int *,double *,double *,double *,double *,double *,double *,double *,int *); - -int rao_tr_rect_ic(int *,double *,double *,double *,double *,double *,double *,int *,double *,double *,double *,double *,double *,double *, - int *,double *,int *,int *,double *,int *,int *,double *,double *,double *,double *,double *,double *,double *,double *, - double *,double *,double *,double *,double *,int *); -int rao_tr_disq(int *,double *,double *,double *,double *,double *,int *,double *,double *,double *,double *,double *,double *,int *,double *,int *,int *,int *, - double *,double *,double *,double *,double *,double *,double *,int *); -int rao_tr_disq_ic(int *,double *,double *,double *,double *,double *,int *,double *,double *,double *,double *,double *,double *,int *,double *,int *,int *,double *, - int *,int *,double *,double *,double *,double *,double *,double *,double *,double *,double *,double *,double *,double *,double *,int *); - -int mimetic_rect(int *,double *,double *, double *,double *,double *,double *,double *,double *, int *,double *,double *,int *,int *,double *,double *,double *,double *,double *,int *); -int mimetic_disq(int *,double *,double *,double *,double *,double *,double *,double *, int *, double *, double *, int *, int *, double *,double *,double *,double *,double *,int *); -int mimetic_tr_rect(int *,double *,double *, double *,double *,double *,double *,double *,int *, double *, double *, double *, double *, double *, double *, - double *, int *, double *, double *, int *, int *, double *,double *, double *,double *,double *,int *); -int mimetic_tr_disq(int *,double *,double *, double *,double *,double *,double *,int *, double *,double *,double *,double *,double *,double *cy, - double *, int *, double *, double *, int *, int *, double *,double *,double *,double *,double *,int *); -double echange_point_rect(int,double *,double *,double,double,double,double,double,double,double,double *,int *,double *,double *,double *); -double echange_point_disq(int,double *,double *,double,double,double,double,double,double,double *,int *,double *,double *,double *); -double echange_point_tr_rect(int,double *,double *,double,double,double,double,int *,double *,double *,double *,double *,double *,double *, - double,double,double,double *,int *,double *,double *,double *); -double echange_point_tr_disq(int,double *,double *,double,double,double,int *,double *,double *,double *,double *,double *,double *, - double,double,double,double *,int *,double *,double *,double *); - -int shen(int *,double *,double *,int *,double *,int *,int *,double *,double *,double *,double *, double *, int *); - -int shen_ic(int *,double *,double *,int *,double *,int *,double *,int *,int *,double *,double *,double *, - double *, double *, double *,double *,double *,double *,double *,double *,int *); -int intertype(int *,double *,double *,int *, double *, double *,int *,double *,double *,double *); +void s_alea_tr_disq(int ,double*,double*,double,double,double,int,double*,double*,double*,double*,double*,double*,double); +void s_alea_tr_rect(int,double*,double*,double,double,double,double,int,double*,double*,double*,double*,double*,double*,double); +int shen_ic(int*,double*,double*,int*,double*,int*,double*,int*,int*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,int*); +int shen(int*,double*,double*,int*,double*,int*,int*,double*,double*,double*,double*,double*,int*); +int shimatani_disq_ic(int*,double*,double*,double*,double*,double*,int*,double*,int*,double*,int*,int*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,int*); +int shimatani_disq(int*,double*,double*,double*,double*,double*,int*,double*,int*,int*,double*,double*,double*,int*); +int shimatani_rect_ic(int*,double*,double*,double*,double*,double*,double*,int*,double*,int*,double*,int*,int*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,int*); +int shimatani_rect(int*,double*,double*,double*,double*,double*,double*,int*,double*,int*,int*,double*,double*,double*,int*); +int shimatani_tr_disq_ic(int*,double*,double*,double*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,int*,double*,int*,int*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,int*); +int shimatani_tr_disq(int*,double*,double*,double*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,int*,int*,double*,double*,double*,int*); +int shimatani_tr_rect_ic(int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,int*,double*,int*,int*,double*,double*,double*,double*,double*,double*,double*,double*,double*,double*,int*); +int shimatani_tr_rect(int*,double*,double*,double*,double*,double*,double*,int*,double*,double*,double*,double*,double*,double*,int*,double*,int*,int*,double*,double*,double*,int*); +double trois_point(double,double,double,double,double,double,double,double,double); +double un_point(double,double,double,double,double,double,double,double,double); +double ununun_point(double,double,double,double,double,double,double,double,double); diff --git a/src/adssub.h b/src/adssub.h index baf0554..85e85ff 100755 --- a/src/adssub.h +++ b/src/adssub.h @@ -3,31 +3,26 @@ #include <limits.h> #include <R_ext/PrtUtil.h> -double Pi(); -void progress(int,int*, int); -/*double alea ();*/ -void freeintvec (int *); -void freetab (double **); -void freevec (double *); -void taballoc (double ***,int,int); -void tabintalloc (int ***,int,int); -void freeinttab (int **); -void vecalloc (double **vec, int n); -void vecintalloc (int **vec, int n); -double bacos(double a); -void decalVal(double *,int,double); -void decalRect(int,double *,double *,double *,double *,double *,double *); -void decalCirc(int,double *,double *,double *,double *,double); -void decalRectTri(int,double *,double *,double *,double *,double *,double *, - int,double *,double *,double *,double *,double *,double *); -void decalCircTri(int,double *,double *,double *,double *,double, - int,double *,double *,double *,double *,double *,double *); -void decalRect2(int,double *,double *,int,double *,double *,double *,double *,double *,double *); -void decalCirc2(int,double *,double *,int,double *,double *,double *,double *,double); -void decalRectTri2(int,double *,double *,int,double *,double *,double *,double *,double *,double *, - int,double *,double *,double *,double *,double *,double *); -void decalCircTri2(int,double *,double *,int,double *,double *,double *,double *,double,int, - double *,double *,double *,double *,double *,double *); -void decalSample(int,double *,double *,double,double); -double** taballoca(int,int *); -void complete_tab(int,double **,double **,int *,int *,int *,double *,double *); +double bacos(double); +void complete_tab(int,double**,double**,int*,int*,int*,double*,double*); +void decalCirc(int,double*,double*,double*,double*,double); +void decalCirc2(int,double*,double*,int,double*,double*,double*,double*,double); +void decalCircTri(int,double*,double*,double*,double*,double, int,double*,double*,double*,double*,double*,double*); +void decalCircTri2(int,double*,double*,int,double*,double*,double*,double*,double,int,double*,double*,double*,double*,double*,double*); +void decalRect(int,double*,double*,double*,double*,double*,double*); +void decalRect2(int,double*,double*,int,double*,double*,double*,double*,double*,double*); +void decalRectTri(int,double*,double*,double*,double*,double*,double*, int,double*,double*,double*,double*,double*,double*); +void decalRectTri2(int,double*,double*,int,double*,double*,double*,double*,double*,double*, int,double*,double*,double*,double*,double*,double*); +void decalSample(int,double*,double*,double,double); +void decalVal(double*,int,double); +void freeinttab(int**); +void freeintvec(int*); +void freetab(double**); +void freevec(double*); +double Pi(); +void progress(int,int*,int); +void taballoc(double***,int,int); +double** taballoca(int,int*); +void tabintalloc(int***,int,int); +void vecalloc(double**,int); +void vecintalloc(int**,int); -- GitLab