HIF: Height Interval Filtering for Efficient Dynamic Points Removal

FreeSpace-based removal methods represent point clouds as voxels and update occupancy states using per-frame scan results based on ray-tracing principles [14, 15]. However, sole reliance on ray-tracing leads to significant misclassification issues, particularly affecting ground points, as shown in Fig. 2. To address this, some methods incorporate post-processing techniques to recover ground points after ray-tracing, effectively reducing misclassification errors [16, 17]. Additionally, improvements in trace efficiency and the introduction of more robust discriminative mechanisms have further enhanced both computational efficiency and accuracy [18, 19]. Despite these advancements, computationally expensive 3D retrieval remains a major bottleneck for real-time performance. Our approach mitigates this by partitioning the space into pillars and transforming 3D retrieval into a more efficient 1D lookup using a hash-mixed operation, significantly reducing computational complexity.

Difference-based methods detect dynamic objects by analyzing visibility or disparity changes. One approach involves converting maps into depth images and comparing depth differences between frames and a global map, as demonstrated in [7, 5]. Other methods prioritize vertical information, detecting changes along the vertical axis to identify dynamic objects [8, 9, 10, 6]. Some techniques employ vertical layering, which can lead to rigid height stratification issues. Most difference-based approaches operate as post-processing techniques, requiring all scans as global priors to ensure accuracy, thereby limiting real-time applicability. Additionally, these methods heavily depend on ground segmentation, increasing computational overhead and system complexity. In contrast, our approach models vertical distributions by constructing flexible height intervals, eliminating reliance on ground segmentation and mitigating rigid height stratification issues.

Learning-based methods have shown strong capabilities across various domains. In dynamic object removal, the task is often reformulated as a point cloud segmentation problem [20, 21, 22]. These approaches transform point clouds into range images to achieve a dense and compact representation [20, 21]. More recent methods introduce residual images to incorporate dynamic information, improving segmentation performance [23, 24, 25, 26, 27]. However, learning-based approaches require extensive training data and suffer from data imbalance issues [28], while generalization and robustness remain significant challenges [18]. Furthermore, their high computational cost limits practical deployment.

Recent studies [8, 9, 10, 6] have demonstrated that analyzing vertical information effectively aids in dynamic object identification. Inspired by these works, we adopt a pillar-based representation to capture vertical structures. Within each pillar, we construct height intervals to model spatial occupancy probabilities. To ensure a globally consistent pillar partitioning, we define a fixed origin point $O$ on the XY plane, preventing misalignment between new scans and the global height intervals.

In our approach, we represent the environment using a 2D pillar-based grid where each pillar has a fixed size of $\Delta x\times\Delta y$ in the XY plane. These pillars are used to discretize the space horizontally, while vertical information is modeled by constructing height intervals within each pillar to maintain spatial occupancy probabilities.

At time $t$, the point cloud scan is denoted as $S_{t}=\{\mathbf{p}_{i}=(x_{i},y_{i},z_{i})\mid i=1,2,...,N\}$, where $N$ represents the total number of points in the scan. Each point $\mathbf{p}_{i}\in S_{t}$ is assigned to a corresponding pillar $P_{mn}$ by computing its offset $(m,n)$ relative to the global origin $O=(x_{o},y_{o})$. The offset is determined as follows:

$$m=\left\lfloor\frac{x_{i}-x_{o}}{\Delta x}\right\rfloor,\quad n=\left\lfloor\frac{y_{i}-y_{o}}{\Delta y}\right\rfloor$$

(1)

This scheme ensures the pillar $P{mn}$ boundaries align with:

$$P_{mn}=\left\{\mathbf{p}_{i}\,\bigg{|}\,\begin{aligned} x_{o}+m\Delta x\leq&x_{i}<x_{o}+(m+1)\Delta x\\ y_{o}+n\Delta y\leq&y_{i}<y_{o}+(n+1)\Delta y\end{aligned}\right\}$$

(2)

