An Analysis of Parallelized Motion Masking Using Dual-Mode Single Gaussian Models

TL;DR

This paper explores parallelized motion masking using dual-mode single Gaussian models, implementing the method with TBB and CUDA, and analyzing speedup improvements for motion detection in video.

Contribution

It introduces a parallel implementation of a dual-mode Gaussian motion masking technique using TBB and CUDA, with a comparative analysis of performance improvements.

Findings

Parallel implementations significantly improve processing speed.

CUDA outperforms TBB in speedup metrics.

Theoretical analysis aligns with experimental results.

Abstract

Motion detection in video is important for a number of applications and fields. In video surveillance, motion detection is an essential accompaniment to activity recognition for early warning systems. Robotics also has much to gain from motion detection and segmentation, particularly in high speed motion tracking for tactile systems. There are a myriad of techniques for detecting and masking motion in an image. Successful systems have used Gaussian Models to discern background from foreground in an image (motion from static imagery). However, particularly in the case of a moving camera or frame of reference, it is necessary to compensate for the motion of the camera when attempting to discern objects moving in the foreground. For example, it is possible to estimate motion of the camera through optical flow methods or temporal differencing and then compensate for this motion in a background subtraction model. We selection a method by Yi et al. using Dual-Mode Single Gaussian Models which does just this. We implement the technique in Intel's Thread Building Blocks (TBB) and NVIDIA's CUDA libraries. We then compare parallelization improvements with a theoretical analysis of speedups based on the characteristics of our selected model and attributes of both TBB and CUDA. We make our implementation available to the public.

Figures

Equations

$$\frac{1}{\sigma \sqrt{2\pi }}{e}^{-{\left(x-\mu \right)}^2/2{\sigma }^2}$$

$$|G_i|=1$$

$${\mu }_i^{(t)}=\frac{{\tilde{\alpha }}_i^{(t-1)}}{{\tilde{\alpha }}_i^{(t-1)}+1}{\tilde{\mu }}_i^{(t-1)}+\frac{1}{{\tilde{\alpha }}_i^{(t-1)}+1}{M}_i^t$$

$$M_i^{(t)}=\frac{1}{|G_i|}\sum _{j\in G_i}{I}_j^{(t)}$$

$${\sigma }_i^{(t)}=\frac{{\tilde{\alpha }}_i^{(t-1)}}{{\tilde{\alpha }}_i^{(t-1)}+1}{\tilde{\sigma }}_i^{(t-1)}+\frac{1}{{\tilde{\alpha }}_i^{(t-1)}+1}{V}_i^t$$

$$V_i^{(t)}=\underset{j\in G_i}{\max }{\left({\mu }_i^{(t)}-I_j^{(t)}\right)}^2$$

$${\alpha }_i^{(t)}={\tilde{\alpha }}_i^{(t-1)}+1$$

$${\left(M_i^t-{\mu }_{A,i}^{(t)}\right)}^2<{\theta }_s{\sigma }_{A,i}^{(t)}$$

$${\left(M_i^t-{\mu }_{C,i}^{(t)}\right)}^2<{\theta }_s{\sigma }_{C,i}^{(t)}$$

$${\alpha }_{A,i}^{(t)}<{\alpha }_{C,i}^{(t)}$$

Full text

1. 1 Introduction

2. 2 Background

3. 2.1 Preprocessing

4. 2.2 Single Gaussian Model

5. 2.3 Dual-mode Single Gaussian Model

6. 2.4 Motion Compensation

7. 3 Methodology

8. 3.1 Serial Implementation

9. 3.1.1 Preprocessing

10. 3.1.2 Dual-Mode SGM

11. 3.1.3 Motion Compensation

12. 3.2 Thread Building Blocks Implementation

13. 3.3 CUDA Implementation

14. 3.3.1 Preprocessing

15. 3.4 Dual-mode SGM

16. 4 Results and Analysis

17. 4.1 Qualitative Algorithm Analysis

18. 4.1.1 Effects of Preprocessing

19. 4.1.2 Comparison to OpenCV MOG

20. 4.1.3 Motion Compensation Results

21. 4.2 Thread Building Blocks

