I want to compute the total surface area of an isosurface in 3D space as
approximated with
`computeContour3d' by adding up the areas of the triangles making up the
contour surface.
problem:
the vertex matrix returned by `computeContour3d' in general seems to
provide the vertices
not in the frame of reference in which the original data are given but
apparently rather
after some linear transformation (scaling + translation (+ rotation?) -- or
I am having some fundamental misconception of what is going on.
I'm interested in the simplest case where the input data are provided as a
3D array
on an equidistant grid (i.e. leaving the x,y,z arguments at their
defaults).
e.g. (slight modification of `example(computeContour3d)'):
library(misc3d)
x <- seq(-1,1,len=11)
g <- expand.grid(x = x, y = x, z = x)
v <- array(g$x^4 + g$y^4 + g$z^4, rep(length(x),3))
con <- computeContour3d(v, max(v), 1)
drawScene(makeTriangles(con))
this is (approximately) a cube with edge length 10 (taking the grid
spacing as the unit of length).
so the expected (approximate) surface area is 600.
indeed,
apply(con, 2, range) yields
[,1] [,2] [,3]
[1,] 1 1 1
[2,] 11 11 11
which might be interpreted as providing the vertices in coordinates
where the grid spacing is used as unit of length. however
I get an area of only about 430 instead of approx. 600 which is already a
much much larger deviation
from the ideal cube surface than I would have expected given the small
amount of smoothing at the
box edges and corners (but I have to double-check whether my triangle area
computation is right, although I believe it is).
choosing instead
x <- seq(-2,2,len=50)
however, the corresponding range of `con' is
[,1] [,2] [,3]
[1,] 13.274 13.274 13.274
[2,] 37.726 37.726 37.726
which cannot be the "grid coordinates" (which should be in the range
[1,50]). adopting this interpretation nevertheless
(vertices are given in grid coordinates)
the sum of the triangle areas only amounts to about 2600 instead of the
expected approx. 49^2*6 = 14406
question 1:
am I making a stupid error (if so which one...)?
if not so:
question 2:
is there a linear transformation from the original grid coordinates (with
range from 1 to dim(v)[n], n=1:3)
involved which yields the reported vertex coordinates?
question 3:
could someone please explain where to find this information (even if
hidden in the source code of the package)
how to convert the vertex coordinates as delivered by `computeContour3D'
to 'grid coordinates' (or true world coordinates
in general (if the x,y,z arguments are specified, too)?
for the wishlist: it would of course be nice if `computeContour3d' would
indeed return the total surface area itself,
e.g. as an attribute of the triangles matrix.
for the devs: there is a typo in the manpage of this function:
Value:
A matrix of three columns representing the triangles making up the
contour surface. Each row represents a vertex and goups of three
rows represent a triangle.
(should be `groups' instead of `goups')
--
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
