MathGroup Archive 2000

[Date Index] [Thread Index] [Author Index]

Search the Archive

Re: Re: Re: 3D graphics

  • To: mathgroup at smc.vnet.net
  • Subject: [mg22822] Re: [mg22799] Re: [mg22720] Re: 3D graphics
  • From: Hartmut Wolf <hwolf at debis.com>
  • Date: Fri, 31 Mar 2000 01:01:24 -0500 (EST)
  • Organization: debis Systemhaus
  • References: <8bhvn1$noc@smc.vnet.net> <200003260758.CAA27554@smc.vnet.net>
  • Sender: owner-wri-mathgroup at wolfram.com

Allan Hayes schrieb:
> 
> Comparison of the following postings for time and picture quality suggests
> that there is some interesting work to be done on this kind of problem.
> ......

Dear Allan, 

I agree that time and picture quality are basic categories for this
problem.
Before commenting on that, let's shortly reflect on the problem as
posed:

> > > Bernhard Adams wrote:
> > > >
> > > > I would like to plot a body in 3D whose volume is defined by a function
> > > > f[x,y,z] which returns 1 if (x,y,z) is inside the body and 0 if outside.
> > > > How to plot this?
> > > >

Jens' proposal certainly gives a quick solution

> > Jens-Peer Kuska schrieb:
> > >
> > > since your function jump from 0 to 1 you can simply use a
> > > ContourPlot3D[] with
> > > contourlevel 1/2
> > >

As function for test I had used 

f[x_, y_, z_] := If[x^2 +  2 y^2 +  3 z^2 - 1 <= 0, 1., 0.]