22. 4.2.1 Test Specifications

23. 4.2.2 Results

24. 4.3 CUDA

25. 4.3.1 Machine Specifications

26. 4.3.2 Preprocessing

27. 4.3.3 Dual-mode SGM and Thread Block Optimization

28. 4.3.4 GPU Comparison

29. 5 Conclusion

30. 5.1 Issues and Usability of TBB

31. 5.2 Issues and Usability of CUDA

32. 5.3 Future Work

33. A Device Specifications: GeForce GT 650M

34. B Device Specifications: GeForce GTX 570

35. C Preprocessing Gaussian Blur Code Snippet

36. D Preprocessing Median Blur Code Snippet

37. E Dual-mode SGM Inner Update Code Snippet

38. F Motion Compensation Code Snippets

1 Introduction

Motion detection in video is important for a number of applications and fields. In video surveillance, motion detection is an essential accompaniment to activity recognition for early warning systems [1, 2, 3]. Robotics also has much to gain from motion detection and segmentation, particularly in high speed motion tracking for tactile systems [4]. Others can clearly benefit from such technologies as well. However, motion detection is a hard problem, particularly with a moving camera reference frame, frequent changes in lighting or dynamic backgrounds. As such, many attempts at tackling this problem have been made. Modern research into motion detection and segmentation systems generally falls into 3 categories [5]: background subtraction, temporal differencing, and optical flow.

Background subtraction involves developing a model in a series of images which differentiates the foreground (i.e. motion) from its background. A background model is typically learned on a per-pixel basis. A large portion of background modeling techniques involve keeping a Gaussian Model (GM), or number of candidate GMs, of the background for each pixel. Wren et al. [6] first used Single Gaussian Models (SGM) (with one GM per pixel) for motion detection. The Gaussian model is updated incrementally per some learning rate depending on whether or not a pixel intensity falls within the model’s variance. While this technique provided acceptable results for certain cases, large objects moving through the image or changes in lighting proves to contaminate the model significantly. Further expansion on this involves keeping several candidate Gaussian models indicating possible backgrounds, and a probabilistic mixture is used to differentiate the foreground from the background [7]. This typically yields better results in most cases, but can still run into the same problems as SGMs.

In temporal differencing, the absolute difference of pixel intensities across two or more frames is averaged and a threshold value is used to determine when motion occurs. Component analysis methods are used to cluster differences into regions of motion [8, 3].

Optical flow methods are widely varied, but there are several methods used often as a basis for many complex systems. Two often discussed algorithms are Lucas-Kanade [9] and Horn-Schunk [10]. Lucas-Kanade find the vector transformation of pixels by using least-squares methods to determine the flow for a neighborhood of pixels, while Horn-Schunk used the smoothness constraint in the aperture problem to find optical flow in an image.

Many newer motion detection systems use an ensemble of these methods to accomplish the complex task of motion detection. Particularly in the case of a moving camera or frame of reference, it is necessary to compensate for the motion of the camera when attempting to discern objects moving in the foreground. For example, it is possible to estimate motion of the camera through optical flow methods or temporal differencing and then compensate for this motion in a background subtraction model [11, 12].

We select a method by Yi et al. [11] to detect motion in a series of images (or video) using one of these ensemble methods (though the focus is on background subtraction methods). We use OpenCV for some basic tools (like image parsing and display) as well as some of the motion compensation and preprocessing. We parallelize the motion detection algorithm using Intel’s Thread Building Blocks (TBB) and NVIDIA’s Compute Unified Device Architecture (CUDA) libraries.

Intel’s TBB is a library which provides abstractions for interaction with operating system provided threading libraries (i.e. POSIX) [13, 14]. Rather than dealing with thread synchronization and handling manually, TBB provides an abstraction which allows for the creation of tasks which are allocated to different machine cores by the library run-time engine. TBB also handles CPU-cache transactions to make processing more efficient. We choose TBB because of the abstractions available to simplify parallelization as well as the optimizations of task scheduling and cache management which the library provides. We also considered OpenMP [15], but ran into compilation issues relating to the Mac OS X gcc compiler’s interaction with OpenCV. OpenMP is not supported by older versions of the gcc tool on Mac. However, newer versions of the tool are not compatible with OpenCV on Mac. We then reverted to using TBB as it seemed like a sensible and powerful alternative.

