"Jaime Vives Piqueres" wrote:
>  If it is not a surface, how can it have a normal?
> I think you don't mean this... but I can't understand it.

Not exactly a normal, but the gradient has a direction. I have no idea how
to find it though...

Rune
--
3D images and anims, include files, tutorials and more:
Rune's World:    http://rsj.mobilixnet.dk (updated Nov 5)
POV-Ray Users:   http://rsj.mobilixnet.dk/povrayusers/
POV-Ray Webring: http://webring.povray.co.uk

Post a reply to this message

"Barron Gillon" wrote:
> Has anyone else done this?  Rune?  Thanks

Sorry, advanced function math is not my domain...

> p.s.  For extra credit, is there a way to determine
> how far a point is from the nearest point on the
> surface of an isosurface?

I think this can only be determined by taking lots of samples, which is
*very* slow...

Rune
--
3D images and anims, include files, tutorials and more:
Rune's World:    http://rsj.mobilixnet.dk (updated Nov 5)
POV-Ray Users:   http://rsj.mobilixnet.dk/povrayusers/
POV-Ray Webring: http://webring.povray.co.uk

Post a reply to this message

To calculate the normal vector of the function you need to know the
derivative of the function with respect to each parameter. The derivative
of a function is too complicated to calculate with pov-script, so you'll
have to provide these derivatives yourself.
  That is, if your function is f(x,y,z), you need to calculate the
derivatives df/dx, df/dy and df/dz.
  After this the normal vector is pretty simple: <df/dx, df/dy, df/dz>

  Here is an example scene which demonstrates this:


camera { location -z*10 look_at 0 angle 35 }
light_source { <100,200,-300>, 1 }

#declare F = function { .5*x^2+y^2+z^2-1 }
#declare dFdx = function { x }
#declare dFdy = function { 2*y }
#declare dFdz = function { 2*z }

#macro nF(P)
  vnormalize(<dFdx(P.x,P.y,P.z), dFdy(P.x,P.y,P.z), dFdz(P.x,P.y,P.z)>)
#end

isosurface
{ function { F(x,y,z) }
  contained_by { box { <-2,-1,-1>,<2,1,1> } }
  max_gradient 2.352
  pigment { rgb x } finish { specular .5 }
}

union
{ #declare Pnt = <-.4,.4,-1>;
  cylinder { Pnt, Pnt+nF(Pnt), .1 }
  #declare Pnt = <.8,.6,-1>;
  cylinder { Pnt, Pnt+nF(Pnt), .1 }
  #declare Pnt = <.1,-.6,-1>;
  cylinder { Pnt, Pnt+nF(Pnt), .1 }
  pigment { rgb y } finish { specular .5 }
}



-- 
#macro N(D)#if(D>99)cylinder{M()#local D=div(D,104);M().5,2pigment{rgb M()}}
N(D)#end#end#macro M()<mod(D,13)-6mod(div(D,13)8)-3,10>#end blob{
N(11117333955)N(4254934330)N(3900569407)N(7382340)N(3358)N(970)}//  - Warp -

Post a reply to this message

Wasn't it Barron Gillon who wrote:
>not trace!  I need to find the normal at an arbitrary point that
>_may_or_may_not be on the surface.  Trace finds an intersection with a
>surface and then returns the normal at the intersection, but I need the
>normal from an arbitrary point given an arbitrary function.  Trace is good,
>but not quite what I'm looking for.  Thanks anyway.

If the point isn't on the surface, then there's obviously no normal at
that point.

Do you mean that you want to know if its possible to determine the
direction of the Grad vector of an arbitrary scalar field?

-- 
Mike Williams
Gentleman of Leisure

Post a reply to this message

Wasn't it Barron Gillon who wrote:
>not trace!  I need to find the normal at an arbitrary point that
>_may_or_may_not be on the surface.  Trace finds an intersection with a
>surface and then returns the normal at the intersection, but I need the
>normal from an arbitrary point given an arbitrary function.  Trace is good,
>but not quite what I'm looking for.  Thanks anyway.

Actually you can still use "trace()". What you have to do is declare a
temporary isosurface which is the same as the original one, but has a
different value of "threshold". The value of "threshold" to use is the
value of the function at the point you are interested in.

E.g. suppose you have an arbitrary function F(x,y,z) and you want to
find the "normal" at the point P, then you can temporarily declare

#declare Temp =
isosurface {function {F(x,y,z)}
        threshold F(P.x,P.y,P.z)
        ...
        }

and perform a "trace()" on that.

[We know that the function evaluates to F(P.x,P.y,P.z) at the point P,
and the Temp isosurface is the set of all points for which the function
evaluates to that value, so the Temp isosurface contains the point P]

Since we don't ever render the Temp isosurfaces, they shouldn't pose too
much of a processing overhead.