To enable efficient pillar retrieval, we employ a hash-mixed operation to generate unique pillar indices:

$$H_{mn}=h(m)\oplus(h(n)+\mathrm{C}+(h(m)\ll 6)+(h(m)\gg 2))$$

(3)

where $h(\cdot)$ represents a standard integer hash function, $\mathrm{C}$ is a constant derived from the golden ratio, commonly used in hash functions to improve distribution uniformity. The operators $\ll$ and $\gg$ denote bitwise left and right shifts, respectively.

This hashing strategy ensures an injective mapping, minimizing the risk of collisions. Furthermore, to optimize storage efficiency, we store only non-empty pillars and construct an index for fast access.

Next, our approach dynamically merges height values within each pillar into adaptive height intervals, enabling finer vertical representation without requiring ground segmentation, as illustrated in Fig. 4 A-1 and Fig. 4 A-2. This design addresses the limitations of existing methods, which suffer from rigid vertical partitioning. The R-POD method in [8] and the approaches in [9, 10] rely on predefined height divisions with binary or symbolic encoding, making them less adaptable to elevation changes and occlusions. Additionally, these methods require auxiliary ground segmentation and lack probabilistic modeling, leading to potential misclassifications in complex scenes.

We define $p_{\text{empty}}$ to maintain the occupancy probability of the empty space in each pillar, which provides reasonable initialization for the observed height interval in new scans. Each pillar $P_{mn}$ is structured as follows:

$$P_{mn}=\left(p_{\text{empty}},\,\mathcal{H}_{mn}=\left\{h_{k}=(b_{k},t_{k},p_{k})\right\}_{k=1}^{K}\right)$$

(4)

where the set $\mathcal{H}_{mn}$ represents the collection of height intervals within a pillar, containing $K$ discrete height intervals $h_{k}$. Each interval $h_{k}$ is defined by three parameters: $b_{k}$ and $t_{k}$, which denote the minimum and maximum height values within the interval, respectively, and $p_{k}$, which represents the occupancy probability of the interval.

We represent each pillar using an index $H_{mn}$ and an associated storage variable $P_{mn}$, forming a key-value pair $H_{mn}:P_{mn}$. To maintain a global mapping of pillars, we define a dictionary $\mathcal{G}=\{H_{mn}^{g}:P_{mn}^{g}\}$, where $H_{mn}^{g}$ denotes the spatial index of a pillar in the global map, and $P_{mn}^{g}$ stores its corresponding height intervals and occupancy probabilities. Upon receiving a new scan $S_{t}$ at time $t$, we perform the Height Interval Construction operation to generate a local pillar set $\mathcal{L}^{t}=\{H_{mn}^{t}:P_{mn}^{t}\}$, where $H_{mn}^{t}$ represents the spatial index of a pillar in the current scan, and $P_{mn}^{t}$ contains its associated height intervals and occupancy probabilities.

For each newly generated local pillar pair $H_{mn}^{t}:P_{mn}^{t}$ in $\mathcal{L}^{t}$, we first check whether a corresponding global pillar $H_{mn}^{g}$ exists in $\mathcal{G}$. If a match is found, we update $P_{mn}^{g}$ using $P_{mn}^{t}$ following the probability filtering operations described in Sections III-B and III-C. If no match is found, the new key-value pair $H_{mn}^{t}:P_{mn}^{t}$ is directly inserted into $\mathcal{G}$.

We employ a binary Bayes filter to update the occupancy probability of height interval. The observation model parameters are determined by the height interval’s observation status in the new scan $S_{t}$ and global mapping of pillars:

$$\mathrm{Bayesfilter}(p,P_{S},P_{D})=\frac{P_{S}\cdot p}{P_{S}\cdot p+P_{D}\cdot(1-p)}$$

(5)

	$\displaystyle P_{S}=\begin{cases}\alpha&\text{if interval observed,}\\ 1-\alpha&\text{otherwise.}\end{cases}$		(6)
	$\displaystyle P_{D}=\begin{cases}\beta&\text{if interval observed,}\\ 1-\beta&\text{otherwise.}\end{cases}$