NVIDIA’s CUDA is a platform which gives developers access to low level GPU instruction sets [16]. References to CUDA often refer to the NVIDIA CUDA C/C++ language and libraries which provide high level control over the GPU. It is this that we use for our project. While this can vary between graphics cards, generally the paradigm for CUDA GPUs is as such (based on the Kepler architecture [17], but applies similarly to others like Fermi [18]). The GPU contains several multiprocessors (MPs) each with many CUDA cores (processing units). Processing functions are divided into blocks of threads. Each block is assigned to one MP. Once a block is assigned, another MP cannot access that block or its shared memory. The block of threads is further assigned to CUDA cores by a warp scheduler. A warp is a group of threads which will be scheduled onto CUDA cores in parallel, typically in blocks of 32, but this depends on the warp size of the card. The threads can either fetch from global or shared memory, but shared memory is typically shared only within a block on an MP.

2 Background

2.1 Preprocessing

To improve the accuracy of the motion detection, as in Yi et al., we first pre-process the incoming frame. Two filters are used for this task: a Gaussian filter and a median filter. A Gaussian filter creates a matrix of size $N\times N$ with a Gaussian distribution of weights across the matrix. The Gaussian distribution function for a one dimensional vector can be seen in Equation 1 where $\sigma$ is the standard deviation, $\mu$ is the mean value and $x$ is the distance parameter desired.

[TABLE]

The image can then be convolved with this normal distribution of weights to set a pixel’s value to the summed normal distribution of its neighbors for some given size and sigma value.

In a median filter, a pixel is set to the median value of its neighbors. This preprocessing is done to smooth the edges of motion detection and ensure that fewer per-pixel lighting anomalies are detected.

2.2 Single Gaussian Model

In SGM background subtraction methods, the mean and variance of every pixel in a series of images is kept. They are updated based on a specified learning rate to develop a model for the background of an image. In their method, Yi et al. use a variable learning rate based on the age of the model to place more emphasis on earlier frames. However, the age is capped at 30 to prevent excessive emphasis and ghosting effects222Ghosting is referred to contamination of the background model such that motion is always observed there because a previously stationary object has moved from the background, leaving a ghosting effect.. On each iteration which affects the background model, the mean, variance and age are updated as follows. We take the equations here directly from Yi et al.’s paper. However, it is important to note that while Yi et al. use an average of pixel intensities in a grid of size size $m\times m$, we calculate these on a per pixel basis, and hence, for each grid the size of our Grid is 1:

[TABLE]

The mean of the SGM is updated according to:

[TABLE]

Here, the mean of the background model $\mu_{i}^{(t)}$ is updated based on the previous mean and the current pixel intensity. While Yi et al. use $M_{i}^{t}$ to be defined as follows:

[TABLE]

We use a per-pixel model, and thus $M_{i}^{t}$ is equivalent to $I_{j}^{(t)}$. Further, the variance is updated according to:

[TABLE]

Here, the variance $\sigma_{i}^{(t)}$ is updated based on the age, the previous variance, and the absolute difference of the mean and the current pixel intensities. Again, though Yi et al., define:

[TABLE]

We use a per-pixel model, and as such this just becomes $\left(\mu_{i}^{(t)}-I_{j}^{(t)}\right)^{2}$. Lastly, the age of the model is updated according to:

[TABLE]

In all of the above, $\tilde{\alpha},\tilde{\mu},\tilde{\sigma}$ represent the previous values after motion compensation, as will be discussed later. While an SGM may be acceptable for a static background and relatively simple background subtraction tasks, large objects moving through the image sequences can contaminate the background and changes in lighting will result in large errors. To overcome this, Yi et al. propose a Dual-mode SGM model.

2.3 Dual-mode Single Gaussian Model