-- 
Mike Williams
Gentleman of Leisure

Post a reply to this message

Mike Williams <mik### [at] nospamplease> wrote:
: If the point isn't on the surface, then there's obviously no normal at
: that point.

  I think that the idea is that if the "threshold value" was set so that the
surface for that function would pass through that point, then what would be
the normal vector. It's the normal vector to the equipotential surface of
that function which passes through that point (or something similar). It is
a well-defined concept.
  (And I answered the original question in another article.)

-- 
#macro M(A,N,D,L)plane{-z,-9pigment{mandel L*9translate N color_map{[0rgb x]
[1rgb 9]}scale<D,D*3D>*1e3}rotate y*A*8}#end M(-3<1.206434.28623>70,7)M(
-1<.7438.1795>1,20)M(1<.77595.13699>30,20)M(3<.75923.07145>80,99)// - Warp -

Post a reply to this message

Barron Gillon wrote:
> 
> Does anyone know (preferably have a macro they could send me) how to find
> the surface normal of an isosurface at an arbitrary <x,y,z>?  If I remember
> correctly, a macro that returns <g'(x),h'(y),i'(z)> where f(x,y,z) =
> <g(x),h(y),i(z)> and x, y, and z are specified will suffice.  The catch is
> of course that f() would be an arbitrary function, and I don't know how to
> find numeric derivatives.  Did I miss something in the standard includes
> that would do this?  Has anyone else done this?  Rune?  Thanks
>...

Warp has shown in another post how to find the gradient
of a function f(x, y, z) for any point in 3D space by
analytical calculation of the 3 partial derivatives of
the function.

But if you want to have POV to estimate this gradient
by applying some numerical methods, (in order to find
the partial derivatives), then you may have a look at
the code I provided in my answer to Jan Walzer's thread
12. Nov. 2001:

"Help on isosurface - distances - other stuff..."

news://news.povray.org/3BF0541A.62BC8C1F%40hotmail.com
http://news.povray.org/povray.advanced-users/19992/


Tor Olav

Post a reply to this message

Barron Gillon wrote:
> 
> Does anyone know (preferably have a macro they could send me) how to find
> the surface normal of an isosurface at an arbitrary <x,y,z>?  If I remember
> correctly, a macro that returns <g'(x),h'(y),i'(z)> where f(x,y,z) =
> <g(x),h(y),i(z)> and x, y, and z are specified will suffice.  The catch is
> of course that f() would be an arbitrary function, and I don't know how to
> find numeric derivatives.
>...

At the end of my code below you'll find such a macro.

There's also an analytical solution provided, so that
you can compare the results from both the methods.
(I have chosen a torus function for this example.)

POV-Ray v3.5 are needed for the code below.
(Because it allows passing of functions as parameters
to macros.)

Note that there are other and more accurate (and
slower) numerical methods to find derivatives for
such functions, but the one used below will suffice
for many applications.


Tor Olav


// ===== 1 ======= 2 ======= 3 ======= 4 ======= 5 ======= 6 ======= 7
// Copyright 2001 by Tor Olav Kristensen
// Email: tor### [at] hotmailcom
// http://www.crosswinds.net/~tok
// ===== 1 ======= 2 ======= 3 ======= 4 ======= 5 ======= 6 ======= 7

#version 3.5;

// ===== 1 ======= 2 ======= 3 ======= 4 ======= 5 ======= 6 ======= 7

#macro PrintVector(v0)

 #debug concat("<", vstr(3, v0, ", ", 0, -1), ">")

#end // macro PrintVector

// ===== 1 ======= 2 ======= 3 ======= 4 ======= 5 ======= 6 ======= 7
// Define a torus function

#declare TorusFn =
function(x, y, z, Rmaj, Rmin) {
  sqrt((sqrt(x^2 + z^2) - Rmaj)^2 + y^2) - Rmin
}

// ===== 1 ======= 2 ======= 3 ======= 4 ======= 5 ======= 6 ======= 7
// Define constants for the torus function that we are investigating.
// Also define a point in space that we are going to find the gradient
// vector at.

#declare Rmajor = 10;
#declare Rminor = 2;
#declare pA = <-3, 2, -5>*3;

// ===== 1 ======= 2 ======= 3 ======= 4 ======= 5 ======= 6 ======= 7
// Find the analytical solution
// Note that the partial derivatives are not dependant of R minor.)

#declare DerxTorusFn =
function(x, y, z, Rmaj) {
  x/sqrt((sqrt(x^2 + z^2) - Rmaj)^2 + y^2)
  *(sqrt(x^2 + z^2) - Rmaj)/sqrt(x^2 + z^2)
}

#declare DeryTorusFn =
function(x, y, z, Rmaj) {
  y/sqrt((sqrt(x^2 + z^2) - Rmaj)^2 + y^2)
}

