[Vennerable-commits] r72 - pkg/Vennerable/R
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Fri Oct 9 00:22:21 CEST 2009
Author: js229
Date: 2009-10-09 00:22:21 +0200 (Fri, 09 Oct 2009)
New Revision: 72
Added:
pkg/Vennerable/R/Euler.R
Modified:
pkg/Vennerable/R/02TissueDrawing.R
pkg/Vennerable/R/Circles.R
Log:
Modified: pkg/Vennerable/R/02TissueDrawing.R
===================================================================
--- pkg/Vennerable/R/02TissueDrawing.R 2009-10-08 19:12:42 UTC (rev 71)
+++ pkg/Vennerable/R/02TissueDrawing.R 2009-10-08 22:22:21 UTC (rev 72)
@@ -380,12 +380,6 @@
sector.to.xy <- function(edge,dx=.05) {
r <- edge at radius;
hand <- edge at hand
-# if (edge at from==edge at to ) {
-# (if (thetafrom - thetato
-# thetato <- theta
-# if (hand>0) { thetafrom <- 2*pi; thetato <- 0 } else { thetafrom <- 0; thetato <- 2*pi}
-# thetadist <- 2 * pi
-# } else {
thetafrom <- edge at fromTheta;thetato <- edge at toTheta
# if hand > 0 we always go anti-clockwise, decreasing thetafrom
# so thetato must be less than thetafrom
@@ -397,9 +391,8 @@
if (thetato<thetafrom) { thetato <- thetato + 2* pi }
thetadist <- thetato - thetafrom
}
-# }
arclength <- (thetadist) * r
- nintervals <- arclength/dx
+ nintervals <- max(3,arclength/dx)
theta <- seq(from=thetafrom,to=thetato,length=nintervals)
xy <- .theta.to.point.xy(theta,r,edge at centre)
}
@@ -997,8 +990,12 @@
xy1 <- all.xy; xy2 <- all.xy[ c(2:nrow(all.xy),1),]
x1 <- xy1[,1];y1 <-xy1[,2];x2<-xy2[,1];y2<- xy2[,2]
area <- .polygon.area(all.xy)
+ if (area>0) {
cx <- sum((x1+x2) * ( x1 * y2 - x2 * y1))/(6* area)
cy <- sum((y1+y2) * ( x1 * y2 - x2 * y1))/(6* area)
+ } else {
+ cx <- mean(x1);cy <- mean(y1)
+ }
centroid.xy <- matrix(c(cx,cy),ncol=2)
centroid.xy
}
@@ -1126,6 +1123,10 @@
.find.triangle.within.face <- function(drawing,faceName) {
# poor mans triangulation... subtracting ear method cf wikipedia polygon triangulation
xy <- .face.toxy(drawing,faceName)
+ A <- .face.area(drawing,faceName)
+ if (A==0) { # collinear
+ return(xy[1:3,])
+ }
xy <- rbind(xy,xy[1:2,])
fix <- NA
for (ix in 2:(nrow(xy)-1)) {
@@ -1606,7 +1607,11 @@
face2IntersectionPoints <- intersect(face2Points,intersectionPoints)
if (length(face2IntersectionPoints)==0) {
- new1 <- .addNonintersectingFace(new1,drawing2,tempface2Name)
+ if (length(newEdges)==0) {
+ new1 <- .addSetWithExistingEdges(new1,drawing2,tempface2Name)
+ } else {
+ new1 <- .addNonintersectingFace(new1,drawing2,tempface2Name)
+ }
if (!inherits(new1,"TissueDrawing")) { return(NA) } # when we can't build an invisible edge joining the two faces
} else {
new1 <- .addIntersectingFace(new1,new2,tempface2Name,face2IntersectionPoints)
@@ -1621,6 +1626,24 @@
}
+
+.addSetWithExistingEdges<- function(drawing1,drawing2,tempface2Name) {
+ # the new Set is solely comprised of existing edges
+ # so there are no new faces
+ # all we need to do is update the signatures
+ # probably an easier way to do this, by tracking which edges came from which sets...
+ # but we just cycle through the faces and update through brute force
+
+
+ for (faceName in setdiff(.faceNames(drawing1),"DarkMatter")) {
+ aPoint <- .find.point.within.face(drawing1,faceName)
+ inFace <- .is.point.within.face(drawing=drawing1,faceName=tempface2Name,point.xy=aPoint,type="set")
+ inFaceSig <- if (inFace) { "1"} else {"0"}
+ drawing1<- updateSignature(drawing1,faceName,inFaceSig )
+ }
+ drawing1
+}
+
.addNonintersectingFace <- function(drawing1,drawing2,tempface2Name) {
# drawing2 contains a single face
#must be inside one of the faces or outside them all
@@ -1857,36 +1880,6 @@
faceEdgeList
}
-if(FALSE){
-.create.edge.joining.faces <- function(drawing,outerFaceName,innerFaceName) {
- outerPoints <- .points.of.face(drawing,outerFaceName)
- innerPoints <- .points.of.face(drawing,innerFaceName)
- foundEdge <- NULL
- for (op in outerPoints) {
- for (ip in innerPoints) {
- xy <- do.call(rbind,drawing at nodeList[c(op,ip)])
- testEdge <- newEdgeLines(from=op,to=ip,xy=xy,visible=FALSE)
- # does it intersect any existing edges?
- if (!.internal.edge.drawing.intersection(drawing,testEdge)) {
- foundEdge <- testEdge
- break
- }
- }
- if (!is.null(foundEdge)) { break }
- }
- if (is.null(foundEdge)) {
-# cat("Try specifying different points for internal face")
- return(list(ok=FALSE))
- # try the xy points themselves, sigh
-
- }
- edgeName <- paste(op,ip,"invisible",sep="|")
- tel <- list(testEdge); names(tel) <- edgeName
- drawing at edgeList <- c(drawing at edgeList,tel)
- return(list(edgeName=edgeName,drawing=drawing,ok=TRUE))
-
-}
-}
@@ -2210,16 +2203,28 @@
}
.merge.faces.invisibly.split <- function(diagram) {
+ doneamerge<- TRUE
+ while (doneamerge) {
+ res <- .try.merge.faces.invisibly.split(diagram)
+ doneamerge <- res$merged
+ diagram <- res$diagram
+ }
+ diagram
+}
+
+.try.merge.faces.invisibly.split <- function(diagram) {
+ # first we identify multople faces with the same signature
fsigs <- data.frame(cbind(Name=unlist(.faceNames(diagram)),Signature=unlist(.faceSignatures(diagram))),stringsAsFactors=FALSE);
rownames(fsigs)<- 1:nrow(fsigs)
nSets <- unique(nchar(setdiff(fsigs$Signature,"DarkMatter")))
fsigs$Signature[fsigs$Signature==paste(rep("0",nSets),collapse="")] <- "DarkMatter"
- w <- lapply(split(fsigs$Name,fsigs$Signature),length)
- w <- w[w>1]
- if (length(w)>0) {
- for (wn in names(w)) {
+ FacesPerSignature <- lapply(split(fsigs$Name,fsigs$Signature),length)
+ FacesPerSignature <- FacesPerSignature [FacesPerSignature >1]
+ doingamerge <- FALSE
+ if (length(FacesPerSignature )>0) {
+ for (wn in names(FacesPerSignature )) {
wnames <- fsigs$Name[fsigs$Signature==wn]
- if (length(wnames)!=2) { stop("Can't merge multiple invisibly split faces")}
+ if (length(wnames)!=2) { warning("Can't merge multiple invisibly split faces, just trying first two")}
faceEdges1 <- .face.to.faceEdges(diagram,wnames[1])
faceVisible1 <- sapply(faceEdges1,function(x)x at visible)
if (all(faceVisible1)) { break; }
@@ -2228,6 +2233,7 @@
if (all(faceVisible2)) { break; }
commonInvisibleEdges <- intersect(sub("^-","",names(faceEdges1)),sub("^-","",names(faceEdges2)))
if (length(commonInvisibleEdges)==0) break;
+ doingamerge <- TRUE
firstVisibleIx <- min(which(faceVisible1))
firstVisibleFrom <- faceEdges1[[firstVisibleIx ]]@from
diagram<- .startFaceAtPoint(diagram,wnames[1],firstVisibleFrom )
@@ -2256,6 +2262,6 @@
diagram at edgeList[lostedges] <- NULL
}
}
- diagram
+ return(list(diagram=diagram,merged=doingamerge))
}
Modified: pkg/Vennerable/R/Circles.R
===================================================================
--- pkg/Vennerable/R/Circles.R 2009-10-08 19:12:42 UTC (rev 71)
+++ pkg/Vennerable/R/Circles.R 2009-10-08 22:22:21 UTC (rev 72)
@@ -268,17 +268,21 @@
if (cp+a < b ) {
b <- (cp+a)*smidge
}
-
- # use the SSS rule to compute angles in the triangle
- CP <- acos( (a^2+b^2-cp^2)/(2*a*b))
- B <- acos( (a^2+cp^2-b^2)/(2*a*cp))
- A <- acos( (cp^2+b^2-a^2)/(2*b*cp))
- # arbitrarily bisect one and calculate xy offsets from first centre
- theta <- B/2
- o21 <- cp * c( sin(theta), - cos(theta))
- o31 <- a * c( -sin(B-theta), -cos(B-theta))
- # pick a centre for the first one
+ # pick a centre for the first one
c1 <- c(0, dList12$r1)
+ # then place the others
+ if (a==0 || b==0 || cp==0) {
+ o21 <- cp *c(1,0) ; o31 <- a * c(0,1)
+ } else {
+ # use the SSS rule to compute angles in the triangle
+ CP <- acos( (a^2+b^2-cp^2)/(2*a*b))
+ B <- acos( (a^2+cp^2-b^2)/(2*a*cp))
+ A <- acos( (cp^2+b^2-a^2)/(2*b*cp))
+ # arbitrarily bisect one and calculate xy offsets from first centre
+ theta <- B/2
+ o21 <- cp * c( sin(theta), - cos(theta))
+ o31 <- a * c( -sin(B-theta), -cos(B-theta))
+ }
c2 <- c1 + o21
c3 <- c1 + o31
C3 <- ThreeCircles(r=c(dList12$r1,dList12$r2,dList31$r1),
Added: pkg/Vennerable/R/Euler.R
===================================================================
--- pkg/Vennerable/R/Euler.R (rev 0)
+++ pkg/Vennerable/R/Euler.R 2009-10-08 22:22:21 UTC (rev 72)
@@ -0,0 +1,53 @@
+
+
+Euler.from.Signature <- function(Signature ) {
+ w <- lapply(Signature,function(x){c(0,1)})
+ names(w) <- Signature
+ Eulers <- do.call(expand.grid,w)
+ Eulers$ESignature <- apply(data.matrix(Eulers),1,paste,collapse="")
+ Eulers
+}
+
+EulerClasses <- function(numberOfSets) {
+ V <- Venn(numberOfSets=numberOfSets)
+ vs <- data.frame(Indicator(V))
+ vs$Signature <- apply(data.matrix(vs),1,paste,collapse="")
+ vs <- subset(vs,Signature!=dark.matter.signature(V))
+ E2 <- Euler.from.Signature (vs$Signature)
+ E2 <- E2[order(E2$ESignature),]
+
+ library(gtools)
+ worder <- permutations(numberOfSets,numberOfSets)
+ worder <- lapply(1:nrow(worder),function(x){worder[x,]})
+ P2 <- lapply(worder,function(x) {
+ vs.order <- vs[,x]
+ vs.order$Signature <- apply(data.matrix(vs.order),1,paste,collapse="")
+ vs.perm <- match(vs.order$Signature,vs$Signature)
+ E2.perm <- E2[,vs.perm]
+ E2.perm$ESignature <- apply(data.matrix(E2.perm),1,paste,collapse="")
+ E2.perm$ESignature
+ })
+
+
+ # now E2 has the indicator strings generated by all possible orderings
+ # with the corresponding row names
+ E3 <- do.call(rbind,P2)
+ F3 <- unique(apply(E3,2,function(x)(unique(sort(x)))))
+ iclasses <- (sapply(F3,paste,collapse=";"))
+ rclasses <- sapply(F3,function(x)x[1])
+ irclasses <- data.frame(ESignature=rclasses,iclasses=iclasses,stringsAsFactors=FALSE)
+ Eclass <- merge(E2,irclasses)
+ rownames(Eclass) <- 1:nrow(Eclass)
+ Eclass <- Eclass[order(Eclass$ESignature),]
+
+ #However some of these (eg 0000010) correspond to
+ #patterns in which every region at least one set is empty.
+ vsnames <- names(E2[,!colnames(E2)%in% c("ESignature","iclasses")])
+
+ vsmat <- do.call(rbind,strsplit(vsnames,split=""))
+ isset <- lapply(1:ncol(vsmat),function(col)vsnames[vsmat[,col]=="1"])
+ haveset <- sapply(isset,function(setsigs)apply(Eclass[,setsigs],1,sum)>0)
+ Eclass$SetsRepresented <- apply(haveset,1,sum)
+ Eclass
+
+}
More information about the Vennerable-commits
mailing list