where $\alpha$ and $\beta$ are empirically determined confidence parameters ($\alpha>0.5$, $\beta<0.5$).

When $P_{mn}^{g}$ is successfully matched with $P_{mn}^{t}$ in the new scan, it confirms that $P_{mn}^{g}$ has been observed. Consequently, we update the occupancy probability of empty space within $P_{mn}^{g}$ using Eq. 7, thereby reducing the likelihood of it being occupied:

$$p_{\text{empty}}\leftarrow\mathrm{Bayesfilter}(p_{\text{empty}},1-\alpha,1-\beta)$$

(7)

We extract all height interval endpoints from $P_{mn}^{t}$ and $P_{mn}^{g}$ to construct the set of interval endpoints. Specifically, we extract $b_{i}^{t}$ and $t_{i}^{t}$ from the height intervals in $P_{mn}^{t}$, where $b_{i}^{t}$ and $t_{i}^{t}$ denote the minimum and maximum height values of the $i$-th interval $h_{i}^{t}$ in $\mathcal{H}_{mn}^{t}$, respectively. Similarly, we take $b_{j}^{g}$ and $t_{j}^{g}$ from all height intervals in $P_{mn}^{g}$, where $b_{j}^{g}$ and $t_{j}^{g}$ represent the minimum and maximum height values of the $j$-th interval $h_{g}^{t}$ in $\mathcal{H}_{mn}^{g}$, respectively. The resulting interval endpoint set is defined as:

$$\mathcal{E}=\left\{b_{i}^{t},t_{i}^{t}\mid i=1,...,I\right\}\cup\left\{b_{j}^{g},t_{j}^{g}\mid j=1,...,J\right\}$$

(8)

After sorting these endpoints into an ordered sequence $\mathcal{E}^{sorted}=\{e_{1},e_{2},...,e_{K}\}$ where $e_{k}<e_{k+1}$, we generate new height intervals through pairwise combination:

$$\mathcal{H}^{new}=\left\{h_{k}^{new}=(e_{k},e_{k+1},p_{k}^{new})\mid k=1,...,K-1\right\}$$

(9)

For each new interval $h_{k}^{\text{new}}\in\mathcal{H}^{\text{new}}$, we determine its temporal correspondence by analyzing overlaps with existing height intervals. Specifically, we define $\Omega_{t}$ and $\Omega_{g}$ to represent the observation states of the new height interval within $\mathcal{L}^{t}$ and $\mathcal{G}$, respectively:

		$\displaystyle\Omega_{t}=\{h_{i}^{t}\mid h_{i}^{t}\cap h_{k}^{new}\neq\emptyset\}$		(10)
		$\displaystyle\Omega_{g}=\{h_{j}^{g}\mid h_{j}^{g}\cap h_{k}^{new}\neq\emptyset\}$

When $\Omega_{t}=\Omega_{g}=\emptyset$, this means that the interval is not covered by any historical or current observations, thus we discard this interval. When $\Omega_{t}\neq\emptyset\land\Omega_{g}\neq\emptyset$, it indicates that we have consistent observations, so we apply the Bayesian filter to increase the occupancy probability based on $p_{j}^{g}$:

$$p_{k}^{new}=\mathrm{Bayesfilter}(p_{j}^{g},\alpha,\beta)$$

(11)

When $\Omega_{t}=\emptyset\land\Omega_{g}\neq\emptyset$, the height interval exists in the global pillar but remains unobserved in the current scan $S_{t}$. Consequently, its occupancy probability is reduced based on $p_{j}^{g}$:

$$p_{k}^{new}=\mathrm{Bayesfilter}(p_{j}^{g},1-\alpha,1-\beta)$$

(12)

When $\Omega_{t}\neq\emptyset\land\Omega_{g}=\emptyset$, we update the probability based on $p_{\text{empty}}$ as the static probability for the new height interval:

$$p_{k}^{\text{new}}=\mathrm{Bayesfilter}(p_{\text{empty}},\alpha,\beta)$$

(13)

In our approach, we propose the Low-Height Preservation (LHP) strategy, which leverages the spatial characteristics of the pillar-based structure to prevent incorrect updates of unknown space in a new scan. Specifically, when matching indices between $P_{mn}^{t}$ and $P_{mn}^{g}$, we traverse each height interval in $P_{mn}^{t}$, denoted as $h_{i}^{t}=(b_{i}^{t},t_{i}^{t},p_{i}^{t})$ for $i=1,2,\dots,I$. We then extract the minimum $b^{t}$ and define it as the base height $b_{\text{base}}$ of the current pillar.

This $b_{\text{base}}$ has two possible physical interpretations:

1. 1.

   It may correspond to the ground level observed within the pillar.

2. 2.

   Alternatively, it could be the lowest detected point, constrained by occlusion from objects in the foreground.

The region below $b_{\text{base}}$ is treated as unknown space in the current scan.

To define the lowest base point, we extract the minimum $b^{t}$ value from all height intervals in $P_{mn}^{t}$, using the following formula:

$$b_{\text{base}}=\min\{b_{i}^{t}\mid i=1,2,\dots,I\}$$

(14)

When updating the global height intervals, we apply the Low-Height Preservation (LHP) strategy, maintaining the occupancy probability for the interval as:

$$p_{k}^{\text{new}}=p_{j}^{g}\quad\text{if}\quad\Omega_{t}=\emptyset\land\Omega_{g}\neq\emptyset,\,t_{k}\leq b_{\text{base}}$$

(15)

To prevent excessive probability convergence, we apply a bounded update:

$$p_{h}^{\text{new}}=\text{clip}(p_{h}^{\text{new}})\triangleq\max(0.1,\min(0.9,p_{h}^{\text{new}}))$$

(16)

Finally, we replace the height intervals in the global pillar $P_{mn}^{g}$ with the updated set $\mathcal{H}^{\text{new}}$:

$$P_{mn}^{g}\leftarrow\mathcal{H}^{\text{new}}\circ\text{clip}(p_{h}^{\text{new}})$$

(17)

where $\circ$ denotes the element-wise probability clipping operation.

We classify the compared methods into online and offline approaches, selecting representative works from each category. For offline methods, we include Removert, ERASOR, and the latest SOTA method, DUFOMap. For online methods, we evaluate the classic OctoMap and Dynablox.

For parameter settings, we use the benchmark configurations from [16] for Removert, ERASOR, and Dynablox. OctoMap is set with a voxel size of 0.1 meters. For DUFOMap, we follow the original paper’s parameters, using a voxel size of 0.1 meters, $d_{s}=0.2$ meters, and a sub-voxel localization error of $d_{p}=1$.

For our method, we set the pillar size to 1 meter along both the X and Y axes. For different datasets, we only fine-tuned the observation model parameters $P_{S}$ and $P_{D}$.

All experiments were conducted on a desktop computer equipped with an Intel Core i9-13900KF processor and 64GB of RAM.

To validate the effectiveness of our algorithm, we conducted a precision evaluation and runtime evaluation using the KITTI dataset, with ground truth labels provided by Semantic KITTI. The dataset is captured using the HDL-64E LiDAR and includes pose and point cloud labels. Following the setup in [8] and [16], we selected sequences 00 (frames 4,390 - 4,530), 01 (frames 150 - 250), 02 (frames 860 - 950), and 05 (frames 2,350 - 2,670). During preprocessing, we discarded points labeled as 0 (unknown) and 1 (unlabeled).

To evaluate the adaptability of our algorithm to different sensors, we conducted supplementary experiments on the Argoverse 2 dataset, which utilizes the VLP-32C LiDAR. The dataset setup follows the configuration in [16].