#declare DerzTorusFn =
function(x, y, z, Rmaj) {
  z/sqrt((sqrt(x^2 + z^2) - Rmaj)^2 + y^2)
  *(sqrt(x^2 + z^2) - Rmaj)/sqrt(x^2 + z^2)
}

#declare vGradAnalytical =
<
  DerxTorusFn(pA.x, pA.y, pA.z, Rmajor),
  DeryTorusFn(pA.x, pA.y, pA.z, Rmajor),
  DerzTorusFn(pA.x, pA.y, pA.z, Rmajor)
>;

#debug "\n"
PrintVector(vGradAnalytical)
#debug "\n"

// ===== 1 ======= 2 ======= 3 ======= 4 ======= 5 ======= 6 ======= 7
// Find a numerical solution

#macro vGradient(Fn, p0)
  
  #local H = 0.00001;

  <Fn(p0.x + H, p0.y, p0.z) - Fn(p0.x - H, p0.y, p0.z),
   Fn(p0.x, p0.y + H, p0.z) - Fn(p0.x, p0.y - H, p0.z),
   Fn(p0.x, p0.y, p0.z + H) - Fn(p0.x, p0.y, p0.z - H)>/2/H;

#end // macro vGradient

#declare vGradNumerical =
vGradient(function { TorusFn(x, y, z, Rmajor, Rminor) }, pA);

#debug "\n"
PrintVector(vGradNumerical)
#debug "\n\n"

// ===== 1 ======= 2 ======= 3 ======= 4 ======= 5 ======= 6 ======= 7

Post a reply to this message

Wasn't it Barron Gillon who wrote:
>not trace!  I need to find the normal at an arbitrary point that
>_may_or_may_not be on the surface.  Trace finds an intersection with a
>surface and then returns the normal at the intersection, but I need the
>normal from an arbitrary point given an arbitrary function.  Trace is good,
>but not quite what I'm looking for.  Thanks anyway.

Here's an example


global_settings {assumed_gamma 1.0}

camera { location  <-1, 1, -4.5> look_at <0, 0, 0>}

sky_sphere { pigment { rgb 1} }

light_source {<-100,200,-500> colour rgb 1}

// ----------------------------------------

// An arbitrary function
#declare F=function { sin(1.5*x)*y }

// An isosurface using this function
isosurface {function {F(x,y,z)}
        threshold (F(0.2,-0.2,0))
        max_gradient 5
        contained_by{box{-2,2}} open
        pigment {rgb 0.9}
        finish {phong 0.5 phong_size 10}
}

// Scan for "normals"
#declare Norm=<0,0,0>;
#declare X=-2;
#while (X<2)
 #declare Y=-2;
 #while (Y<2)
  #declare P=<X,Y,-0>;//Current point of interest
  // Declare a temporary isosurface 
  #declare Temp = 
    isosurface {function {F(x,y,z)}
      threshold F(X,Y,0)
      max_gradient 5000
      contained_by{box{-2.1,2.1}} open
      }
  // Trace it
  #declare Q=trace(Temp, P-y*0.1, y, Norm);
  // Just in case P-y*0.1 happens to be in the surface
  #if (vlength(Norm)=0)
    #declare Q=trace(Temp, P-x*0.1, x, Norm);
  #end
  #if (vlength(Norm)=0)
    #declare Q=trace(Temp, P-z*0.1, z, Norm);
  #end

  // A cylinder pointing in the direction of the "normal"        
  cylinder {P,P+Norm*0.2,0.02 pigment {rgb x} no_shadow}
 #declare Y=Y+0.2;
 #end
#declare X=X+0.2;
#end

Post a reply to this message

Tor Olav Kristensen wrote:
>...
> 
> #macro vGradient(Fn, p0)
> 
>   #local H = 0.00001;
> 
>   <Fn(p0.x + H, p0.y, p0.z) - Fn(p0.x - H, p0.y, p0.z),
>    Fn(p0.x, p0.y + H, p0.z) - Fn(p0.x, p0.y - H, p0.z),
>    Fn(p0.x, p0.y, p0.z + H) - Fn(p0.x, p0.y, p0.z - H)>/2/H;
> 
> #end // macro vGradient
>...

I just noticed that this macro had a little error.
Here's a better version of it:


#macro vGradient(Fn, p0)

  #local H = 0.00001;

  (<Fn(p0.x + H, p0.y, p0.z) - Fn(p0.x - H, p0.y, p0.z),
    Fn(p0.x, p0.y + H, p0.z) - Fn(p0.x, p0.y - H, p0.z),
    Fn(p0.x, p0.y, p0.z + H) - Fn(p0.x, p0.y, p0.z - H)>/2/H)

#end // macro vGradient


Tor Olav

Post a reply to this message