{'fraction_non_alphanumeric': 0.029728407144604844, 'fraction_numerical': 0.015170051382432102, 'mean_word_length': 4.591313269493844, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Optimization and Performance Evaluation of Stereo-Matching Software on Many-core Processors

The high-performance Geospatial Information System (GIS) is expected to provide an innovative infrastructure for Earth sciences, enabling near-real-time data and processing services. Stereo-matching software has often been used in generating a Digital Elevation Model (DEM) from a pair of satellite imagery data sets to compute height from a parallax views using two photographic images. There is a need to reduce the computation time required for processing large images. We optimize stereo-matching software on multi-core/many-core processors, including Xeon, Cell and GPGPU. We describe optimization approaches of the correlation calculation part which occupied about 55% of the ovarall computation time. After porting and optimizing software for multi-core/many-core processors, we achieved processing time of 4.79 second (Xeon), 2.28 second (Cell) and 0.97 second (GPGPU) on each platform.


Introduction
GEO Grid [1] is a system and a concept that enables users to handle integrated processing and analysis of geographical information system data as well as various other observation data located in distributed environments, by the large-scale archiving and high-speed processing of Earth observation satellite data using grid technology. Its objective is to provide an environment in which research communities and companies can safely and securely use diverse data and make the required computations for Earth observations. One of the most important and interesting research issues is to implement the high-performance geospatial information system (GIS), which is expected to provide an innovative infrastructure for Earth sciences, enabling near-real-time data and processing services. For example, a Digital Elevation Model (DEM), which is created from a pair of satellite imagery data sets, is often used by Earth scientists to understand changes in topology. However, it normally takes from 10-15 minutes to as long as 4-6 hours to create a DEM of 4000 × 4000 pixels or 19000 × 19000 pixels, respectively, on a commodity PC. Reducing the computation time enables implementation of a near-real-time IT infrastructure based on Service Oriented Architecture (SOA). Optimization of signal processing for radio astronomy on multi-core/many-core processors has been presented [2]. In this paper, we optimize stereo-matching software, which has often been used in generating a DEM on multi-core/many-core processors, including Intel Xeon X5500 (Nehalem), IBM PowerXCell 8i (Cell), NVIDIA Tesla C1060 and C2050. Table 1 shows the platforms used in the experiments. Details of each platform are described as follows.

Intel Xeon 5500 series (Nehalem)
The Xeon 5500 series made by the Intel Corporation has an Intel x86 architecture CPU comprising a Nehalem core. It is a homogeneous multi-core processor that integrates four cores having the same performance. Because the same program operates in all cores for parallel processing, development is comparatively easy. SSE4.2 has instructions that are advantageous for the census transform method described in Section 4.2. The HP Z800 Workstation (hereafter "Z800") is used for evaluation in this paper.

IBM PowerXCell 8i
The PowerXCell 8i made by the IBM corporation is a processor with Cell architecture, and it is suitable for media processing. Its double precision floating-point arithmetic performance is five times faster than the Cell Broadband Engine (Cell/B.E.). Because double-precision arithmetic is necessary in a DEM generator program, it is more suitable than Cell/B.E. Moreover, because Cell/B.E uses a special memory (XDR DRAM), it has a limitation in terms of memory capacity. As the PowerXCell 8i uses general-purpose DDR2 DRAM, its memory capacity can be enlarged, which makes it suitable for large-scale computation.
The PowerXCell 8i is a heterogeneous processor. It comprises one PPE core that executes general purpose instructions and eight SPE cores that execute SIMD instructions at high speed. Because the PPE is not that fast, the processing is distributed to SPE to speed it up. DMA transfer to 256KB of local storage is required for parallel programming in SPE. Optimizing on this processor requires a lot of effort compared with a homogeneous multi-core processor such as Xeon. In this paper, IBM BladeCenter QS22 (hereafter "QS22") equipped with two CPUs is used for the evaluation.

NVIDIA Tesla C1060 and C2050
Recently, a GPGPU that uses a GPU for a general-purpose programming is attracting attention because it obtains good performance. A NVIDIA Tesla C1060 and C2050 (hereafter "C1060"and "C2050", respectively) are installed in the Z800, and it used for the evaluation in this paper. The C1060 comprises 30 SMs (streaming multiprocessors) equipped with eight SPs (streaming processors) that operate the scalar, and has a total of 240 operation cores. The C2050 comprises of 14 SMs equipped with 32 CUDA cores that operate the scalar, and has a total of 448 operation cores.
CUDA [8] is used for programming. It divides a process into a large number of threads. It processes an operate in a unit called a grid, which collects blocks, and blocks that collects threads. It can execute two or more blocks at the same time if conditions are suitable. The number of threads that can be executed in one Block is decided by the resource that each thread uses. Moreover, it is necessary to use shared memory in each SM to increase processing speed.

Stereo matching
Stereo matching is processing that searches for the same points in two images using different view points, and restructures depth information using the principle of triangulation. In the case of satellite image processing, it is processing that searches for the same points in under and rear sensor images, and converts them into altitude information.
The method used to speed up the stereo matching program in this paper consists of the preprocessing part, the main processing part and the postprocessing part. The preprocessing part consists of data input, memory allocation and initializing. The postprocessing part consists of data output. Stereo matching in the main processing part consists a correlation calculation, a median filter, and interpolation of the unexpected value. Fig. 1 shows the outline of the loop structure in the main processing part. The first loop employs coarse-tofine method. It roughly searches for a correspondence point by using the coarse image, which is a low resolution image. The parameter of the range is changed using the result of the previous loop. It gradually searches for the correspondence point for the switch to a fine image, which is a high resolution image. The substage loop searches again for low correlation points. There is a loop-carried dependence in the coarse-to-fine method loop and the substage loop. There is no loop-carried dependence in the y axis loop or the x axis loop, and these can operated in parallel. However, when a speed increase by the recursive calculation described in Section 4.1 is implemented, the x axis loop has loop-carried dependence.
The image data used to evaluate the optimization of the stereo matching program is shown in Fig. 2. The size of these images is 3962 × 4200 pixels. These images are provided by the ASTER satellite sensor [3]. They are images around the Mauna Kea on the Island of Hawaii. Because a range from the sea level is covered from 0m to 4,205m in the same images, the amount of the parallax is very suitable as bench-mark data.
We measured the processing time using the Z800 platform and the program before it was optimized. The correlation calculation occupied about 55% of the ovarall processing time. Then, we assumed that we had emphatically optimized the correlation calculation part.

Correlation methods
The search processing for the same point in stereo matching using the correlated calculation is described here. ZNCC is often generally used for the stereo matching of satellite images as a correlation calculation method. However, there is a problem with a high computational load. Methods such as using a low-load formula with preprocessing [4], and the recursive calculation of NCC [5] are reported to resolve this problem. In this paper, we describes a method of speeding up ZNCC using recursive calculation.
Another correlation calculation method is the census transform method [6] used for the roof detection of buildings [7]. Because an increase in speed of this method using SIMD etc. is obtaind, the outline is described here.
One image is expressed as I and another image is expressed as J. When the point in the image is expressed by vector q, I(q) express the brightness of image I in point q. The template is the rectangle region of image I. The program searches parallax vector p by measuring the correlation between the template and image J. X, Y is the rectangular size of the template. r = (x r , y r ) is coordinate on the top left of the template. Set W conprises coordinates in the template as follows: (1)

Recursive calculation in ZNCC
ZNCC C ZNCC (p, W) is as follows: I,J are the averages of the rectangle region as follows: |W| (= X × Y) is the number of elements in set W. The computational load of ZNCC is high because of the repetition of the sum of products. NCC (Normalized Cross Correlation) is an adjusting of the averages of ZNCC to 0. It is reported that NCC can greatly decrease computational complexity by recursive calculation [5]. NCC C NCC (p, W) is as follows: N, Q I , Q J are following: Q J (p, W) = q∈W J 2 (p + q).
Q I does not depend on p. The computational load to calculate N, Q J is high because it has to be calculated every time.
Here, we think about W as shown in Fig. 3. W \ W express a difference set between W and W . The value of N, Q J at region W is preserved in the memory, N, Q J at new region W can be calculated by calculating only the difference value as follows: A similar calculation can be done also for the y axis. As a result, the operation frequency can be decreased to 1/|W|. The recursive calculation in ZNCC is computable as well as the case of NCC. Eq. (2) represents M, D I , D J are following: This calculation is well known as a formula of covariance. The calculation of the averagesĪ,J can use recursive calculation. Elements M, D I , D J are calculated from N, Q I , Q J ,Ī,J which can use recursive calculation.

Census transform
The census transform method [6] uses the bit string of 0 and 1 that shows the bigness and smallness of the brightness value with a pixel at the center, and uses the Hamming distance. It is known to be robustness to brightness change and occlusion. Details of the method are as follows: Hamming{a, b} means Hamming distance between a and b. s is a vector of center coordinates on the template; s = p + (X/2, Y/2). It is said that census transform is faster than ZNCC which requires a product-sum operation, because it can be implemented only by big and small comparing and operating the bit string. However, because a value of center coordinates is used, a recursive calculation like ZNCC cannot be used.

Optimization on Xeon
In optimization on Xeon, the entire operation is parallelized using OpenMP first. Afterward, we optimize it by algorithm improvement or SIMD instructions. Fig. 4 shows the effectiveness of parallelization using OpenMP. The horizontal axis shows the number of threads and the vertical axis shows speed increase rates : ( speed increase rate ) = ( processing time in one thread ) / ( processing time in each thread ). It is increased around 6.5 times by using eight cores.
Next, we optimize the correlation methods using the recursive calculation of ZNCC, and a census transform using SIMD instructions. The processing time is shown in Table 2. ZNCC was greatly sped up using recursive calculation. A census transform is faster than ZNCC before optimization. We optimize a census transform using SIMD instructions. The POPCNT instruction that counts the bit is added from the Nehalem core of SSE4.2. It is efficient in computing a census transform.

Optimization on CELL
The optimization methods in Cell (IBM PowerXCell 8i) and their effectiveness are described as follows. In Cell, Processing speed of PPE which is the core for general processing is not fast. In paticular, the random access rate is very slow. However, the speed improves dramatically by the distribution of processing to SPEs and by using local storage that enables very fast access. However, it is necessary to perform optimization individually, and a lot of work is required compared with a Xeon. We implement optimization in the order shown in Table 3. The table shows the processing time and the speed-up rate in the process of each optimization phase. First, the performance improves greatly by paralleling the calculation part of ZNCC from PPE to SPEs. No too much of an improvement in performance is seen though processing SIMD instructions. In this optimization process, the DMA transfer time is a problem. The method of dividing the image in a block of constant width is used for latency hiding of a DMA transfer. Next, recursive calculation of ZNCC is implemented. The implementation of all these optimizations increases speed by a factor of about 765.

Optimization on GPGPU
We port the stereo matching program to the GPGPU using CUDA. The technique and the effectiveness of optimization are described as follows. The CUDA programing divides a process into many blocks. The image is divided as shown in Fig. 5, and we allocate one divided image in one block. The width of the partition is 32 pixels that is an efficiency value of the GPU programing. Height h decides the best value through measurement of actual performance.
Optimization includes loop-unrolling and using shared memory. The recursive calculation of ZNCC is implemented. The effectiveness of optimization on GPGPU is shown in Table 4. The processing time includes about 0.75 seconds for data transfer between the main memory on the host PC (Z800) and the global memory on the GPGPU board. The increase in speed obtained by the optimization is 2.8 to 4.28 times. In ZNCC, the recursive calculation increases speed by 6.12 times. The effect of the recursive calculation is less than the case on Cell and XEON, because the number of threads is decreased from 384 to 256.

Overall study
The final optimization results are shown in Table 5. With respect to the increased processing rate using the same processor, Cell is the highest because the original program porting on PPE was very slow. C2050 had the  fastest processing, but the speed increase rate was lowest because the original program was already parallelized using CUDA. The original C2050 program is very fast because C2050 includes a cache memory. Although, optimized program using shared memory is faster than original program. QS22 equipped with two processors showed almost equal performance to C1060 equipped with one processor. Z800 had about half the performance of QS22. A summary of the optimization approaches and the results is shown in Table 6. The findings show optimization in CELL requires of the most work, while because OpenMP can be used for automatic parallelization, the work required for Xeon is the least. The work required for GPGPU lies in the middle of these two.

Conclusion
We ported and optimized a stereo matching program in DEM generation using many-core processor platforms, and evaluated the differences and effectiveness of the optimization techniques on each platform. Two correlation methods, ZNCC and census transform, were evaluated in this paper. Moreover, a high-speed recursive calculation method was implemented in ZNCC. In optimization on Xeon, we used OpenMP, a recursive calculation on ZNCC, and SIMD instructions on census transform. In the optimization on Cell, we ported from a low-speed PPE to SPE, and optimizing the DMA transfer.
In optimization on GPGPU, the optimum size of block was properly requested. The performance of C2050 was 2.35 times (ZNCC) and 2.95 times (census) faster than C1060. The performance of two Cells was almost the same as C1060. The performance of two Xeons was half of two Cells.
Our future work will include clustering that assumes the use of high-resolution satellite images.