We follow the evaluation approach of [16], using point-level static accuracy (SA) and dynamic accuracy (DA) to directly measure the algorithm’s performance, avoiding the need for ground truth downsampling. Additionally, we employ Associated Accuracy (AA) as a comprehensive metric, defined in Eq. 18, to ensure sensitivity to both static and dynamic points:

$$AA=\sqrt{SA\times DA}$$

(18)

For runtime evaluation, we measure the average runtime per frame and FPS to assess computational efficiency, while the standard deviation of runtime is used to evaluate the algorithm’s stability.

We evaluated our method on the KITTI and AV2 datasets, with results presented in Table I. Removert achieves high static accuracy (SA) but struggles with dynamic accuracy (DA), indicating difficulty in detecting dynamic objects, whereas ERASOR and OctoMap perform well in dynamic object detection but suffer from significant misclassification, leading to lower static accuracy. Dynablox, DUFOMap, and our method HIF maintain a good balance between static and dynamic accuracy, with HIF ranking just below the SOTA method DUFOMap on KITTI while achieving the best performance on KITTI 02, and on AV2, the primary gap between HIF and DUFOMap lies in static accuracy, which we further analyze in the following visualization.

First, we selected the representative KITTI 05 sequence from the KITTI dataset for visualization, as shown in Fig. 5. OctoMap and ERASOR exhibit significant misclassification issues: the former frequently misclassifies ground points, while the latter suffers from sparse ground point retention due to insufficient preservation in the segmentation module, leading to entire pillars being misclassified. Dynablox excels in retaining static points but struggles with detecting dynamic objects. In contrast, our method, HIF, achieves strong performance for both static and dynamic points, with only minor misclassification observed at the top of the canopy.

Then, as shown in Fig. 6, Removert effectively retains static points but struggles with dynamic object removal. Compared to the SOTA method DUFOMap, our approach shows weaker performance in high-rise buildings and canopy areas, with more frequent misclassifications. This is primarily due to the fact that in actual scan data, higher regions appear less frequently within their respective pillars due to occlusions or radar scanning angle limitations, leading to an incorrect reduction in static probabilities.

We evaluated the runtime of all methods on the KITTI sequences, with the results presented in Table II. Our method demonstrated a significant advantage across all sequences, achieving a speedup of several times compared to other approaches. Notably, its runtime remains consistent across different sequences, indicating stable performance across varying environments, whether open or narrow. Furthermore, the per-frame runtime variation is minimal, with low variance, highlighting the method’s robustness. The combination of fast execution and high stability makes our algorithm well-suited for integration into existing SLAM systems, enhancing performance without introducing instability.

To analyze the influence of pillar width on both accuracy and computational efficiency, we conducted a series of experiments (Fig. 7). When the pillar width is small, a high miss rate is observed in dynamic object detection. As the width increases, dynamic accuracy (DA) improves and converges to 95.5%, while static accuracy (SA) remains nearly unchanged.

Furthermore, when the pillar width is small, the algorithm exhibits higher runtime and greater variance. In contrast, at 1 meter and 1.5 meters, the algorithm achieves the best runtime performance.

Additionally, we found that regardless of the pillar width, the algorithm’s peak memory usage remained stable at 801 ± 1MB, indicating high memory efficiency. The maximum memory consumption is primarily determined by the scale of the point cloud in the scene rather than the pillar width.

Since Height Interval Construction and Interval Probability Filtering form the core components of our approach, we conduct ablation experiments exclusively on the Low-Height Preservation (LHP) module. To evaluate its effectiveness, we test this module across all data sequences.

Results presented in Table III, showing that while LHP has a negligible impact on dynamic accuracy (DA), it significantly improves static accuracy (SA); to further investigate this effect, we visualized the results for the KITTI 02 sequence in Fig. 8, where the complex scene—featuring parked vehicles, bushes, and occluded objects—leads to relatively low static accuracy for most algorithms, as shown in Table I, whereas our LHP strategy effectively identifies these occlusions and preserves static points in occluded areas.