Search Unity

  1. Unity support for visionOS is now available. Learn more in our blog post.
    Dismiss Notice

Infinite parallax hair volume using wrapping grid tracing (attempt)

Discussion in 'Shaders' started by neoshaman, Dec 5, 2018.

  1. neoshaman


    Feb 11, 2011
    There is this thing I need to get done, it's an extension and end goal of experiment with relief function shader. Basically I want to trace in fragment shader hair that had parallax property. Since I'm having problem with simple math due to a past burn out, my brain blanking out on thing I should know, I figure out that chronicling the attempt publicly and verbalizing it to other people might help me get thing fall into place, and maybe some people might be interested and drop some clue.

    Hair rendering is notoriously hard and costly, due to the high number of occluded elements, which may or not contribute to the final image. I propose a technique that use a simplified model of hair, as parallel line in an extruded grid, and a mathematical representation that would allow to find the first visible intersection to a line without having to pay the cost for the occluded many, and potentially infer how many hair has been traversed with the same operation at the exit of the volume.

    The main idea is to realize that the simplified hair is "coherent", ie all are in the same position in the grid (idea to shuffle it will be explored later) and therefore all cells are indistinguishable from another, they can be superposed. This mean that a ray entering at one cell boundary will have a predictable trajectory wrapping in the superposed space. The challenge will be to figure out the math of intersection in that wrapping space, which would be the equivalent of determining occlusion in one hit. To simplify further, we will move the hair position to the origin corner of the cell so it's 0,0. We will adapt the DDA algorithm to get to that result.
    Further analysis show we can simplified the problem more by simply projecting to 2d position as it's extruded, therefore we can align the grid with one axis, then project to 1d, where it become a problem similar to the " how many timed clock hand overlap in a day" or planetary alignment puzzle

    The end goal is to have a version of the shader that can also represent curly hair (intersection with helix line).

    Given I found all the mathematical foundation it should have been trivial from there ... I'll add image later to illustrate the technique.
    AlanMattano likes this.
  2. neoshaman


    Feb 11, 2011
    Oups, end of years + a loss in the family + a new occupation didn't make me live up to those statement. Slowly coming back to it...

    a. the basic
    I was trying to model afro textured hair cut, but wasn't satisfied with the result. Afro textured is inherently volumetric and I wanted to reflect that fact.
    There is multiple way to implement that idea but my first instinct was to try a parallax shader, to model the superposition of strand, basically representing multiple layer displaced virtually by the view vector.

    Going through the math I realized there was opportunity to expand on that. But let's review the math of layer parallax (I hope there is no mistake, I'm not the most competent mathematician here).

    Let's consider an UV surface, with layer at regular interval d, and a view normal vector v. Let's simplify the representation of a 2d cross section.

    hair parallax.png
    Due to thales theorem, we can see than the offset is a ratio k of distance to the plane. If we set d=1 then k<1 because v is a unit vector as such k is equivalent to the depth when v is looking from the top.

    Once we see the result we understand that the sample are regularly space in a predictable way, we might be able to derive mathematically an intersection with a fixed position, and skip sampling multiple time an occlusion assuming it happen at the same regular interval.

    If we assume hair strands are regularly space line in uv space, we can approximate the hair volume occlusion mathematically. While strands have thickness, we will assume they are mathematical line and try to find the basis of the math from there.

    The insight is that multiple plane at the same distance is equivalent to a wrapping plane, also UV wrap around the 0-1 range. SO if we can port the math to wrapping space ... it would be equivalent to tracing through an infinite regular grid. Also since we have wrapping, we can realize it's equivalent to looping, stepping at regular interval is therefore equivalent to the clock needle overlap problem.

    The other trick is to realize that the support (the triangle face) can bias the visual such has the straight line become curved while the math stay the same (operating on UV space not in visual space). Of course it has limitation (such has curving beyond the horizon wouldn't work, as the ray would need to intersect twice the same strand.

    It's not really a new idea, some people have experimented with modulo space, though I need to find one that give interesting result for hair strands occlusion.

    I'll explore that next time.
  3. neoshaman


    Feb 11, 2011
    Wrapping the space

    First let’s assume the hair strands are all situated at position 1 in x in all layers, and the view rays start at position 0 in x on the first. We will trace rays through all layers such as they intercept the strand position at each layer.

    If we trace back to the first layer, every intersections with a layer before hitting a strand position, we can observe a pattern. Because the ray slope is constant, and each layer distance is the same, the projected intersection to the first layer are equidistant. We can see there is a correlation between the depth of the layer, where the ray hit the strand is, and the number of projected intersection. Therefore the intersections number give us the depth. We then know we can just divide the distance, between the ray x origin and the strand x position, by the number of intersection to get the depth.

    The corollary is that we can guess which depth a strand is hit, by looking at the number of step it takes, which is a function of the slope. Figuring the slope at the first layer give you the delta of the step, you then just to divide the distance.
    test hair drawing.png
    Wrapping horizontally is similar, the difference is that the number of step happen vertically. But really we just need to observe, that the distance to the targeted strand horizontally, is a multiple of the distance to the reference strand on the first later.

    Those simplified observation are here to demonstrate how we can infer positions of virtual strands (ie strand beyond the first layer and horizontal wrapping), based on the rays slopes vs the reference strand position, in a kind of wrapping trigonometry. And also how we can basicaly simplify the math to a linear convergence using division or modulo.

    However to generalize further, we still need to demonstrate how it works when we start changing the slopes and the offset of the ray, instead of assuming all ray point to a strand from the zero position. And we need to verify when the math allow convergence to a virtual position of a strands and when it does not. Doing so would allow us to guess the first hit in a virtual volume.

    Once we have done that, we have effectively found a way to trace an infinite volume of primitive, however hair has thickness, we will need to find a way to intercept more than an infinitely small position. Ideally we would also generalized to different kind of wrapping shape, such has helix, which would allow us to emulate curly hair.
  4. neoshaman


    Feb 11, 2011
    Up until we cheated by locking the ray slopes and position to the a single place into the receiving surface and only to slopes guaranteed to hit a target point. It was done to demonstrate how it would translate to a 1d problem of interval overlapping. But now we need to introduce two new parameter:
    1. arbitrary slopes
    2. arbitrary starting position of the ray.

    1. variation of slopes
    As we have shown earlier, the slope of the ray translate into discrete steps, with same size, on the "surface", and there is collision when they are basically a multiple of the wrapping size, and that multiple give you the depth or distance of the collision within the infinite space. But what happen when we have an arbitrary slope, that give us an arbitrary step, and we want to find the multiple to the wrapping size, does that multiple even exist?

    Basically we can reformulate the problem as beat synchronicity, let say you have two beat, one that goes on twos and the other on three, when will the beat pulse together? In this simple set up it's easy, they will be synchrone every 6 beats, that's the "collision depth". But what does that mean in mathematical form, how can we generalize it? seems simple though isn't it?

    2. variation of position
    Up until we cheated (again), by fixing the position at the start we made it easier, since the space is wrapping the starting and the target are basically the same position, and since we normalized the wrapping distance, it also mean we had neat number to play with. But what happen when we move the ray to anywhere within the wrap space, as it should happen in a real case? Well let's see the beat analogy again, let's say we have two beat on twos, but one start half way after the first beat, they will never be in sync again, ie similar beat with an offset, they will never pulse, can beat we translate that into mathematical form, to check that two beat will never converge to the same time?

    I use beat analogy, but we can basically use other one, since the space is wrapping, we can also see the surface "segment" as a circle, since we have stuff moving at different rate on that circle is analogous to a clock needles, which mean we can relate to the clock needle superposition problem (itself analogous to the planet alignement problem), with the added complication that we are moving by "discrete steps" rather than a "continuous rate".

    So basically it's this, given two elements moving at discrete step, can they overlap and when? The minimal distance will gave us the depth of the collision, and probably we can define the periodicity of the collision after that.

    So here is where my competence end, I need to translate those insight into mathematical form to resolve.
  5. neoshaman


    Feb 11, 2011
    Up until now I showed that the problem could express as a 1d line with interval overlapping. But since the space is wrapping on itself, that is the first line and teh last line is basically the same, instead of an infinite line with repeated "beat", we can represent it as a circle with interval stop around it:
    In this example the black tick overshoot the target after one wrap, the red tick is the resulting position of overshooting. We need to calculate how many wrap before it can actually overlap the target.

    Using this representation allow us to find equivalence, it looks like a clock, it also look like circular planet orbital. Turns out we can probably get inspire by that to find mathematical equivalence. If you think a bit about it, it turns out that's a similar problem than planets alignment, or the clock needle superposition riddle. The hair target

    Let's use planetary alignment.

    We find this googling on the net:

    which gave us this:
    If it's correct we have this:
    let say we have two point moving at two rate,
    - the position of each point = rate * t.
    - Therefore position1 = position2 it mean that rate1 * t = rate2 * t,
    - that is (rate1-rate2) * t .
    What's missing is the bigger period in which rate1 and rate2 are inscribed,
    - we can get it by having rate1 * rate2.
    - So rate1 and rate2 are sync for t= k(rate1*rate2) / (rate1 -rate2) with k being the number of period.

    However we have a specific case where rate2 = the normalized wraping size, rate2 = 1, for k = 1 we have t = rate1 / (rate1-1).

    But we only assume no rate as an offset,
    - so if we check with offset we have position1+ offset = position2,
    - which translate (rate1 *t) +offset = rate2 *t
    - which lead to t * (rate1 - rate2) + offset
    and we can probably write the final result as t = k(rate1*rate2) / (rate1-rate2+offset)
    which finalized as
    t = k(rate1) / (rate1+offset-1) which has no solution when rate1+offset = 1

    Probably? Any mistake?

    Extruding 2d intersections points
    More math (help given by someone else in another forum), this would help us define the line in 3d, since hair is axis align with the support geometry UV space (the deformation/orientation of that support geometry being the approximation of the hair direction) once we have the virtual depth, we can compute the virtual position as a function of the origin and the slope.

    Jittering the depth alignment
    Anyway, we can have some more insight, we have assumed the wrapping itself was fixed, if the wrapping windows itself move, it's equivalent to find 3 planets alignment ... What's the use? well if we stay in the two beat system, we get all hair perfectly aligned in regular interval, if the wrapping move too, we can jitter the depth predictably by multiplying by clever number.

    However this isn't that easy, if we go back to the 2d representation, let say we have an angle where the beat is bigger than the wrapping size, we would wrap horizontally first, so the jittering position will not match a ray sampling the nearby cells position, because the order of traversal would matter, so we need to define the wrapping displacement as a 2d problem. Which I haven't figure out yet.

    More works is needed
    Up until then we have assume mathematical object, ie point that have no size. In order to properly compute hair, we need thickness, we need to find a way to detect when our beat is within a distance of the target to get a hit, not just right on the target. Once we can do that, we can probably extend to a 3d volume cross section (the circle shape of the hair), the problem still being finding a single point, so the extrusion should remain the same.

    We also need to find a way to intersect helicoidal shape, they are too repeating, so a wrapping function is fine, the challenge is that there isn't cheap formula to intersect sin wave, which is needed to get the helix.

    The end game of this exercice, infinite field of infinite helix shape!

    If we were to just do infinite straight hair, probably a beat base solution is overkill, on shadertoy, there is already wrapping tracing shader, we could just extend them to get hair approximation.
    Last edited: Jul 12, 2019
    Invertex likes this.
  6. neoshaman


    Feb 11, 2011
  7. neoshaman


    Feb 11, 2011
    Doing recursive interception on wrapping helix is hard, the equation of the helix is x = r cos t, y =r sin t, z=a t ... all those sin and cos create problem, there is no simple equation to intercept them. And it's a curve that is not convex, it's hard to find simplification.

    I try simplifying by boxing the shape abusing the inherent symmetry, and then approximate it with a x² curves that we can reliably intercept but I wasn't able to see anything that will put me closer to the goal.

    Up until I decided to step back and instead of trying to intercept the coil, intercept the support cylinder, it can be reduce to a 2D problem of line circle intersection, which doesn't have any sin or cos, but then that lead to a realization that the cylinder unroll into a 2D wrapping space where the coil is basically just the diagonal line ...

    Now the idea is to try and translate the cylinder intersection equation into an equation of all potential wrapping hit points, given that all points are in the form of (x+n, y+m, z+o), with n,m,o being the wrapping period. Convert that into cylinder space, hope that no sincos appear and trying to intercept the diagonal for hit on the helix...
    Last edited: Oct 1, 2019
  8. neoshaman


    Feb 11, 2011
    I spend some time with math equation ... that lead nowhere :oops: Trying to express thing into varying variable, couldn't do it, maybe I wasn't using the proper visualization:


    Oh ...
    You can turn the whole 2D cylinder projection into a single dimension where the line move at integer wrapping step ... where I have seen that before?

    That also open the possibility of doing 3D sampling of pseudo noise integral, ie doing cloud without raymarching a noise?

    Anyway, now I can cleanly project the point back to the cylinder with a cos and have a definitive way to measure hit and miss of rays. Also it solve the vertical on the unwrapping of the cylinder, but now we need to find the delta ie where on that vertical the infinite set of point fall and how to find that intersection with the helix ... We have the horizontal delta but we need the vertical delta, if possible as a predictable pattern. My guess is that it will be nasty SIN, my hope is that it is quadric instead ...

    Anybody has any idea?
    Last edited: Oct 1, 2019
  9. neoshaman


    Feb 11, 2011
  10. neoshaman


    Feb 11, 2011
    Code (CSharp):
    1. Hairtech flat only
    2. [code]
    3. #define DARKEN_OVER_DISTANCE 1  // This makes it easier to see the different layers in a static image, not to hide a max distance.
    4. #define SHOW_2D_SHAPE        1  // This makes the 2d shape be shown in upper left corner
    5. const float c_camearaDistance    = 6.0;
    6. const float c_cameraViewWidth    = 24.0;
    7. //============================================================
    8. //============================================================
    9. bool BooleanFunction_Square (vec2 current)
    10. {
    11.     return (current.x < 0.25 || current.x > 0.75 || current.x < 0.25 || current.x > 0.75);
    12. }
    13. //============================================================
    14. float NumberStepsFunction_Square (vec2 current, vec2 stepValue)
    15. {
    16.     // if it's already in the shape, no steps to take
    17.     if (BooleanFunction_Square(current))
    18.         return 0.0;
    20.     float stepsX;
    21.     if (stepValue.x < 0.0)
    22.         stepsX = (current.x - 0.25) / -stepValue.x;
    23.     else
    24.         stepsX = (0.75 - current.x) / stepValue.x;
    26.     return ceil(stepsX);
    27. }
    28. //============================================================
    29. void mainImage( out vec4 fragColor, in vec2 fragCoord )
    30. {
    31.     // set up the camera
    32.     vec3 cameraPos;
    33.     vec3 rayDir;
    34.     {
    35.         vec2 percent = (fragCoord / iResolution.xy) - vec2(0.5,0.5);
    36.         vec3 offset = vec3(0.5, 0.5, iTime+0.01);
    37.         float angleX = 0.0;
    38.         float angleY = 0.0;
    39.         if (iMouse.z > 0.0) {
    40.             vec2 mouse = iMouse.xy / iResolution.xy;
    41.             angleX = 3.14 + 6.28 * mouse.x;
    42.             angleY = (mouse.y - 0.5) * 3.14;//(mouse.y * 3.90) - 0.4;
    43.         }
    44.         vec3 cameraFwd    = (vec3(sin(angleX)*cos(angleY), sin(angleY), cos(angleX)*cos(angleY)));      
    45.         vec3 cameraRight = normalize(cross(vec3(0.0,1.0,0.0),cameraFwd));
    46.         vec3 cameraUp = normalize(cross(cameraFwd, cameraRight));
    47.         cameraPos = vec3(0.0, 0.0, -1.0) + offset;
    48.         vec3 cameraTarget = vec3(0.0, 0.0, 0.0) + offset;
    49.         float cameraViewHeight    = c_cameraViewWidth * iResolution.y / iResolution.x;
    50.         vec3 rayTarget = cameraPos +  cameraFwd * c_cameraDistance + cameraRight * c_cameraViewWidth * percent.x + cameraUp * cameraViewHeight * percent.y;
    51.         rayDir = normalize(rayTarget - cameraPos);
    52.     }
    53.     // keep the camera in a unit cube
    54.     cameraPos = fract(cameraPos);
    56.     // If ray facing negative on z axis, just flip direction and invert where we are in the cube on the z axis.
    57.     // Now we only have to deal with positive z directions.
    58.     if (rayDir.z < 0.0) {
    59.         rayDir *= -1.0;
    60.         cameraPos.z = 1.0 - cameraPos.z;
    61.     }
    63.     // calculate the 3d position of the next two plane intersections
    64.     float intersection1Distance = (1.0 - cameraPos.z) / rayDir.z;
    65.     float intersection2Distance = (2.0 - cameraPos.z) / rayDir.z;
    66.     vec3 intersection1 = fract(cameraPos + rayDir * intersection1Distance);
    67.     vec3 intersection2 = fract(cameraPos + rayDir * intersection2Distance);
    69.     // Calculate how much the uv changes from intersection1 to intersection2.
    70.     // Convert it from [0,1] to [-0.5, 0.5].
    71.     // We need to know this to know if the uvs are going positive or negative and by how much, on each axis.
    72.     vec2 uvStep = intersection2.xy - intersection1.xy;
    73.     uvStep = fract(uvStep + 0.5) - 0.5;
    75.     // calculate how many steps it takes to hit something on the X and Y axis and take whichever hits first.
    76.     float steps = 0.0;
    77.     steps = NumberStepsFunction_Square(intersection1.xy, uvStep);
    78.     // calculate how far it is to the intersection we found
    79.     float dist = (1.0 - cameraPos.z) / rayDir.z + steps / rayDir.z;
    81.     #if DARKEN_OVER_DISTANCE
    82.     float tint = clamp(1.0 - dist / 5.0, 0.0, 1.0);
    83.     #else
    84.     float tint = 1.0;
    85.     #endif
    87.     // calculate the hit point
    88.     vec3 hitPoint = cameraPos + rayDir * dist;
    89.     vec2 uv = hitPoint.xy;
    91.     // sample the texture
    92.     fragColor = vec4(texture(iChannel0, uv).rgb * tint, 1.0);
    93. }


    Code (CSharp):
    1. //============================================================
    2. float binarySign (float v)
    3. {
    4.     return step(0.0, v) * 2.0 - 1.0;
    5. }
    7. //============================================================
    8. // returns t
    9. // circle xy = position, z = radius
    10. // Adapted from "real time collision detection" IntersectRaySphere()
    11. float RayIntersectCircle (in vec2 rayPos, in vec2 rayDir, in vec3 circle)
    12. {
    13.     // rayDir isn't normalized, so normalize it but remember it's length
    14.     float rayLen = length(rayDir);
    15.     rayDir = normalize(rayDir);
    17.     vec2 m = rayPos - circle.xy;
    18.     float b = dot(m, rayDir);
    19.     float c = dot(m, m) - circle.z*circle.z;
    21.     // Exit if the ray is outside the circle and pointing away from the circle
    22.     if (c > 0.0 && b > 0.0)
    23.         return -1.0;
    25.     float discr = b*b - c;
    27.     // A negative discriminant means it missed the sphere
    28.     if (discr < 0.0)
    29.         return -1.0;
    31.     float t = -b - sqrt(discr);
    32.     if (t < 0.0)
    33.         t = -b + sqrt(discr);
    35.     return t / rayLen;
    36. }
    Code (CSharp):
    2. //============================================================
    3. float NumberStepsFunction_L_OneAxis (float current, float stepValue)
    4. {
    5.     float steps;
    6.     if (stepValue < 0.0)
    7.         steps = (current - 0.5) / -stepValue;
    8.     else
    9.         steps = (1.0 - current) / stepValue;
    10.     return steps;
    11. }
    13. //============================================================
    14. bool BooleanFunction_L (vec2 current)
    15. {
    16.     return (current.x < 0.5 || current.y < 0.5);
    17. }
    19. //============================================================
    20. float NumberStepsFunction_L (vec2 current, vec2 stepValue)
    21. {
    22.     if (BooleanFunction_L(current))
    23.         return 0.0;
    25.     float stepsX = NumberStepsFunction_L_OneAxis(current.x, stepValue.x);
    26.     float stepsY = NumberStepsFunction_L_OneAxis(current.y, stepValue.y);
    27.     return ceil(min(stepsX, stepsY));
    28. }
    Code (CSharp):
    2. //============================================================
    3. bool BooleanFunction_Checker (vec2 current)
    4. {
    5.     if (current.x >= 0.5)
    6.     {
    7.         return current.y < 0.5;
    8.     }
    9.     else
    10.     {
    11.         return current.y >= 0.5;
    12.     }
    13. }
    15. //============================================================
    16. float NumberStepsFunction_Checker (vec2 current, vec2 stepValue)
    17. {
    18.     // if it's already in the shape, no steps to take
    19.     if (BooleanFunction_Checker(current))
    20.         return 0.0;
    22.     // there are different values to reach based on where we are in the pattern.
    23.     float lower = step(0.5, current.x) * 0.5;
    24.     float upper = lower + 0.5;
    26.     // see how long to escape the box on each axis, take sooner event
    27.     float stepsX;
    28.     if (stepValue.x < 0.0)
    29.         stepsX = (current.x - lower) / -stepValue.x;
    30.     else
    31.         stepsX = (upper - current.x) / stepValue.x;
    33.     float stepsY;
    34.     if (stepValue.y < 0.0)
    35.         stepsY = (current.y - lower) / -stepValue.y;
    36.     else
    37.         stepsY = (upper - current.y) / stepValue.y;    
    39.     return ceil(min(stepsX, stepsY));  
    40. }
    Code (CSharp):
    2. //============================================================
    3. bool BooleanFunction_Square (vec2 current)
    4. {
    5.     return (current.x < 0.25 || current.x > 0.75 || current.y < 0.25 || current.y > 0.75);
    6. }
    8. //============================================================
    9. float NumberStepsFunction_Square (vec2 current, vec2 stepValue)
    10. {
    11.     // if it's already in the shape, no steps to take
    12.     if (BooleanFunction_Square(current))
    13.         return 0.0;
    15.     float stepsX;
    16.     if (stepValue.x < 0.0)
    17.         stepsX = (current.x - 0.25) / -stepValue.x;
    18.     else
    19.         stepsX = (0.75 - current.x) / stepValue.x;
    21.     float stepsY;
    22.     if (stepValue.y < 0.0)
    23.         stepsY = (current.y - 0.25) / -stepValue.y;
    24.     else
    25.         stepsY = (0.75 - current.y) / stepValue.y;  
    27.     return ceil(min(stepsX, stepsY));
    28. }
    Code (CSharp):
    2. //============================================================
    3. bool BooleanFunction_Circle (vec2 current)
    4. {
    5.     const float radius = 0.25;
    6.     const float radiusSq = radius * radius;
    7.     vec2 rel = current - vec2(0.5, 0.5);
    8.     return rel.x*rel.x + rel.y*rel.y > radiusSq;
    9. }
    11. //============================================================
    12. float NumberStepsFunction_Circle (vec2 current, vec2 stepValue)
    13. {
    14.     // if it's already in the shape, no steps to take
    15.     if (BooleanFunction_Circle(current))
    16.         return 0.0;
    18.     float steps = RayIntersectCircle(current, stepValue, vec3(0.5, 0.5, 0.25));
    19.     return ceil(steps);
    20. }
    Code (CSharp):
    2. //============================================================
    3. bool BooleanFunction_ThinSquare (vec2 current)
    4. {
    5.     return (current.x < 0.05 || current.x > 0.95 || current.y < 0.05 || current.y > 0.95);
    6. }
    8. //============================================================
    9. float NumberStepsFunction_ThinSquare (vec2 current, vec2 stepValue)
    10. {
    11.     // if it's already in the shape, no steps to take
    12.     if (BooleanFunction_ThinSquare(current))
    13.         return 0.0;
    15.     float stepsX;
    16.     if (stepValue.x < 0.0)
    17.         stepsX = (current.x - 0.05) / -stepValue.x;
    18.     else
    19.         stepsX = (0.95 - current.x) / stepValue.x;
    21.     float stepsY;
    22.     if (stepValue.y < 0.0)
    23.         stepsY = (current.y - 0.05) / -stepValue.y;
    24.     else
    25.         stepsY = (0.95 - current.y) / stepValue.y;  
    27.     return ceil(min(stepsX, stepsY));
    28. }
    Last edited: May 22, 2020
    Invertex likes this.
  11. neoshaman


    Feb 11, 2011
    Here is the comment of the source above

    I'm trying to solve specifically this,
    This limitation only allow for negative shapes (hole in a plane) instead of positive...
    However an observation is that they do a lossy modulo, they don't keep the original vector, I don't know if that help, but maybe using a division and a modulo can reconstruct the missing data?

    About my solution on earlier post:
    Turns out the math I did, I realized was an extension of the DDA algorithm to corner of grid, so the final algorithm is likely to be a généralization of the DDA.

    DDA main part is in the form of

    cell size / normal vector component

    and return the closest plane

    My algorithm is

    common period / delta rate

    but can both can be interpreted as

    period / interception velocity

    which make sense, it mean how long before we catch up to a point for a given velocity, the main difference is that to intercept corner we scale the grid to the common denominator of the grid and velocity, that is when they overlap. DDA can be understood as a function to find the periodic intersection of planes.

    Now the main problem is how do I handle interval testing?

    The problem is we have 3 cases:
    - the stepping is bigger than the interval and jump over it
    - the stepping hit one time in any point inside the interval
    - the stepping hit multiple time inside the interval.

    So If we consider stepping as interval too, it's about finding recurrent interval overlap.

    Some random search that may or may not help:, the greatest common,8 and 12 is 4.
    Last edited: Jun 13, 2020
    PROE_ likes this.
  12. neoshaman


    Feb 11, 2011

    I feel like i'm so close on that case but missing one insight to make it gel ... :confused:

  13. neoshaman


    Feb 11, 2011

    That's a visual that's not math, how do I turn that into math?

    The general idea is that I find a way to turn the problem into a 1d interval problem, ALMOST ...
    The gist:
    - intersection with circle is done by projecting the center to the ray line (using a simple dot product)
    - that give you an intersection point with the perpendicular of the ray that goes through the center
    - because all space are basically the same, the line (the projection line) that goes through the intersection and the center of the circle is the same in all space, it's constant as it has the same direction and goes through the same point
    - by only wrapping one axis (here the vertical axis) I make evident all interval, (here projected on teh horizontal axis), the cell size , the radius and the vertical intercept of the projection line are all nested constant.
    - but the vertical rate of the ray casted have differing offset, isn't tied to the frequency of the cell size BUT the ray vertical rate is constant.

    Which mean that the recurrent intersection of the ray with the fixed projection line is the key to solving the problem, IF the delta between intersection is constant, it's basically just a distance equation from the first intersection to the center of the circle minus the radius and divided by the delta step to get how many step before first hit.

    The magic is to find that delta formula, given we have a lot of triangle rectangle, it's basic trigonometry or Thales theorem ... should be ... if someone figure it out
    JoeStrout likes this.
  14. neoshaman


    Feb 11, 2011
    Studying the representation above I realized that while I was abusing symmetry to stay in a single gradient I could go on further, and only have rising rays (black angled lines), since the circle axis (pink lines) are always perpendicular we know they are intersected at regular rate by the rays, at least in a single span (in between two projection of the vertical wrapping), the thing is we can see that we kinda can predict where the "bubble" are, ie basically place where the ray will hit or miss in succession. That make sense, in some way the ray is kind of scanning the circle axis at a regular rates, we are jumping through the axis at fixed interval, which mean base on the first intersection we can "march" analytically on the circle axis of a single wrapping interval to define when there is a hit and how many, basically dividing the distance of the intersection to the radius divided by the step size, with an extra check to see if we jump over the primitive interval.

    If we can generalize that to get the set of all interval, we are set, but the main observation is that it looks like, once projected on the horizontal axis, to a problem similar to interference pattern:é_pattern there might be some math useful there
  15. lenneth4


    Jun 20, 2013
    Cant help bro, but hope you managed to find something !
  16. neoshaman


    Feb 11, 2011
    I forgot to update:

    lilacsky824 and Invertex like this.
  17. neoshaman


    Feb 11, 2011
    Some javascript to test a theory that hit me while trying to simplify further.


    - I find a new way to predict the recurrence of primitive intersection, using "empty skip" to solve the first hit problem before one wrapping of the ray.
    - Now I'm looking to solve the wrapping occurrence. It's also potential to be generalize to any primitives and at least greatly accelerate any recurring shape. It might also help for "moving" positions based on the grid position.
  18. neoshaman


    Feb 11, 2011
    SO turning this into a 1d problem I have:

    - Range B > Range S
    - B = offset + kS + (S/d)
    - rn = S-(S/dn-1)
    - offsetn = rn-1

    - offset0 = 0

    - B represent the ray interval along the minor axis trick. That is I wrapped vertically, with horizontal representing the axis of the largest component, aka the shallowest slope.
    - S represent the recurring grid interval.

    I'm trying to determine occurrences as the interval drift relative to each other. At 0 they start at the same, but end with a slight difference, which is compounded at each respective recurrence, until they eventually sync back, if they ever do.

    - k is the number of time S repeat within B, so basically B/S as an integer division.
    - S/d = the fragment of S needed to add to kS to reach size of B, that at time 0 B = kS + (S/d)
    - r is the reminder of B - (k+1)S, and basically the start of a new recurrence offset.

    I need to predict offsetn and (S/dn), let's call the later dsn, and I cannot guarantee than k is fixed for all n.

    Now things looks funny, is that offsetn and dsn are dependent on the previous state of each others. Worse there seems to be funny behavior depending on whether ds is above or below 0.5.

    I'm not good enough about math to recognize a pattern, or known conjoncture or anything. More research are needed.

    - r is teh reminder
  19. neoshaman


    Feb 11, 2011
    Imagine you have:

    case 1
    - a car A going at 30km/h
    - a car B going at 60km/h 1km behind
    how long for car B to reach car A?

    case 2
    - a series of cars A going at 30km/h but spaced at regular interval
    - a car B going at a speed S in between two car A
    For what speed car B reach any car A?

    case 3
    - a series of cars A going at 30km/h but spaced at regular interval
    - Car B is not a car, but a frog that can jump over the cars A, with jump of size J
    How many jump before the frog land on a car for any position?

    That has been my logic progression to try to solve this, I struggle a bit on part 3 though
    Code (CSharp):
    1. function drawFastHits(){
    2.         const scale = GlobalStruct.gridScale
    4.         const start = collisionData.fasthit.start
    5.         const end   = collisionData.fasthit.end
    7.         // const intEnd = floor(end)
    8.         const intEnd   = ceil(end)
    9.         const intStart = ceil(start)
    10.         //loop fastHitList
    11.             //draw markers
    13.         const overlap = intEnd - intStart
    15.         //taking box size
    16.         const rayboxSize    = collisionData.raybox.x2 - collisionData.raybox.x1
    17.         const circleboxSize = collisionData.circlebox.x2 - collisionData.circlebox.x1
    18.         //expending the box size together to test against single point
    19.         const minkowski     = rayboxSize + circleboxSize
    20.         const halfCirclebox = circleboxSize/2
    22.         const marker1 = {x:intStart*scale,y:50}
    23.         const marker2 = {x:intEnd*scale,y:50}
    24.         const marker3 = {x:start*scale,y:75}
    25.         const marker4 = {x:end*scale,y:75}
    29.         color("black")
    30.         vertical(start *scale)
    31.         vertical(end   *scale)
    33.         color("red")
    34.         horizontalLine(50, intStart*scale, intEnd *scale )
    36.         line(marker1.x,marker1.y, marker3.x,marker3.y)
    37.         line(marker2.x, marker2.y, marker4.x,marker4.y)
    39.         drawMarker(marker1)
    40.         drawMarker(marker2)
    42.         drawMarker(marker3)
    43.         drawMarker(marker4)
    45.         text(100,200,end)
    46.         text(10,200,intEnd)
    48.         text(100,210,start)
    49.         text(10,210,intStart)
    51.         text(10,220,overlap)
    53.         text(10,240,"ray: "+rayData.length)
    54.         text(10,250,"raybox: "+rayboxSize)
    55.         text(10,260,"minkowski: "+minkowski)
    57.         // const rayboxSize = collisionData.raybox.x2 - collisionData.raybox.x1
    58.         const gap = collisionData.raybox.x2 + rayData.length - rayboxSize
    59.         //verticalLine(gap *scale,0,100)
    61.         //----
    62.         //minkowski box
    63.         color("pink")
    64.         //draw the minkowski as a slab
    65.         //draw the minkowski as interger offset
    67.         //compute start and end of minkowski box
    68.         const minkowski_start = collisionData.raybox.x1 - halfCirclebox
    69.         const minkowski_end   = collisionData.raybox.x2 + halfCirclebox
    70.         const minkowski_size  = minkowski_end - minkowski_start
    72.         //draw minkowski line
    73.         horizontalLine(60, minkowski_start*scale, minkowski_end *scale )
    75.         //find next offset
    76.         const next = minkowski_end + rayData.length - minkowski_size +halfCirclebox
    77.         //draw next offset
    78.         verticalLine(next *scale,0,100)
    80.         //minkowski delta (snap to integer)
    81.         //if minkowski < 1 (if bigger always overlap circle box)
    82.         if (minkowski_size < 1){
    83.             color("green")
    85.             //this delta seems wrong
    87.             //I want the end of minkowski at 0
    88.             //vs the end of minkowski a 1
    89.             //that is 0+1 length
    91.             //0 = min start (gap before minkowski)
    92.             //+ minkowski size
    93.             const delta_start   = (minkowski_start + minkowski_size)%1
    94.             //1 = gap after minkowski + gap avant minkowski
    95.             //+ minkowski (basically the whole length)
    96.             const delta_end     = (minkowski_end + rayData.length)%1
    97.             const delta         = (delta_end - delta_start)
    99.             // const delta_start   = (minkowski_start - 0.5)%1
    100.             // const delta_end     = (minkowski_end   - 0.5)%1
    101.             // const delta
    104.             const sign = sgn(delta)
    106.             // distance to next circle box //if 0 no hit
    107.             const distance = (rayData.length-minkowski_size)/abs(delta)
    108.             //(how many time in integer
    109.             //there is delta (aka space to cross))
    110.             text(10,280,"delta: "+delta + " s: "+1/delta)
    111.             text(10,290,"sign: " +sign)
    112.             text(10,300,"distance delta: "+distance)
    114.             //every delta, the ray box get closer
    115.             //to the next relevant circle box
    117.             //actually not raydata length, just gap 1+2
    118.             //aka length - minkowski size
    119.             const next_hit = minkowski_end+floor(distance * (rayData.length-minkowski_size))//+0.5 -1
    120.             text(10,310,"next: "+next_hit)
    121.             text(10,320,"length: "+rayData.length)
    123.             //draw next hit
    124.             verticalLine(next_hit *scale,0,1 *scale)
    126.             //next hit should be delta * times in distance
    128.             //draw raylength filling delta
    130.             //---------
    131.             sign > 0 ? color("magenta") : color("aquamarine")
    132.             const inter = horizontalIntercept(mice,1)
    133.             const markerx = {x:inter*scale, y:100}
    134.             drawMarker(markerx)
    135.             //----------
    136.         }
    137.         // const gap2 = collisionData.raybox.x2 + rayData.length - rayboxSize
    138.         // verticalLine(gap2 *scale,0,100)
    140.         //rectangle(x1*scale,y1*scale, x2*scale,y2*scale, false)
    142.     }
    Is there anyone whom can help? :eek:

    And I still have to solve the sin interception on the cylinder,
    which mean tracing the hit trajectories on the cylinder,
    where probably the sin curve reappear but in distorted form
    I just hope that form can be solved by taking advantage of the diagonal trick (ie x=y)

    Meanwhile this can always be used as a starting point too.

    Last edited: Sep 8, 2023