In the Dual-mode SGM model, an absolute background model and a candidate background model are updated in concert. The absolute model represents the current model used for differentiating foreground from background. The candidate model represents a possible new background after, say, a shift in lighting or camera position. These models are updated as follows. Per pixel, if the squared absolute difference of the absolute model’s mean and the pixel intensity is less than some weighted value of the variance, then the absolute model is updated as per the aforementioned SGM rules. This can be formulated as:

[TABLE]

If the absolute model is not updated, then the candidate model will be updated if it falls within the similarly defined threshold as follows:

[TABLE]

If neither model is updated, the candidate model is reset with the current mean set to the current pixel intensity, the variance set to some initial value (255 for the purposes of our implementation), and the age to 1.

If the age of the candidate model exceeds the age of the absolute model as in:

[TABLE]

the absolute model is set to the candidate model and the candidate model is reset.

2.4 Motion Compensation

The motion compensation described in Yi et al. is as follows. Kanade-Lucas-Tomasi’s (KLT) feature tracking algorithm333KLT is a method of tracking the movement of pixels in an image sequence. In it, features (pixel sequences) are identified which are unique in an image. Methods based on Lucas-Kanade’s optical flow are used to track these pixels through a series of images. Other methods similar to this include SURF [19]. is used to find feature points and the vector shifts of the pixels [9, 20]. RANSAC [21] is then used to create a homography matrix from these vectors. Because in their method, Yi et al. use grid blocks rather than single pixels, they compute the mixture between grid blocks based on this homography matrix. However, because we compute SGMs on a per-pixel basis, we can simply use the homography matrix to map the SGM values to the appropriate pixels based on the overall motion of the image.

3 Methodology

3.1 Serial Implementation

3.1.1 Preprocessing

For Gaussian and median blurring, we use two methods of implementation. First, we use the current OpenCV implementation of both. Second, we implement our own unoptimized methods for both. In the serial Gaussian blur, for each pixel, we multiply a 2D Gaussian kernel (created by the Gaussian equation) by the surrounding pixels. The sum of this is then placed in the current pixel’s position in the updated image. For median blurring, we take the surrounding 8 pixels (though this can be changed, we generally stayed with 8), as well as the current pixel, and use a simple sorting algorithm to order their intensities. The median value is then used as the new value of the pixel. A snippet of the core serial code of these methods can be found in Appendices D and C.

3.1.2 Dual-Mode SGM

As per the aforementioned equations described by Yi et al., we implement the Dual-mode SGM model as a core loop, keeping the age, variance, and mean of the two models as matrices. The core code iterated for each pixel can be seen in Appendix E.

3.1.3 Motion Compensation

We implement motion compensation using OpenCV’s implementation of KLT and RANSAC. We first use the findGoodFeatures function to generate a list of possible points. We then use KLT (by the method calcOpticalFlowPyrLK) to find the movement vectors associated with these points. Finally we feed these vectors into OpenCV’s RANSAC (through the findHomography method). Lastly, we transform the current frame to match the previous frame’s keypoints via this homography matrix to attempt to reduce error from overall camera motion. Snippets of this code can be found in Appendix F. It is important to note that we did not achieve significant improvements from this method, which we will discuss in the results section. Because of the poor results and due to the efficiency of the existing OpenCV implementation and the complexity of the algorithms in comparison with the amount of time available for the project, we determined it to be out of the scope of the project to parallelize motion compensation as well as all other aspects.

3.2 Thread Building Blocks Implementation

Generally, for both the preprocessing parallelization and the Dual-mode SGM parallelization, we separate the inner process functions from the loops iterating over the pixels and formulate them as tasks with the overloaded TBB functions. We specify the number of tasks desired using the OpenCV TBB extension (cv::parallel_for method). The image matrix is decomposed by rows based on the number of tasks allotted. TBB then parallelizes the tasks among multiple threads and processors.

3.3 CUDA Implementation

With CUDA, we used a general scheme of decomposition for the current image frame. We choose a number of threads per block to form square blocks ($N\times N$) and then we create a grid of blocks to match the size of the image. We provide CUDA with a simple kernel such that each thread will process a given pixel. We generally place the entire image into the global memory of the GPU, since past experimentation [22] seems to indicate that implementing a shared mechanism for the size of images we were using would have more overhead than a global memory implementation. The pixel that a thread worked on is determined according to standard CUDA conventions as:

