kiste <- function(x1, x2, y1, y2, ...)
{
	polygon( c(x1,x2,x2,x1,x1), c(y1,y1,y2,y2,y1) , ...)
}	

new.mos <- function(...)
{
	oldp <- par(mai=c(1,1,0,0))
	plot( c(0,1), c(0,1), type="n", axes=F, xlab="", ylab="")
	tab <- table(...)
	posvec <- vector("numeric",length=length(dim(tab)))
	draw.mosaic(tab, c(0,1,0,1), sum(tab), "x", posvec)
	par(oldp)
}

draw.mosaic <- function( tab, coord, n, dir, posvec)
{
	if( sum(tab) == 0)
		return()
	level <- length(dim(tab))
	if( level == 0 )
		level <- 1
	eps <- 0.001+0.001*level
	if( length(posvec) < 3 ||
	    (all(rev(posvec)[(1:(length(posvec)/2))*2]<=1)   && dir == "x" ) || 
	    (all(rev(posvec)[(1:(length(posvec)/2))*2-1]<=1) && dir == "y" ) )
		label <- T
	else
		label <- F
	if( level >= 2 )
	{
		cum  <- c(0,cumsum(apply( tab, 1, sum)))
		ncum  <- cum / n
		for (i in 1:dim(tab)[1])
		{
			posvec[level] <- i
			switch(length(dim(tab))-1,
				{newtab <- as.array(tab[i,])},
				{newtab <- as.array(tab[i,,])},
				{newtab <- as.array(tab[i,,,])},
				{newtab <- as.array(tab[i,,,,])},
				{newtab <- as.array(tab[i,,,,,])})
			#
			if( dir == "x" )
			{	draw.mosaic(newtab,
					c(coord[1]+ncum[i]*(coord[2]-coord[1])+eps,
				        coord[1]+ncum[i+1]*(coord[2]-coord[1])-eps,
				 	coord[3:4]+c(eps,-eps)),
					cum[i+1]-cum[i],
					"y", posvec)
				if( label == T )
				      text(coord[1]+(ncum[i]+ncum[i+1])/2*(coord[2]-coord[1]),
				        -0.015*level,
					 dimnames(tab)[[1]][i])}
			else
			{
				draw.mosaic(newtab,
					c(coord[1:2]+c(eps,-eps),
					coord[3]+ncum[i]*(coord[4]-coord[3])+eps,
				        coord[3]+ncum[i+1]*(coord[4]-coord[3])-eps),
					cum[i+1]-cum[i],
					"x", posvec)
				if( label == T )
				      text(-0.015*level,
        				coord[3]+(ncum[i]+ncum[i+1])/2*(coord[4]-coord[3]),
					   dimnames(tab)[[1]][i])}
		}
	}
	else
	{
		#
		cum <- c(0,cumsum( tab ))
		ncum  <- cum / rev(cum)[1]		
		for (i in 1:length(tab))
		{
			posvec[1]   <- i
			shapevec    <- append(posvec,c(1,1,1,1))
			if( dir == "x" )
			{
				kiste(coord[1]+ncum[i]*(coord[2]-coord[1])+eps,
				      coord[1]+ncum[i+1]*(coord[2]-coord[1])-eps,
			      	      coord[3]+eps, coord[4]-eps,
				      col=shapevec[1], 	    density=shapevec[2]*10, 
				      angle=35*shapevec[3], lty=shapevec[4])
				if( label == T )
				      text(coord[1]+(ncum[i]+ncum[i+1])/2*(coord[2]-coord[1]),
				        -0.015*level,
					 dimnames(tab)[[1]][i])
			}
			else
			{
				kiste(coord[1]+eps, coord[2]-eps,
				      coord[3]+ncum[i]*(coord[4]-coord[3])+eps,
				      coord[3]+ncum[i+1]*(coord[4]-coord[3])-eps,
				      col=shapevec[1], 	    density=shapevec[2]*10, 
				      angle=35*shapevec[3], lty=shapevec[4])
				if( label == T )
				      text(-0.015*level,
        				coord[3]+(ncum[i]+ncum[i+1])/2*(coord[4]-coord[3]),
					   dimnames(tab)[[1]][i])
			}
		}

	}		
}	 

new.mos( c(rep("black",108),
 	  rep("brown",286),
	  rep("red",71),
 	  rep("blond",127)),
 	c(rep("brown",68),
 	  rep("blue",20),
 	  rep("hazel",15),
	  rep("green",5),
   	  rep("brown",119),
 	  rep("blue",84),
 	  rep("hazel",54),
 	  rep("green",29),
           rep("brown",26),
 	  rep("blue",17),
 	  rep("hazel",14),
 	  rep("green",14),
           rep("brown",7),
 	  rep("blue",94),
 	  rep("hazel",10),
 	  rep("green",16)) ,
   				as.integer(runif(592)+0.5), 
   				as.integer(runif(592)+0.5), 
#  				as.integer(runif(592)+0.5),
# 			 	as.integer(runif(592)+0.5)
 			 )

tmp <- t(matrix(c( rep( c("male","n-pr","n-ex","div"), 68),
	  rep( c("male","n-pr","n-ex","marr"), 130), 
	  rep( c("male","n-pr","extra","div"), 17), 
	  rep( c("male","n-pr","extra","marr"), 4), 
	  rep( c("male","pre","n-ex","div"), 60), 
	  rep( c("male","pre","n-ex","marr"), 42), 
	  rep( c("male","pre","extra","div"), 28), 
	  rep( c("male","pre","extra","marr"), 11),
	  rep( c("female","n-pr","n-ex","div"), 214),
	  rep( c("female","n-pr","n-ex","marr"), 322), 
	  rep( c("female","n-pr","extra","div"), 36), 
	  rep( c("female","n-pr","extra","marr"), 4), 
	  rep( c("female","pre","n-ex","div"), 54), 
	  rep( c("female","pre","n-ex","marr"), 25), 
	  rep( c("female","pre","extra","div"), 17), 
	  rep( c("female","pre","extra","marr"), 4)), 4, 1036))


new.mos( tmp[,1], tmp[,2], tmp[,3], tmp[,4])