I said: > myfun <- function(x, y, ints) { > fx <- x %/% (1/ints) > fy <- y %/% (1/ints) > txy <- hist(fx + ints*fy+ 1, breaks=0:(ints*ints), plot=FALSE)$counts > dim(fxy) <- c(ints, ints) ^^^ > return(txy) > } Of course it should be: dim(txy) <- c(ints, ints) ^^^ Sorry about that, Ray