! This code is adapted from the C++ implementation by ! [Lester Hedges](http://lesterhedges.net). The original ! C++ implementation can be found [here](https://github.com/lohedges/aabbcc#readme). ! ! The zlib/libpng License (Zlib) ! ! Copyright (c) 2009 Erin Catto http://www.box2d.org ! Copyright (c) 2016 Lester Hedges <lester.hedges+aabbcc@gmail.com> ! ! This software is provided 'as-is', without any express or implied ! warranty. In no event will the authors be held liable for any damages ! arising from the use of this software. ! ! Permission is granted to anyone to use this software for any purpose, ! including commercial applications, and to alter it and redistribute it ! freely, subject to the following restrictions: ! ! 1. The origin of this software must not be misrepresented; you must not ! claim that you wrote the original software. If you use this software ! in a product, an acknowledgment in the product documentation would be ! appreciated but is not required. ! ! 2. Altered source versions must be plainly marked as such, and must not be ! misrepresented as being the original software. ! ! 3. This notice may not be removed or altered from any source distribution. ! ! This code was adapted from parts of the Box2D Physics Engine, ! http://www.box2d.org #:include 'asserts.fypp' module aabb_m use constants_m use strings_m implicit none type aabb_t real(rp), dimension(3) :: lbnd real(rp), dimension(3) :: ubnd real(rp), dimension(3) :: center real(rp) :: srfarea contains procedure :: init => aabb_init procedure :: print => aabb_print procedure :: clear => aabb_clear procedure :: get_extent => aabb_get_extent procedure :: update => aabb_update procedure :: fatten => aabb_fatten procedure :: includes => aabb_includes procedure :: overlaps => aabb_overlaps procedure, private :: calc_center => aabb_calc_center procedure, private :: calc_srfarea => aabb_calc_srfarea end type aabb_t interface operator(+) module procedure :: aabb_union end interface contains !******************************************************************************* subroutine aabb_init(this, lbnd, ubnd) !! Initializes an *aabb_t* instance from lower and upper bounds. class(aabb_t), intent(out) :: this !! An *aabb_t* instance. real(rp), dimension(3), intent(in) :: lbnd !! Lower bound real(rp), dimension(3), intent(in) :: ubnd !! Upper bound @:ASSERT(all(lbnd <= ubnd)) this%lbnd = lbnd; this%ubnd = ubnd call this%calc_center() call this%calc_srfarea() end subroutine !******************************************************************************* subroutine aabb_print(this, frmt, str) !! Prints an AABB. class(aabb_t), intent(in) :: this !! An *aabb_t* instance. character(len=*), intent(in), optional :: frmt !! Fortran-style format string for a real number. Default: *(g0.6)*. character(len=:), allocatable, intent(out), optional :: str !! If present, the output is printed to this string instead of STDOUT. character(len=:), allocatable :: frmt_, buf frmt_ = '(g0.6)' if (present(frmt)) frmt_ = trim(adjustl(frmt)) buf = 'lbnd: ['//str_from_num(this%lbnd(1), frmt_)// ', ' & &//str_from_num(this%lbnd(2), frmt_) // ', ' & &//str_from_num(this%lbnd(3), frmt_) // ']' buf = buf // ', ubnd: ['//str_from_num(this%ubnd(1), frmt_)// ', ' & &//str_from_num(this%ubnd(2), frmt_) // ', ' & &//str_from_num(this%ubnd(3), frmt_) // ']' buf = buf // ', center: ['//str_from_num(this%center(1), frmt_)// ', ' & &//str_from_num(this%center(2), frmt_) // ', '& &//str_from_num(this%center(3), frmt_) // ']' buf = buf // ', srfarea: ' // str_from_num(this%srfarea, frmt_) if (present(str)) then str = buf else write(*,*) buf end if end subroutine !******************************************************************************* subroutine aabb_clear(this) !! Clears all attributes of an AABB and sets them to zero. class(aabb_t), intent(in out) :: this !! An *aabb_t* instance. this%lbnd = 0.0_rp; this%ubnd = 0.0_rp this%center = 0.0_rp; this%srfarea = 0.0_rp end subroutine !******************************************************************************* subroutine aabb_get_extent(this, extent) !! Calculates the extent of an *aabb*. The *extent* of an AABB is defined as !! the difference between its upper and lower bounds. class(aabb_t), intent(in) :: this !! An *aabb_t* instance. real(rp), dimension(3), intent(out) :: extent !! Extent of an AABB. extent = this%ubnd - this%lbnd end subroutine !******************************************************************************* subroutine aabb_update(this, lbnd, ubnd) !! Updates an AABB with new bounds. class(aabb_t), intent(in out) :: this !! An *aabb_t* instance. real(rp), dimension(3), intent(in), optional :: lbnd !! Lower bound real(rp), dimension(3), intent(in), optional :: ubnd !! Upper bound !Nothing to do if no bounds are present if ( (.not. present(lbnd)) .and. (.not. present(ubnd)) ) return !Update bounds based on input if (present(lbnd)) this%lbnd = lbnd if (present(ubnd)) this%ubnd = ubnd @:ASSERT( all(this%lbnd <= this%ubnd) ) call this%calc_center() call this%calc_srfarea() end subroutine !******************************************************************************* subroutine aabb_fatten(this, frac) !! Fattens an AABB by a fraction of its base extent. class(aabb_t), intent(in out) :: this !! An *aabb_t* instance. real(rp), intent(in) :: frac !! Fraction of AABB base extent. real(rp), dimension(3) :: extent extent = this%ubnd - this%lbnd !New bounds this%lbnd = this%lbnd - frac*extent; this%ubnd = this%ubnd + frac*extent @:ASSERT( all(this%lbnd <= this%ubnd) ) call this%calc_center() call this%calc_srfarea() end subroutine !******************************************************************************* function aabb_includes(this, other) result(res) !! Returns *true* if *this* includes *other*, *false* otherwise. Inclusion !! is considered in a strict sense. class(aabb_t), intent(in) :: this !! An *aabb_t* instance. type(aabb_t), intent(in) :: other !! An *aabb_t* instance. logical :: res res = ( all(this%lbnd < other%lbnd) .and. all(this%ubnd > other%ubnd) ) end function !******************************************************************************* function aabb_overlaps(this, other) result(res) !! Returns *true* if *this* overlaps *other*, *false* otherwise. Touching !! does not count as an overlap. class(aabb_t), intent(in) :: this !! An *aabb_t* instance. type(aabb_t), intent(in) :: other !! An *aabb_t* instance. logical :: res res = ( all(this%ubnd > other%lbnd) .and. all(other%ubnd > this%lbnd) ) end function !******************************************************************************* function aabb_union(x, y) result(z) !! Combines AABBs *x* and *y* to return a new AABB *z*. type(aabb_t), intent(in) :: x !! An *aabb_t* instance. type(aabb_t), intent(in) :: y !! An *aabb_t* instance. type(aabb_t) :: z z%lbnd = merge( x%lbnd, y%lbnd, (x%lbnd < y%lbnd) ) z%ubnd = merge( x%ubnd, y%ubnd, (x%ubnd > y%ubnd) ) call z%calc_center() call z%calc_srfarea() end function !******************************************************************************* subroutine aabb_calc_center(this) !! Calculates the center of an *aabb*. class(aabb_t), intent(in out) :: this !! An *aabb_t* instance. this%center = 0.5_rp*(this%ubnd + this%lbnd) end subroutine !******************************************************************************* subroutine aabb_calc_srfarea(this) !! Calculates the surface area of an *aabb*. class(aabb_t), intent(in out) :: this !! An *aabb_t* instance. real(rp) :: dx, dy, dz dx = this%ubnd(1) - this%lbnd(1) dy = this%ubnd(2) - this%lbnd(2) dz = this%ubnd(3) - this%lbnd(3) this%srfarea = 2.0_rp*(dx*dy + dx*dz + dy*dz) end subroutine !******************************************************************************* end module aabb_m