ray_ellipsoid Function

public function ray_ellipsoid(origin, dir, center, axes, radii, distance, point, local) result(hit)

Where the ray (origin, unit dir) first meets the ellipsoid centered at center with semi-axes radii along the columns of axes (unit vectors): the distance along the ray, the point, and the point in the ellipsoid's axes (from its center). False if it misses (or the ellipsoid is behind).

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: origin(3)

the ray (dir: a unit vector)
real(kind=wp), intent(in) :: dir(3)

the ray (dir: a unit vector)
real(kind=wp), intent(in) :: center(3)

the ellipsoid's center
real(kind=wp), intent(in) :: axes(3,3)

its axes (columns: unit vectors)
real(kind=wp), intent(in) :: radii(3)

its semi-axes along them
real(kind=wp), intent(out) :: distance

along the ray (huge() if missed)
real(kind=wp), intent(out) :: point(3)

the point met
real(kind=wp), intent(out) :: local(3)

the point in the ellipsoid's axes, from its center

Return Value logical


Called by

proc~~ray_ellipsoid~~CalledByGraph proc~ray_ellipsoid ray_ellipsoid proc~pick_ellipsoid pick_ellipsoid proc~pick_ellipsoid->proc~ray_ellipsoid proc~model_pick Model%model_pick proc~model_pick->proc~pick_ellipsoid proc~sphere_pick Sphere%sphere_pick proc~sphere_pick->proc~pick_ellipsoid proc~terrain_pick TerrainBody%terrain_pick proc~terrain_pick->proc~pick_ellipsoid

Source Code

    logical function ray_ellipsoid(origin, dir, center, axes, radii, distance, point, local) result(hit)
        real(wp), intent(in) :: origin(3), dir(3)           !! the ray (dir: a unit vector)
        real(wp), intent(in) :: center(3)                   !! the ellipsoid's center
        real(wp), intent(in) :: axes(3,3)                   !! its axes (columns: unit vectors)
        real(wp), intent(in) :: radii(3)                    !! its semi-axes along them
        real(wp), intent(out) :: distance                   !! along the ray (huge() if missed)
        real(wp), intent(out) :: point(3)                   !! the point met
        real(wp), intent(out) :: local(3)                   !! the point in the ellipsoid's axes, from its center
        real(wp) :: o(3), d(3), qa, qb, qc, disc, t
        hit = .false.
        distance = huge(1.0_wp)
        point = 0.0_wp
        local = 0.0_wp
        if (any(radii <= 0.0_wp)) return
        o = matmul(transpose(axes), origin - center) / radii    ! (the unit sphere)
        d = matmul(transpose(axes), dir) / radii
        qa = dot_product(d, d)
        qb = 2.0_wp * dot_product(o, d)
        qc = dot_product(o, o) - 1.0_wp
        disc = qb * qb - 4.0_wp * qa * qc
        if (disc < 0.0_wp .or. qa <= 0.0_wp) return
        t = (-qb - sqrt(disc)) / (2.0_wp * qa)
        if (t < 0.0_wp) t = (-qb + sqrt(disc)) / (2.0_wp * qa)  ! (the eye inside it: where the ray leaves)
        if (t < 0.0_wp) return
        distance = t
        point = origin + t * dir
        local = (o + t * d) * radii
        hit = .true.
    end function ray_ellipsoid