(following the wisdom as to use un-symmetric test cases) but with 
Contours -> 1 (does work!), 0.5, and 0.1 (0 doesn't work) I had got
(it's 
interesting to see their differences) ugly, tank-shaped objects
contrasting 
the beauty of 

ContourPlot3D[x^2 +  2 y^2 +  3 z^2 - 1, {x, -1, 1}, {y, -1, 1}, {z, -1,
1}, 
    PlotPoints -> 5] // Timing

{36.943*Second, Graphics3D[]}

and as observed by
 
> "Bojan Bistrovic" <bojanb at python.physics.odu.edu> wrote in message
> news:8bhvn1$noc at smc.vnet.net...
> >
> > One that looks good is
> >
> > In[14]:= ContourPlot3D[fun[x, y, z],
> >   {x, -1.1, 1.1}, {y, -1.1, 1.1}, {z, -1.1, 1.1}, Contours -> {0.5},
> >     PlotPoints->{7,9}]//Timing
> >
> > Out[13]= {278.31 Second, , - Graphics3D - }
> >

it takes considerable time to get something pleasing. (And, Allan, me
too, 
I don't share the MegaHertzes of Bojan's.)

Now it is not natural to get a function with two discrete binary values
from
an analytical expression, so the first attempt should be to rework it as 
to get a smooth function to ContourPlot3D. However this is not always
possible,
e.g. when the "function" is experimental data from probing a space
region 
(not unplausible at that domain "desy" -- I asked Bernhard, but he
hasn't yet
got the opportunity to reply.) Using ListContourPlot3D then again gives
those
(slightly improved) tanks.

Having the function values tabulated

In[58]:=
(tt2 = Table[f[x, y, z], {x, -1, 1, #}, {y, -1, 1, #}, {z, -1, 1, #}] &[
          0.05]); // Timing
Out[58]= {25.126 Second, Null}

a simple idea would be to do

In[59]:=
(gg2 = Flatten[MapIndexed[If[#1 == 1, Cuboid[#2], {}] &, tt2, {-1}]]);
// 
    Timing // First
Out[59]= 5.198 Second

In[60]:= Show[Graphics3D[gg2]] // Timing
Out[60]= {862.269*Second, Graphics3D[]}

Now its the rendering algorithm that takes it's time for

In[61]:= Count[tt2, 1., 3]
Out[61]= 13725

Cuboids. Most of them are inside (for this example), so we try to get
rid of them:

In[96]:=
d = Function[{mat, trans}, Transpose[
    Map[Drop[Flatten[{#1, 0.}] - Flatten[{0., #1}], -1]&, 
       Transpose[mat, trans], {2}], trans]]; 

In[97]:= tt2x1 = d[tt2, {1, 2, 3}]; // Timing // First
Out[97]= 3.695 Second
In[99]:= tt2x2 = d[tt2, {1, 3, 2}]; // Timing // First
Out[99]= 3.505 Second
In[100]:= tt2x3 = d[tt2, {3, 2, 1}]; // Timing // First
Out[100]= 3.756 Second
In[101]:=
tt2xxx = MapThread[Plus, {tt2x1, tt2x2, tt2x3}, 3]; // Timing // First
Out[101]= 2.293 Second
In[104]:=
(gg2x = Flatten[MapIndexed[If[#1 =!= 0., Cuboid[#2], {}] &, tt2xxx,
{-1}]]); //
     Timing // First
Out[104]= 5.268 Second
In[105]:= Show[Graphics3D[gg2x]] // Timing
Out[105]= {47.418*Second, Graphics3D[]}

Times for display sum up to 65.935 Second, compared to 867.467 Second
above. It's
interesting to compare the resulting picture with Bojan's "looks good"
(which for 
this function, on my machine, was 554.417 Second)

Otherwise, as Allen had stated, you need an interpolating or smoothing
function
to ContourPlot3D. His idea 
> 
> Allan Hayes:
> 
> f1[p__] := (Plus @@
>         f @@@ Transpose[{p} +
>               Transpose[{{0, 0, 0}, {d, 0, 0}, {-d, 0, 0}, {0, d, 0},
> {0, -d,
>                     0}, {0, 0, d}, {0, 0, -d}}]])/7.0
> 
On my machine, with my function:
In[13]:=
ContourPlot3D[f1[x, y, z], {x, -1, 1}, {y, -1, 1}, {z, -1, 1}, 
    Contours -> {0.5}, PlotPoints -> 5] // Timing
Out[13]= {39.848*Second, Graphics3D[]}

gives surprisingly good results! I had previously tried another method
to mollify the data:

ftt = Plus @@ Flatten[MapIndexed[
      If[#1 === 1, With[{r = {x, y, z} - #2}, 
         E^(-r . r)], {}] & , Transpose[tt, {3, 2, 1}], 
      {-1}]]; 

Dimensions[tt]
{21, 21, 21}

The plotting however takes it's time

ContourPlot3D[Evaluate[ftt], {x, -1, 23}, {y, -1, 23}, {z, -1, 23}, 
    Contours -> {3.}, PlotPoints -> 5] // Timing
{927.073*Second, Graphics3D[]}

Of course there's the same problem as with the Cuboids above (and that 
could be fixed with similar methods) and the quality is midway between 
Allan's Graphic and the analytical ellipsoid (for same PlotPoints). The 
other proposal of Allan to use FunctionInterpolation (on my machine)

In[20]:= ContourPlot3D[f2[x, y, z], {x, -1, 1}, {y, -1, 1}, {z, -1, 1}, 
    Contours -> {0.5}, PlotPoints -> 5] // Timing
Out[20]= {44.554*Second, Graphics3D[]}

however gives comparably ugly results (I don't understand why).

Finally I had tried another approach, namely to directly generate
Polygons 
which touch the surface. For that I considered the elementary cells
(cubes) 
with values from {0, 1} at the corners. There are 256 possible
elementary 
cells, however only 22 fundamental cases. I analyzed these (effectively
only 
10 or so) and generated all 256 from these through rotations and
inversion. 
(This perhaps is an exercise in group theory.) I defined a _double_
layer one 
containing the body of the '1's and complimentary one at the innermost
'0's. 
That then will also contain isolated points (of '1's). For instance

cell = {1, 0, 0, 0, 1, 1, 0, 1};
showperm[cell, cccornrot2, pstyle[0], Polygon[{p4, p3, p7}], 
  Polygon[{p2, p4, p3}], pstyle[1], Polygon[{p1, p5, p8}], 
  Polygon[{p1, p6, p8}]]

You see the generated graphical objects in the arguments, cccornrot2 are
the 24 
rotations to generate all cases of this type (makes the graphics for
visual 
checking and simultaniously the rules for the programm to be built). 

If you want to see the program then just click here 
mailto:hwolf at debis.com?subject-xrules and send, and I'll fix up my
notebook 
and return it to you.

The resulting graphs were, though coarse quite pleasing (to me) -- no
tanks.
And perhaps, this could be the base for a more efficient ContourPlot3D.

-- Hartmut


  • Prev by Date: Re: Trying to define: Fractional Derivatives & Leibniz' display form for output and templates
  • Previous by thread: Re: Re: Re: 3D graphics
  • Next by thread: RE: Re: Re: 3D graphics