Fast online high-order time lacunarity for characterizing active sonar echographs of harbor environment Open Access

A challenge for the application of active sonar systems is to detect low-strength threats such as divers and unmanned underwater vehicles in harbor environments. The sonar echo energy from different ranges and bearings is normally displayed as an echograph in brightness levels (Yang et al., 2012). However, the performance of active sonar in shallow water is often limited by false alarms arising from multipath-induced high-level clutter, which tends to be confused with real targets (Fialkowski and Gauss, 2010; Kessel and Hollett, 2006). Thus, this Letter focuses on effective target detection in active sonar echographs in practice.

Since the matched-filter envelope in shallow water has been observed to depart from the traditional Rayleigh assumption, a variety of statistical methods have been investigated for characterizing clutter probability distributions functions (PDFs). Some models attempt to fit the observed amplitude distributions, including log-normal (Chotiros et al., 1985) and Weibull (Gensane, 1989) distributions. Others, such as K-distribution (Abraham and Lyons, 2002) with a more physical basis, relate clutter PDFs to sonar systems and environmental mechanisms. However, with the drawbacks of the difficulty of parameter estimation and lack of generality, it is not robust enough under complex scenes.

With growing interest in automatic target hunting on sonar imagery, some unsupervised texture features worked well in distinguishing background areas from potential targets. The feature of anisotropy is helpful for detecting dense linear objects, and the feature of complexity reflects the amount of background clutter in sonar images (Fakiris et al., 2013; Picard et al., 2018). More recently, a more robust feature of lacunarity, a measure of pixel-intensity variation, has been effectively applied on hunting static mines in high-resolution sonar imagery (Myers, 2002; Williams, 2015). On this basis, Zhao et al. (2020) extends a new high-order time lacunarity (HOT-Lac), a feature of highlighting non-stationarity, to detect moving targets from high-level clutters of active sonar echographs. However, during the surveillance of moving targets, it is necessary to ensure real-time and reliable detections for threat defense. If the processing time of the target detection and tracking is too long, namely, longer than the transmitting cycle of the signal (generally no more than 1 s), then the threats might have escaped from the attack center. Therefore, it raises a higher demand of computation efficiency of HOT-Lac for fast detection of underwater moving targets.

An integral image (Viola and Jones, 2004) is a special data structure that makes Haar-like features computed very rapidly at any scale or location. Using the integral image any rectangular sum can be computed in most four array references. Since this representation was effective for real-time detections in the computer-vision field, it has been extended to sonar image classification (Daniell et al., 2012) and underwater target detection (Williams and Groen, 2011). Moreover, Williams (2015) has successfully applied it to the fast computation of lacunarity on high-resolution synthetic-aperture-sonar imagery.

Based on the idea of integral image, this Letter proposes the echo-intensity integral sequence of filtering sonar echographs along the time domain recursion. With this representation, online calculation methods of HOT-Lac are performed quickly in a constant number of operations. Furthermore, the feature map can be obtained frame-by-frame by using the current sonar data to correct the previous result. The contribution of this work lies in fast online methods of computing HOT-Lac for its promising real-time application of harbor surveillance.

The HOT-Lac method is fully described in Zhao et al. (2020). Given a sequence of active sonar echographs ${Pxy,i}i∈N$ in a set N, the HOT-Lac is defined as

where k$(≥2)$ represents the order of HOT-Lac, n is the number of frames of calculated echographs (i.e., the length of the time window), and $μ=(1/n)∑i∈NPxy,i$ is the corresponding mean value. Sliding a time window over the sonar echograph sequence, a feature map of the HOT-Lac is obtained at any specified frame interval. However, when the time window sliding from frames P₁-P_n to P₂-$Pn+1$⁠, according to Eq. (1), most frames need to be stored in advance and reused to update of HOT-Lac. However, it not only wastes storage space but also affects the processing speed because of the loop calculation. Considering on these problems, we attempt to update the HOT-Lac online using the current frame of sonar echograph efficiently.

The integral image (Crow, 1984) as an intermediate representation allows very fast summation of image values at any scale or location in constant time. Since this representation was popularized by Viola and Jones (2004) for real-time detections in the computer-vision field, more advantage of it has been taken for the purpose of rapidly computing approximations. As shown in Fig. 1(a), an integral image I(x, y) at location $x,y$ corresponds to the sum of the pixels above and to the left of (x, y), in the original image, f(x, y), inclusive of