const size_t r = blockIdx.x * blockDim.x + threadIdx.x;  const size_t c = blockIdx.y * blockDim.y + threadIdx.y;  const size_t index = r * numCols + c;

Here, r is the row of the pixel, c is the column of the pixel, and index becomes the position of the current pixel in an array based on the current Thread ID, Block Dimension, and Block ID (provided by CUDA libraries).

3.3.1 Preprocessing

To parallelize the Gaussian blur, we first attempted to simply write a CUDA kernel which would take as input the image and Gaussian kernel. The kernel would determine which pixel the current thread was assigned to through previously described means. Then the surrounding pixels would be convolved with this kernel. This however, is a naive implementation. In this scenario $O(nmk^{2})$ operations are needed, where $n\times m$ is the size of the image, and $k\times k$ is the size of the kernel. To significantly reduce this cost, we can use a separable convolution. In this, a 1-Dimensional Gaussian kernel is used to convolve first the rows and then the columns [23]. This results in an equivalent output while requiring only $O(nm2k)$ operations. For the 2-Dimensional implementation of the Gaussian filter, we generally use the same code as the serial implementation in the CUDA kernel (obviously without looping over the pixels). For the separable convolution, we use similar code, but rather than running through a loop over a 2D Gaussian kernel for convolution we first run a CUDA kernel over the rows with a 1D Gaussian kernel, then we run it over the columns using the same code, with a flag raised. The core inner code remains fairly similar to the serial functions and can be found in the cuda_kernels.cu file

For median blurring, we simply take the surrounding pixels, order them through a simple swap sorting mechanism and then choose the median value just as in the serial implementation. However, we divide this per kernel such that each thread calculates the values for only one pixel. Generally, the core of this code remains the same as the serial implementation.

3.4 Dual-mode SGM

To parallelize the Dual-mode SGM, we take the same inner code as the serial version and similarly remove the two outer loops which iterate over the pixels in favor of the CUDA kernel thread model described above.

4 Results and Analysis

For both qualitative and quantitative analysis, we use the University of Sherbrooke CDNet database [24]. The database provides image sequences and videos for use in motion detection in different scenarios (including with motion of the camera).

4.1 Qualitative Algorithm Analysis

We first analyze the qualitative results of the algorithm. We note that all parallelized versions were verified against the serial version to ensure that no errors were introduced into the algorithm by parallelization. As such, this qualitative comparison is not necessary.

4.1.1 Effects of Preprocessing

We determined that preprocessing has is an important aspect in filtering out noise during motion detection. This is clearly visible in Figures 1,2,3, and 4. Without preprocessing, the motion detection segmentation mask becomes grainy as more small anomalies in lighting are picked up on a per pixel basis. After Gaussian and median filtering, the anomalies disappear as the image smooths out.

4.1.2 Comparison to OpenCV MOG

Next, we wanted to compare our implementation of the Yi et al. algorithm against the Mixture of Gaussian model in OpenCV. MoG can potentially yield better results as it has a number of background models which can be used to develop a more robust representation of the background. In contrast to our expected results, we found that the MoG implementation of OpenCV did not produce readily identifiable objects in the image as seen in Figure 5. Most of the man is considered to be part of the background, while the man is more identifiable in Figure 6 for our implementation. However, it is worth noting that any place where the color of the man’s clothes blends in with his background, the algorithm cannot identify him as part of the foreground (as in his pants blending into the couch).

4.1.3 Motion Compensation Results

While we implemented motion compensation using OpenCV’s LKT, RANSAC, and transform functions, our results were not ideal and as such we did not focus on it for parallelization or other aspects. In a dataset with significant jitter as seen in Figures 8, 9, and 10, we saw minimal benefit to having motion compensation versus not. As can be seen in Figure 11, the MoG algorithm has a much cleaner image in this case. However, this is to be expected as it is both better suited to adaptable backgrounds and is much less sensitive than our implementation.

4.2 Thread Building Blocks

