Recently I've needed a Binary Space Partitioning (BSP) tree and I was surprised that there was not a "C++ container-ish" implementation available. I've decided to make my own with custom allocator support and to work in any number of dimensions alongside. This is my first ever code dealing with custom allocators. Is there anything I missed in terms of C++ best practices? Or even potential UB? Is there anything I've missed to be a proper STL container? (I've mostly looked at `std::map`, as that's usually a tree)

```
#define FWD(...) ::std::forward<decltype(__VA_ARGS__)>(__VA_ARGS__)

/**
* A view for BSP tree to simplify iterations.
* Basically an iterator begin and end pair.
*/
template <typename It>
class bsptree_view {
private:
    It from;
    It to;

public:
    constexpr bsptree_view(It from, It to) noexcept
        : from(from), to(to) {}

    [[nodiscard]] constexpr auto begin() noexcept {
        return from;
    }

    [[nodiscard]] constexpr auto end() noexcept {
        return to;
    }
};

/**
* An N dimensional point in space.
*/
template <std::size_t Dim, typename T>
struct point {
    std::array<T, Dim> coordinates;
};

/**
* An N dimensional axis-aligned bounds type, consisting of the offset and size in each dimension.
*/
template <std::size_t Dim, typename T>
struct bounds {
    std::array<T, Dim> offset;
    std::array<T, Dim> size;

    /**
    * Checks if a point is contained within this bounds.
    */
    [[nodiscard]]
    constexpr bool contains(point<Dim, T> const& p) const noexcept {
        for (std::size_t i = 0; i < Dim; ++i) {
            if (p.coordinates[i] < offset[i] || p.coordinates[i] >= offset[i] + size[i]) return false;
        }
        return true;
    }

    /**
    * Checks if two bounds intersect.
    */
    [[nodiscard]]
    constexpr bool intersects(bounds const& o) const noexcept {
        for (std::size_t i = 0; i < Dim; ++i) {
            if (offset[i] >= o.offset[i] + o.size[i] || o.offset[i] >= offset[i] + size[i]) return false;
        }
        return true;
    }
};

/**
* A generic, N-dimensional Binary Space Partitioning tree.
*/
template <
    std::size_t Dim,
    typename Key,
    typename T,
    typename Allocator = std::allocator<std::pair<const point<Dim, Key>, T>>
>
class bsptree {
public:
    // Common type definitions for std::map

    using key_type = point<Dim, Key>;
    using bounds_type = bounds<Dim, Key>;
    using mapped_type = T;
    using value_type = std::pair<const key_type, T>;
    using size_type = std::size_t;
    using difference_type = std::ptrdiff_t;
    using allocator_type = Allocator;
    using reference = value_type&;
    using const_reference = value_type const&;
    using pointer = typename std::allocator_traits<Allocator>::pointer;
    using const_pointer = typename std::allocator_traits<Allocator>::const_pointer;

    inline static constexpr size_type subdivision_count = 1 << Dim;

private:
    using storage_type = std::aligned_storage_t<sizeof(value_type), alignof(value_type)>;

public:
    // The node structure type
    struct node_type {
        // Parent for iteration
        node_type* parent;
        // The bounds of the node
        bounds_type bounds;
        // Stored data
        storage_type* data;
        size_type length = 0;
        // Subnodes
        std::array<node_type*, subdivision_count> subnodes;

        node_type(node_type* parent, bounds_type const& bounds, storage_type* data)
            : parent(parent), bounds(bounds), data(data) {
            for (size_type i = 0; i < subdivision_count; ++i) subnodes[i] = nullptr;
        }

        [[nodiscard]] value_type* values() noexcept {
            return std::launder(reinterpret_cast<value_type*>(data));
        }

        [[nodiscard]] value_type const* values() const noexcept {
            return std::launder(reinterpret_cast<value_type const*>(data));
        }

        [[nodiscard]] value_type& value(size_type i) noexcept {
            return values()[i];
        }

        [[nodiscard]] value_type const& value(size_type i) const noexcept {
            return values()[i];
        }
    };

private:
    // A simple helper for iterating through the node structure
    template <typename Derived>
    class node_stepper {
    protected:
        node_type* current;

        explicit constexpr node_stepper(node_type* current) noexcept
            : current(current) {}

        [[nodiscard]] Derived& underlying() { return static_cast<Derived&>(*this); }
        [[nodiscard]] Derived const& underlying() const { return static_cast<Derived const&>(*this); }

        [[nodiscard]] constexpr bool accepts_node(node_type* n) {
            return n != nullptr && underlying().intersects_bounds(n->bounds);
        }

        void descend_first_node() {
            if (current == nullptr) return;
            while (true) {
            start:
                for (size_type i = 0; i < subdivision_count; ++i) {
                    if (accepts_node(current->subnodes[i])) {
                        current = current->subnodes[i];
                        goto start;
                    }
                }
                break;
            }
        }

        bool move_next_node() {
            if (current == nullptr) return false;
        begin:
            auto parent = current->parent;
            if (parent == nullptr) return false;

            for (size_type i = 0; i + 1 < subdivision_count; ++i) {
                if (current == parent->subnodes[i]) {
                    for (size_type j = i + 1; j < subdivision_count; ++j) {
                        if (accepts_node(parent->subnodes[j])) {
                            current = parent->subnodes[j];
                            descend_first_node();
                            return true;
                        }
                    }
                    break;
                }
            }

            // This was the last node we visited, now stay at the parent
            current = parent;
            if (!accepts_node(current)) goto begin;
            return true;
        }

        bool move_prev_node() {
            // Step down in reverse order
            for (size_type i = subdivision_count; i > 0; --i) {
                if (accepts_node(current->subnodes[i - 1])) {
                    current = current->subnodes[i - 1];
                    return true;
                }
            }
            // No more children, step up
            // NOTE: We want to keep the child if there are no more nodes, so we don't copy back until found
            auto next = current;
            while (true) {
                auto parent = next->parent;
                if (parent == nullptr) return false; // We hit the start
                // Now we need to step back in reverse order
                for (size_type i = subdivision_count; i > 1; --i) {
                    if (next == parent->subnodes[i - 1]) {
                        for (size_type j = i - 1; j > 0; --j) {
                            if (accepts_node(parent->subnodes[j - 1])) {
                                current = parent->subnodes[j - 1];
                                return true;
                            }
                        }
                        break;
                    }
                }
                // We gotta step up more
                next = parent;
            }
            return false;
        }
    };

    // Base for iterators to factor out stepping
    template <typename Derived>
    class value_iterator_base : public node_stepper<Derived> {
    private:
        using base_type = node_stepper<Derived>;

        friend class bsptree;

        template <typename>
        friend class value_iterator_base;

    public:
        using difference_type = typename bsptree::difference_type;
        using iterator_category = std::bidirectional_iterator_tag;

    protected:
        size_type index;

        constexpr value_iterator_base(node_type* current, size_type index) noexcept
            : base_type(current), index(index) {
        }

        [[nodiscard]] constexpr bool is_end() const noexcept {
            return this->current == nullptr
                || (this->current->parent == nullptr && index >= this->current->length);
        }

        [[nodiscard]] bool accepts_value() const noexcept {
            return this->current != nullptr
                && index < this->current->length
                && this->underlying().contains_point(value().first);
        }

    public:
        void descend_first() {
            this->descend_first_node();
            if (!accepts_value()) move_next();
        }

        void move_next() {
        begin:
            if (this->current == nullptr) return;
            if (index < this->current->length) {
                index += 1;
                if (!accepts_value()) goto begin;
                return;
            }
            if (this->move_next_node()) {
                index = 0;
                if (!accepts_value()) goto begin;
            }
        }

        void move_prev() {
        begin:
            if (this->current == nullptr) return;
            if (index > 0) {
                --index;
                if (!accepts_value()) goto begin;
                return;
            }
            if (this->move_prev_node()) {
                index = this->current->length - 1;
                if (!accepts_value()) goto begin;
            }
        }

        [[nodiscard]] typename bsptree::value_type& value() noexcept {
            return this->current->value(index);
        }

        [[nodiscard]] typename bsptree::value_type const& value() const noexcept {
            return this->current->value(index);
        }

        Derived& operator++() {
            this->move_next();
            return this->underlying();
        }

        Derived operator++(int) {
            auto cpy = this->underlying();
            this->operator++();
            return cpy;
        }

        Derived& operator--() {
            this->move_prev();
            return this->underlying();
        }

        Derived operator--(int) {
            auto cpy = this->underlying();
            this->operator--();
            return cpy;
        }

    public:
        template <typename OtherDerived>
        [[nodiscard]]
        constexpr bool operator==(value_iterator_base<OtherDerived> const& r) noexcept {
            return (is_end() && r.is_end())
                || (this->current == r.current && index == r.index);
        }

        template <typename OtherDerived>
        [[nodiscard]]
        constexpr bool operator!=(value_iterator_base<OtherDerived> const& r) noexcept {
            return !(*this == r);
        }
    };

    using node_allocator = typename std::allocator_traits<Allocator>::template rebind_alloc<node_type>;

    bounds_type root_bounds;
    Allocator alloc;
    size_type bucket_size;
    node_type* root;
    size_type count;

    // Initializes all members
    bsptree(
        bounds_type const& bounds, Allocator const& alloc, size_type bucket, node_type* root, size_type count)
        : root_bounds(bounds), alloc(alloc), bucket_size(bucket), root(root), count(count) {}

    // Allocation

    node_type* allocate_node(node_type* parent, bounds_type const& bounds) {
        // Allocate the node
        node_allocator node_alloc = node_allocator(std::move(alloc));
        node_type* result = node_alloc.allocate(1);
        // Allocate bucket
        alloc = Allocator(std::move(node_alloc));
        value_type* data = alloc.allocate(bucket_size);
        // Initialize the node
        new (result) node_type(parent, bounds, std::launder(reinterpret_cast<storage_type*>(data)));
        return result;
    }

    void deallocate_node(node_type* n) {
        // First destruct every element
        for (size_type i = 0; i < n->length; ++i) {
            n->value(i).~value_type();
        }
        // Free the space for the data
        alloc.deallocate(std::launder(reinterpret_cast<value_type*>(n->data)), bucket_size);
        auto parent = n->parent;
        if (parent != nullptr) {
            // Null out the proper subnode of parent
            for (size_type i = 0; i < subdivision_count; ++i) {
                if (parent->subnodes[i] == n) {
                    parent->subnodes[i] = nullptr;
                    break;
                }
            }
        }
        // Now we can free the node itself
        node_allocator node_alloc = node_allocator(std::move(alloc));
        n->~node_type();
        node_alloc.deallocate(n, 1);
        alloc = Allocator(std::move(node_alloc));
    }

    void destruct() {
        auto current = root;
        while (current) {
            for (size_type i = 0; i < subdivision_count; ++i) {
                if (current->subnodes[i] != nullptr) { 
                    current = current->subnodes[i]; 
                    continue;
                }
            }
            // This is a leaf, we can free it
            auto parent = current->parent;
            deallocate_node(current);
            current = parent;
        }
        root = nullptr;
        count = 0;
    }

    node_type* clone_node(node_type* other, node_type* parent) {
        auto n = allocate_node(parent, other->bounds);
        for (; n->length < other->length; ++n->length) {
            new (&n->data[n->length]) value_type(other->value(n->length));
        }
        return n;
    }

    void clone(bsptree const& other) {
        if (other.root == nullptr) return;
        // We at least have a root
        root = clone_node(other.root, nullptr);

        auto current = root;
        auto other_current = other.root;

        while (true) {
            // Clone the first non-null child
            for (size_type i = 0; i < subdivision_count; ++i) {
                if (current->subnodes[i] == nullptr && other_current->subnodes[i] != nullptr) {
                    current->subnodes[i] = clone_node(other_current->subnodes[i], current);
                    current = current->subnodes[i];
                    other_current = other_current->subnodes[i];
                    continue;
                }
            }
            // Either a leaf or cloned everything, step up
            if (current->parent == nullptr) break;
            current = current->parent;
            other_current = other_current->parent;
        }
    }

public:
    class iterator : public value_iterator_base<iterator> {
    private:
        using base_type = value_iterator_base<iterator>;

    public:
        using base_type::difference_type;
        using base_type::iterator_category;
        using reference = typename bsptree::reference;
        using pointer = typename bsptree::pointer;
        using value_type = typename bsptree::value_type;

        using base_type::base_type;

        using base_type::operator++;
        using base_type::operator--;

        [[nodiscard]] reference operator*() noexcept {
            return this->value();
        }

        [[nodiscard]] pointer operator->() noexcept {
            return &this->value();
        }

        [[nodiscard]] constexpr bool intersects_bounds(bounds_type const& bounds) const noexcept {
            return true;
        }

        [[nodiscard]] constexpr bool contains_point(key_type const& point) const noexcept {
            return true;
        }
    };

    class const_iterator : public value_iterator_base<const_iterator> {
    private:
        using base_type = value_iterator_base<const_iterator>;

    public:
        using base_type::difference_type;
        using base_type::iterator_category;
        using reference = typename bsptree::const_reference;
        using pointer = typename bsptree::const_pointer;
        using value_type = typename bsptree::value_type;

        using base_type::base_type;

        constexpr const_iterator(iterator it) noexcept
            : const_iterator(it.current, it.index) {}

        using base_type::operator++;
        using base_type::operator--;

        [[nodiscard]] reference operator*() noexcept {
            return *this->value();
        }

        [[nodiscard]] pointer operator->() noexcept {
            return this->value();
        }

        [[nodiscard]] constexpr bool intersects_bounds(bounds_type const& bounds) const noexcept {
            return true;
        }

        [[nodiscard]] constexpr bool contains_point(key_type const& point) const noexcept {
            return true;
        }
    };

    template <typename Shape>
    class query_iterator : public value_iterator_base<query_iterator<Shape>> {
    private:
        using base_type = value_iterator_base<query_iterator<Shape>>;

        Shape shape;

    public:
        using base_type::difference_type;
        using base_type::iterator_category;
        using reference = typename bsptree::reference;
        using pointer = typename bsptree::pointer;
        using value_type = typename bsptree::value_type;

        constexpr query_iterator(node_type* current, size_type index, Shape const& shape) noexcept
            : base_type(current, index), shape(shape) {}

        using base_type::operator++;
        using base_type::operator--;

        [[nodiscard]] reference operator*() noexcept {
            return this->value();
        }

        [[nodiscard]] pointer operator->() noexcept {
            return &this->value();
        }

        [[nodiscard]] constexpr bool intersects_bounds(bounds_type const& bounds) const noexcept {
            return shape.intersects(bounds);
        }

        [[nodiscard]] constexpr bool contains_point(key_type const& point) const noexcept {
            return shape.contains(point);
        }
    };

    template <typename Shape>
    class const_query_iterator : public value_iterator_base<const_query_iterator<Shape>> {
    private:
        using base_type = value_iterator_base<const_query_iterator<Shape>>;

        Shape shape;

    public:
        using base_type::difference_type;
        using base_type::iterator_category;
        using reference = typename bsptree::const_reference;
        using pointer = typename bsptree::const_pointer;
        using value_type = typename bsptree::value_type;

        constexpr const_query_iterator(node_type* current, size_type index, Shape const& shape) noexcept
            : base_type(current, index), shape(shape) {}

        constexpr const_query_iterator(query_iterator<Shape> it) noexcept
            : const_query_iterator(it.current, it.index, it.shape) {}

        using base_type::operator++;
        using base_type::operator--;

        [[nodiscard]] reference operator*() noexcept {
            return this->value();
        }

        [[nodiscard]] pointer operator->() noexcept {
            return &this->value();
        }

        [[nodiscard]] constexpr bool intersects_bounds(bounds_type const& bounds) const noexcept {
            return shape.intersects(bounds);
        }

        [[nodiscard]] constexpr bool contains_point(key_type const& point) const noexcept {
            return shape.contains(point);
        }
    };

    class node_iterator : public node_stepper<node_iterator> {
    private:
        friend class bsptree;

        using base_type = node_stepper<node_iterator>;

    public:
        using difference_type = typename bsptree::difference_type;
        using iterator_category = std::bidirectional_iterator_tag;
        using value_type = node_type;
        using reference = value_type const&;
        using pointer = value_type const*;

        explicit constexpr node_iterator(node_type* current) noexcept
            : base_type(current) {}

        [[nodiscard]] constexpr bool intersects_bounds(bounds_type const& bounds) const noexcept {
            return true;
        }

        [[nodiscard]] reference operator*() noexcept {
            return *this->current;
        }

        [[nodiscard]] pointer operator->() noexcept {
            return this->current;
        }

        node_iterator& operator++() {
            this->move_next_node();
            return *this;
        }

        node_iterator operator++(int) {
            auto cpy = *this;
            this->operator++();
            return cpy;
        }

        node_iterator& operator--() {
            this->move_prev_node();
            return *this;
        }

        node_iterator operator--(int) {
            auto cpy = *this;
            this->operator--();
            return cpy;
        }

    public:
        [[nodiscard]]
        constexpr bool operator==(node_iterator const& r) noexcept {
            return this->current == r.current;
        }

        [[nodiscard]]
        constexpr bool operator!=(node_iterator const& r) noexcept {
            return !(*this == r);
        }
    };

    explicit bsptree(bounds_type const& bounds, size_type bucket = 1)
        : bsptree(bounds, bucket, Allocator()) {}

    bsptree(bounds_type const& bounds, Allocator const& alloc)
        : bsptree(bounds, 1, alloc) {}

    bsptree(bounds_type const& bounds, size_type bucket, Allocator const& alloc)
        : bsptree(bounds, alloc, bucket, nullptr, 0) {}

    bsptree(bsptree&& other, Allocator const& alloc)
        : bsptree(other.root_bounds, alloc, other.bucket_size, other.root, other.count) {
        other.root = nullptr;
        other.count = 0;
    }

    bsptree(bsptree&& other)
        : bsptree(std::move(other), other.data_alloc) {}

    ~bsptree() {
        destruct();
    }

    bsptree& operator=(bsptree&& other) noexcept {
        destruct();
        root_bounds = other.root_bounds;
        alloc = other.alloc;
        bucket_size = other.bucket_size;
        root = other.root;
        count = other.count;
        other.root = nullptr;
        other.count = 0;
        return *this;
    }

    bsptree(bsptree const& other)
        : bsptree(other, other.alloc) {}

    bsptree(bsptree const& other, Allocator const& alloc)
        : bsptree(other.root_bounds, alloc, other.bucket_size, nullptr, other.count) {
        clone(other);
    }

    bsptree& operator=(bsptree const& other) {
        destruct();
        root_bounds = other.root_bounds;
        alloc = other.alloc;
        bucket_size = other.bucket_size;
        count = other.count;
        clone(other);
        return *this;
    }

    [[nodiscard]] allocator_type get_allocator() const noexcept { return alloc; }
    [[nodiscard]] node_type* root_node() const noexcept { return root; }

    // TODO: Element access?

    // Iterators

    [[nodiscard]] iterator begin() noexcept {
        auto it = iterator(root, 0);
        if (root != nullptr) it.descend_first();
        return it;
    }

    [[nodiscard]] iterator end() noexcept {
        if (root == nullptr) return iterator(root, 0);
        return iterator(root, root->length);
    }

    [[nodiscard]] const_iterator begin() const noexcept {
        return cbegin();
    }

    [[nodiscard]] const_iterator end() const noexcept {
        return cend();
    }

    [[nodiscard]] const_iterator cbegin() const noexcept {
        auto it = const_iterator(root, 0);
        if (root != nullptr) it.descend_first();
        return it;
    }

    [[nodiscard]] const_iterator cend() const noexcept {
        if (root == nullptr) return const_iterator(root, 0);
        return const_iterator(root, root->length);
    }

    template <typename Shape>
    [[nodiscard]] query_iterator<Shape> query_begin(Shape const& shape) noexcept {
        auto it = query_iterator<Shape>(root, 0, shape);
        if (root != nullptr) it.descend_first();
        return it;
    }

    template <typename Shape>
    [[nodiscard]] query_iterator<Shape> query_end(Shape const& shape) noexcept {
        if (root == nullptr) return query_iterator<Shape>(root, 0, shape);
        return query_iterator<Shape>(root, root->length, shape);
    }

    template <typename Shape>
    [[nodiscard]] const_query_iterator<Shape> query_begin(Shape const& shape) const noexcept {
        return cquery_begin(shape);
    }

    template <typename Shape>
    [[nodiscard]] const_query_iterator<Shape> query_end(Shape const& shape) const noexcept {
        return cquery_end(shape);
    }

    template <typename Shape>
    [[nodiscard]] const_query_iterator<Shape> cquery_begin(Shape const& shape) const noexcept {
        auto it = const_query_iterator<Shape>(root, 0, shape);
        if (root != nullptr) it.descend_first();
        return it;
    }

    template <typename Shape>
    [[nodiscard]] const_query_iterator<Shape> cquery_end(Shape const& shape) const noexcept {
        if (root == nullptr) return const_query_iterator<Shape>(root, 0, shape);
        return const_query_iterator<Shape>(root, root->length, shape);
    }

    [[nodiscard]] node_iterator nodes_begin() const noexcept {
        auto it = node_iterator(root);
        if (root != nullptr) it.descend_first_node();
        return it;
    }

    [[nodiscard]] node_iterator nodes_end() const noexcept {
        return node_iterator(root);
    }

    // Views

    template <typename Shape>
    [[nodiscard]] bsptree_view<query_iterator<Shape>> query(Shape const& shape) & noexcept {
        return { query_begin(shape), query_end(shape) };
    }

    template <typename Shape>
    [[nodiscard]] bsptree_view<const_query_iterator<Shape>> query(Shape const& shape) const& noexcept {
        return cquery(shape);
    }

    template <typename Shape>
    [[nodiscard]] bsptree_view<const_query_iterator<Shape>> cquery(Shape const& shape) const& noexcept {
        return { cquery_begin(shape), cquery_end(shape) };
    }

    [[nodiscard]] bsptree_view<node_iterator> nodes() const& noexcept {
        return { nodes_begin(), nodes_end() };
    }

    // Capacity

    [[nodiscard]] bool empty() const noexcept { return count == 0; }
    [[nodiscard]] size_type size() const noexcept { return count; }
    [[nodiscard]] size_type max_size() const noexcept { return std::numeric_limits<difference_type>::max(); }
    [[nodiscard]] bounds_type const& bounds() const noexcept { return root_bounds; }

    // Modifiers

    // TODO: Remaining?

    void clear() noexcept {
        destruct();
    }

    iterator insert(value_type const& value) {
        return emplace(value);
    }

    template <typename... Args>
    iterator emplace(Args&&... args) {
        // NOTE: Cheeky emplace, as I don't really know how to extract the key from here
        return insert(value_type(FWD(args)...));
    }

    iterator insert(value_type&& value) {
        if (root == nullptr) root = allocate_node(nullptr, root_bounds);
        // Search the first free subnode
        auto subnode = root;
        while (true) {
            if (subnode->length < bucket_size) {
                // Free space in subnode
                new (&subnode->data[subnode->length++]) value_type(std::move(value));
                ++count;
                return iterator(subnode, subnode->length - 1);
            }
            // No free space in current node, find where to traverse
            auto const& bounds = subnode->bounds;
            // Calculate center
            auto center = std::array<Key, Dim>();
            for (size_type i = 0; i < Dim; ++i) center[i] = bounds.offset[i] + bounds.size[i] / Key{2};
            // Search for the next subnode
            // We search by constructing an index bit-wise
            node_type** next_subnode = nullptr;
            bounds_type next_subbounds;
            size_type next_subbounds_index = 0;
            for (size_type i = 0; i < Dim; ++i) {
                if (value.first.coordinates[i] < center[i]) {
                    next_subbounds.offset[i] = bounds.offset[i];
                    next_subbounds.size[i] = center[i] - bounds.offset[i];
                }
                else {
                    next_subbounds_index |= (1 << i);
                    next_subbounds.offset[i] = center[i];
                    next_subbounds.size[i] = bounds.offset[i] + bounds.size[i] - center[i];
                }
            }
            next_subnode = &subnode->subnodes[next_subbounds_index];
            if (*next_subnode == nullptr) {
                // The next subnode is unallocated
                *next_subnode = allocate_node(subnode, next_subbounds);
            }
            subnode = *next_subnode;
        }
    }

    iterator erase(iterator pos) {
        return erase_impl(pos);
    }

    iterator erase(const_iterator pos) {
        auto it = erase_impl(pos);
        return iterator(it.current, it.index);
    }

    template <typename Shape>
    query_iterator<Shape> erase(query_iterator<Shape> pos) {
        return erase_impl(pos);
    }

    template <typename Shape>
    query_iterator<Shape> erase(const_query_iterator<Shape> pos) {
        auto it = erase_impl(pos);
        return query_iterator<Shape>(it.current, it.index, it.shape);
    }

    size_type erase(iterator first, iterator last) {
        return erase_impl(first, last);
    }

    size_type erase(const_iterator first, const_iterator last) {
        return erase_impl(first, last);
    }

    template <typename Shape>
    size_type erase(query_iterator<Shape> first, query_iterator<Shape> last) {
        return erase_impl(first, last);
    }

    template <typename Shape>
    size_type erase(const_query_iterator<Shape> first, const_query_iterator<Shape> last) {
        return erase_impl(first, last);
    }

    template <typename It>
    size_type erase(bsptree_view<It> view) {
        return erase_impl(view.begin(), view.end());
    }

    // TODO: Lookup

private:
    template <typename It>
    size_type erase_impl(It start, It end) {
        size_type removed = 0;
        for (; start != end; ++removed) start = erase_impl(start);
        return removed;
    }

    template <typename It>
    It erase_impl(It pos) {
        // Remove element
        auto node = pos.current;
        auto idx = pos.index;
        pos.value()->~value_type();
        for (size_type i = idx + 1; i < node->length; ++i) {
            new (&node->data[i - 1]) value_type(std::move(node->value(i)));
        }
        --count;
        --node->length;
        if (idx < node->length) {
            // The position is still valid, contains the next element
            return pos;
        }
        // The position is no longer valid, gotta step
        pos.move_next();
        return pos;
    }
};

// Some aliases

template <
    typename Key,
    typename T,
    typename Allocator = std::allocator<std::pair<const point<2, Key>, T>>
>
using quadtree = bsptree<2, Key, T, Allocator>;

template <
    typename Key,
    typename T,
    typename Allocator = std::allocator<std::pair<const point<3, Key>, T>>
>
using octree = bsptree<3, Key, T, Allocator>;

```

Some sample usage code is provided here:

```
// Create a quadtree with float coordinates, with an area from (0,0) to (100, 100), bucket size of 4
// Associated values are Particle pointers
auto qtree = quadtree<float, Particle*>{ {0.0f, 0.0f, 100.0f, 100.0f}, 4 };
// Add a few particles at specific positions
qtree.insert({ { 10.0f, 10.0f }, &p1 });
qtree.insert({ { 10.0f, 30.0f }, &p2 });
qtree.insert({ { 30.0f, 30.0f }, &p3 });
// Let's list all the particles intersecting the bounds starting (5, 5) and dimensions (40, 35)
// Note that any shape can be used that implements contains(point) and intersects(bounds)
for (auto& [pos, particle] : qtree.query(bounds<2, float>{ {5, 5}, {40, 35} })) {
    // TODO: Do something with the particles
}

```