Inspired by this, we extend an echo-intensity integral sequence by filtering sonar echographs along time domain recursion. As shown in Fig. 1(b), for the sequence ${Pxy,i}i∈N$⁠, we define that the echo-intensity integral sequence I(t) is constructed by the corresponding sum of each pixel position ahead of the tth frame, namely,

With the recursive relation as $I(t)=I(t−1)+Pxy,t$⁠, the echo-intensity integral sequence is generated in a 3D data structure along the time dimension. Moreover, with only two array accesses of it, the sum of sonar echographs in a time window n about a given time is easily expressed as $S(t)=I(t)−I(t−n)$⁠, which is important for online HOT-Lac described below.

In order to avoid the waste of computing resources caused by repeated use of sonar echographs for loop summation in the previous HOT-Lac calculation (Zhao et al., 2020), we propose an online method just using the current active sonar echograph to complete the intermediate calculation of HOT-Lac efficiently. The method of online HOT-Lac1 focuses on making tractable computations through the intermediate representation of echo-intensity integral sequences. Through sliding a time window with one frame step across the sequence of active sonar echographs, the HOT-Lac feature map can be obtained continuously with the definition that the updated feature map corresponds to the latest input data.

First, according to the binomial theorem to expand the powers of a binomial, i.e., $(a+b)k=Ck0ak+Ck1ak−1b+⋯+Ckmak−mbm+⋯+Ckkbk$ (⁠$Ckm=k!/m!(k−m)!$ is the combinatorial number), the HOT-Lac of Eq. (1) corresponding to the active sonar echograph t can be expressed as