Originally, OpenMP was pursued for easy parallelization. However, on Mac OS X, compilation difficulties when building with OpenCV and OpenMP simultaneously were insurmountable. Since Thread Building Blocks (TBB) is readily supported on Mac OS X systems in combination with OpenCV and the compatible version of gcc, we use this.

4.2.1 Test Specifications

Combinations of using the OpenCV preprocessing filter versus an implemented version as well as including motion compensation in the computation yields four total variations: OpenCV preprocessing with Motion Compensation, OpenCV preprocessing without Motion Compensation, our preprocessing with Motion Compensation, and our preprocessing without Motion Compensation. Furthermore, serial and parallelized implementations of each variation were implemented.

An example set, showing a man enter the frame, sit on a sofa, and finally leave, was used. We use the first 500 frames of the set which has a resolution of 320x240. An original raw frame can be seen in Figure 7. The test machine consisted of an Intel Core i5 at 2.4Ghz with two cores and four hardware threads.

4.2.2 Results

Figure 12 displays the speedup of each the four parallel variations as compared to their respective serial versions with a varying number of TBB tasks. Figures 13, 14, 15 and 16 reveals the time breakdown spent in each of the major parts of the complete algorithm for all four variations. It is important to note that while we can specify tasks to be parallelized TBB by default creates threads optimized for the number of processor cores and hardware threads. TBB then schedules these tasks onto the cores and hardware threads [25]. As such, we speak of tasks, not threads.

Across all the variations, there are a number of general trends. The optimal number of TBB tasks for running the algorithm plateaus between four and sixteen. Since the test machine has four hardware threads, this is a consistent result. TBB can schedule all the tasks successfully among the hardware threads. Comparing the runtime of the serial version to a single TBB task reveals the parallelization overhead associated with TBB since the serial version is typically slight faster. Furthermore, specifying TBB to use a number of tasks greater than 256 results in times comparable to only using one task. An unreasonably high number of tasks results in the framework not being able to schedule the task among multiple processors or hardware threads, and hence the result is the same as having only one thread. At this point the tasks become to small so the first core is complete with a task before the library can schedule the next one on a different processor so it schedules it on the same thread on the same core [25]. The addition of motion compensation into an algorithm delivers a near constant time addition to the serial portion of the runtime with little variation.

With regards to using the OpenCV Gaussian and median filters for preprocessing as compared to the serial and parallelized filters implemented in this project, the OpenCV implementations are consistently much faster. Despite any TBB parallelization used to speed up filtering, the OpenCV image processing library is optimized based on the architecture of the processor used 444The optimizations within OpenCV were revealed by examining the source code found on GitHub. Many preprocessing flags for specific processor architectures are used, indicating hardware specific optimizations. In particular: https://github.com/Itseez/opencv/blob/f50f249f804b35ad6f39459ec6df9a66fb0825f8/modules/imgproc/src/smooth.cpp. Ultimately, any attempt to re-implement and parallelize parts of OpenCV would most likely be significantly slower.

4.3 CUDA

4.3.1 Machine Specifications

The main graphics card we use for testing is an NVIDIA GeForce GT 650M graphics card in conjunction with a 2.6 GHz Intel Core i7 processor. Though full specifications for the card from CUDA’s deviceQuery program can be seen in Appendix A, the main points to note are the following. It is a Kepler architecture card with two multiprocessors with 192 CUDA cores per multiprocessor. It has a warp size of 32 and a quad-warp scheduler (it schedules 4 warps at once for a total of 128 threads). Additionally it can handle a maximum of 1024 threads per block.

4.3.2 Preprocessing

We first compared our separable Gaussian filter implementation to the 2-Dimensional implementation. We ran this test on the NVIDIA GeForce GT 650M graphics card. We excluded motion compensation, but retained median blurring and the CUDA-parallelized Dual-mode SGM background detection for all versions that we ran. We found as we previously did with TBB that OpenCV Gaussian and Median filters are significantly faster overall due the heavy optimizations done within the code. Additionally, though the CUDA parallelized version of separable Gaussian filtering comes closer to comparing to the OpenCV method, but still cannot compete with its optimizations. As such, for further tests we retain the OpenCV pre-processing methods.