where $t0=t−n+1$⁠. Moreover, based on the echo-intensity integral sequence, let summations of $Pxy,i$ in the numerator of Eq. (4) be described as $Sk(t)=∑i=t0tPxy,ik=Ik(t)−Ik(t−n)$ (the number of created integral sequence is equal to the HOT-Lac's order k); then, we have

Further, a matrix form of computing online HOT-Lac1 is conducted as

where W, U, S are all column vectors of length k + 1. Each element in vector U and S is a 2D-matrix whose size is same as the original sonar echograph, and 1 is a matrix of ones. As a weight coefficient, each element in vector W is a scalar, specifically

It is important to point out that A is a matrix obtained by sums in rows from vectors W, U, and S by the operation of Hadamard product which multiplies corresponding elements of each vector. The same way is applied for the division between matrix A and $S12(t)$⁠. After putting $μ=S1(t)/n$ into Eq. (6), the final online HOT-Lac1 is a function of $S(t)={S1(t), S2(t), … ,Sk(t)}$ constructed as

The function illustrates that when the precalculated S(t) is obtained by new echo-intensity integral sequences, the HOT-Lac feature map can be updated synchronously in a constant number of operations.

Although the online HOT-Lac1 implements continuous update by a sliding time window, calculations start only when the number of accumulated sonar echographs reaches the window size. For this reason, we consider a more extreme case for updating the feature map in a stream (i.e., the feature map is generated with the arrival of each new sonar echograph), and then propose the method of online HOT-Lac2. Even better, the online HOT-Lac2 could directly update the feature map using the latest frame of sonar echograph.

Specifically, we directly make an approximate calculation for Eq. (1), namely,

It is worth noting that unlike HOT-Lac1, where the mean value μ is computed within a window of defined length, the mean value μ_i of HOT-Lac2 is calculated in the time window with varying length. In Eq. (8), the mean value μ_i is just subtracted by the current sonar echograph $Pxy,i$⁠. With the help of an echo-intensity integral sequence $I1(t)$⁠, the online HOT-Lac2 can be implemented in a simple form as

To make it clear, $ϕt=∑i=1t−1(Pxy,i−I1(i)/i)k+(Pxy,t−I1(t)/t)k$ is an iterative process. The advantage of this change is that it avoids subtracting the new mean μ_t from previous t – 1 frames of sonar echographs, the current HOT-Lac2 is only related to the latest input sonar echograph $Pxy,t$ [i.e., $I1(t)$ is derived from $Pxy,t$] which can be discarded after a look. Therefore, the online HOT-Lac2 not only improves computational efficiency but also reduces data storage space.

The proposed two online methods for HOT-Lac calculation are tested by using four data sets of measured sequences of active sonar echographs with four types of moving targets from Zhao et al. (2020), and Fig. 1(c) shows one example from the last frame of each data set.

By repeating the feature extraction experiment many times with different parameters of HOT-Lac, we can obtain a robust estimate of the performance for our proposed methods. Specifically, the order k of HOT-Lac is chosen to be 2, 3, and 4. And for the fixed time window, the size n is selected as 10, 20, and 30 frames (considering the properties of online methods on time-effectiveness application, the window size is not set too large.). In addition, the moving step size of all methods is one frame.

All algorithms were executed on a desktop computer with a 2.3-Ghz Intel Core i5 CPU and 8 GB of RAM. The average computation time per HOT-Lac map is shown in Table 1. As can be seen from the table, both of the two proposed online methods are several times faster than the original method HOT-Lac, and with the increase of order k or window size n, such advantage is more obvious. Moreover, since the method of online HOT-Lac2 does not have a fixed-length time window, there is only one result corresponding to each order of it.

Then, in Fig. 2, we show the real computation time of each feature map in order of appearance. It can be seen that like the method HOT-Lac, the online HOT-Lac1 consumes much more time at the beginning of the process. Because it has to consume a period of time to accumulate enough frames of sonar echographs (we show the result of n = 20) for the start of calculation. Then, stating with the second feature map, the computation time of the two proposed online methods is similar, and is about an order of magnitude faster than the method HOT-Lac with the increasing computational complexity in higher order k.

In general, there is often a trade-off between accuracy and speed of the algorithms, therefore, we also evaluate the detection performance of these two online HOT-Lac methods with the original HOT-Lac by the area under the receiver operating characteristic curves (AUC), which is the same approach used in Zhao et al. (2020). Figure 3 shows the AUC values corresponding to each frame of feature map, where the order k = 3 and the fixed window size n = 20, and it accompanies the result of the last feature map of each method for a more intuitive comparison. It can be seen that the detection performance of the proposed online HOT-Lac1 is basically the same as the method HOT-Lac. Moreover, compared to the original active sonar echograph in Fig. 1(c), the feature map of HOT-Lac1 is pretty much the same as HOT-Lac's with effective clutter suppression.

For the method of online HOT-Lac2, the AUC value is incremented from zero in the first ten frames or so and then basically remains stable, which is similar to the HOT-Lac method. It should be noted that, as shown in Figs. 3(a)–3(b) and Fig. 3(d), when the frame number gradually increases, the AUC value of online HOT-Lac2 begins to decrease slightly. Since it needs to use the sonar echographs from the first frame to the last for HOT-Lac calculation, the target shown on the feature map is the cumulative result of the entire process. If the observation time is too long, then targets that appeared once earlier in some pixels have disappeared for a while, which reduces the intensity of these target pixels on the feature map. This is obvious in the example feature maps of online HOT-Lac2. Moreover, the high-level interference accumulated by dynamic clutters will be retained on more frames of feature maps, which is one of the reasons for the degradation of detection performance for online HOT-Lac2.

Correspondingly, we also compare the average AUC values of all methods under different HOT-Lac parameters in Table 2. For the two algorithms we proposed, the overall detection performance was not influenced significantly while the speed increased, especially for the online HOT-Lac1. Actually, the calculation of HOT-Lac2 is conducted under the variable length of time window. When the window length of HOT-Lac2 is less than the fixed time window of HOT-Lac, the corresponding AUC value is much smaller. [It is also pointed out in Zhao et al. (2020) that when the window length is less than 10 frames, the AUC value is extremely small, which results in the smaller average AUC value.]

For the case of order $k≥5$⁠, which is not discussed in this paper, even though the result of AUC is different under different orders, the relationship of the performance for these three algorithms is the same under one specified order (i.e., the average AUC value of the HOT-Lac2 is lower than that of HOT-Lac and HOT-Lac1). Moreover, Zhao et al. (2020) pointed out that an optimal order for HOT-Lac exists in the case of detection performance, so there is no need to increase the order k indefinitely.

In this Letter, based on the technique of the integral image, we extended an echo-intensity integral sequence to filter images along time dimension, resulting in two kinds of fast and accurate online methods of HOT-Lac update. HOT-Lac computations can be made tractable in a constant number of operations by use of intermediate representations for reducing the high computational complexity of loop summations. Using the active sonar echograph data sets collected in real-world harbor environments in the South China Sea, performances were evaluated for different HOT-Lac parameters. The results showed that the computation time of the proposed methods decreased significantly with several times even an order of magnitude faster than before. In addition, the detection performance represented by AUC values still maintains high rates. These are good indicators of the potential of online HOT-Lac methods for real-time detection and tracking in harbor environments.

This work was supported in part by the National Key R&D Program of China under Grant No. 2016YFC1400200, in part by the National Natural Science Foundation of China under Grant No. 61671388, and in part by the National Engineering Laboratory for Integrated Aero-Space-Ground-Ocean Big Data Application Technology.