Abdulkadir Akin
email: abdulkadir.akin@epfl.ch

Luis Manuel Gaemperle
email: luis.gaemperle@epfl.ch

Halima Najibi
email: halima.najibi@epfl.ch

Alexandre Schmid
email: alexandre.schmid@epfl.ch

Yusuf Leblebici
email: yusuf.leblebici@epfl.ch

Enhanced Compressed Look-Up-Table Based Real-Time Rectification Hardware

Keywords: Stereo Matching, Image Rectification, Compression, Real-Time, Hardware Implementation, FPGA

Real-time disparity estimation requires real-time rectification which involves solving the models of lens distortions, image translations and rotations. Low complexity look-up-table based rectification algorithms usually require an external memory to store large look-up-tables. In this chapter, we present an implementation of the look-up-table based approach which compresses the rectification information to fit the look-up-table into the on-chip memory of a Virtex-5 FPGA. First, a very low complexity compressed look-up-table based rectification algorithm (CLUTR) and its real-time hardware are presented. The implemented CLUTR hardware rectifies stereo images with moderate lens distortion and camera misalignment. Moreover, an enhanced version of the compressed look-uptable based rectification algorithm (E-CLUTR) and its novel real-time hardware are presented. E-CLUTR solves more extreme camera alignment and distortion issues than CLUTR while maintaining the low complexity architecture.

Introduction

Disparity estimation (DE) is an algorithmic step that is applied in a variety of applications such as autonomous navigation, robot and driving systems, 3D geographic information systems, object detection and tracking, medical imaging, computer games, 3D television, stereoscopic video compression, and disparity-based rendering.

The stereo matching process compares the pixels in the left and right images and provides the disparity value corresponding to each pixel. If the cameras could be aligned perfectly parallel and the lenses were perfect, without any distortion, the matching pixels would be located in the same row of the right and left images. However, providing a perfect set-up is virtually impossible. Lens distortion and camera misalignments should be modeled and removed by internal and external stereo camera calibration and image rectification processes [START_REF] Bouguet | Camera Calibration Toolbox for Matlab[END_REF].

Image rectification is one of the most essential pre-processing parts of DE. Nevertheless, many real-time stereo-matching hardware implementations [START_REF] Chang | Algorithm and Architecture of Disparity Estimation with Mini-Census Adaptive Support Weight[END_REF][START_REF] Georgoulas | A Real-time Occlusion Aware Hardware Structure for Disparity Map Computation[END_REF][START_REF] Akin | A Hardware-Oriented Dynamically Adaptive Disparity Estimation Algorithm and Its Real-Time Hardware[END_REF] prove their DE efficiency using already calibrated and rectified benchmarks of the Middlebury evaluation set [START_REF] Scharstein | A Taxonomy and Evaluation of Dense Two-Frame Stereo Correspondence Algorithms[END_REF], while some do not provide detailed information related to the rectification of the original input images [START_REF] Greisen | An FPGA-based Processing Pipeline for High-Definition Stereo Video[END_REF].

In a system that processes the disparity estimation in real-time, image rectification should also be performed in real-time. The rectification hardware implementation presented in [START_REF] Son | A Rectification Hardware Architecture for an Adaptive Multiple-Baseline Stereo Vision System[END_REF] solves the complex equations that model distortion, and consumes a significant amount of hardware resources.

A look-up-table based approach is a straightforward solution that consumes a low amount of hardware resources in an FPGA or ASIC [START_REF] Vancea | LUT-based Image Rectification Module Implemented in FPGA[END_REF][START_REF] Gribbon | A Real-time FPGA Implementation of a Barrel Distortion Correction Algorithm with Bilinear Interpolation[END_REF][START_REF] Park | Real Time Rectification Using Differentially Encoded Lookup Table[END_REF]. In [START_REF] Vancea | LUT-based Image Rectification Module Implemented in FPGA[END_REF][START_REF] Gribbon | A Real-time FPGA Implementation of a Barrel Distortion Correction Algorithm with Bilinear Interpolation[END_REF][START_REF] Park | Real Time Rectification Using Differentially Encoded Lookup Table[END_REF], the mappings between original image pixel coordinates and rectified image pixel coordinates are precomputed and then used as look-up-tables. Due to the significant amount of generated data, these tables are stored in an external memory such as a DDR or SRAM [START_REF] Vancea | LUT-based Image Rectification Module Implemented in FPGA[END_REF][START_REF] Gribbon | A Real-time FPGA Implementation of a Barrel Distortion Correction Algorithm with Bilinear Interpolation[END_REF]. Using external storage for the image rectification process may increase the cost of the disparity estimation hardware system or impose additional external memory bandwidth limitations on the system. In [START_REF] Park | Real Time Rectification Using Differentially Encoded Lookup Table[END_REF], the look-up-tables are encoded to consume 1.3 MB data for 1280×720 size stereo images with a low-complexity compression scheme. This amount of data requires at least 295 Block RAMs (BRAM) without considering pixel buffers, thus it can only be supported by the largest Virtex-5 FPGAs or other recent high-end FPGAs.

In this chapter, we present a novel compressed look-up-table based image rectification (CLUTR) algorithm and its real-time hardware. A preliminary description of CLUTR has been presented in [START_REF] Akin | Compressed look-up-table based real-time rectification hardware[END_REF]. In this chapter, the CLUTR algorithm and its realtime hardware implementation are explained with further details. In addition, an enhanced version of the CLUTR (E-CLUTR) algorithm and its real-time hardware implementation are presented. The real-time hardware of E-CLUTR maintains the low complexity architecture of CLUTR, while it is able to rectify images under excessive mechanical misalignment of the cameras and lens distortions. Moreover, the Caltech rectification algorithm [START_REF] Bouguet | Camera Calibration Toolbox for Matlab[END_REF] which does not benefit from look-up-tables is implemented in hardware, and its hardware resource consumption results are presented to improve the hardware comparison and to evidence the efficiency of CLUTR and E-CLUTR much fairly.

This chapter is structured as follows. A typical look-up-table based rectification process is introduced in Section 2. The compression scheme that is used by the CLUTR algorithm is presented in Section 3. The proposed hardware implementation of the decompression process of CLUTR is presented in Section 4. The analysis of the challenges of CLUTR to rectify images under extreme conditions is presented in Section 5. The compression scheme of E-CLUTR is presented in Section 6. The hardware implementation of the decompression process of E-CLUTR is presented in Section 7. The implementation results of the CLUTR and E-CLUTR algorithms and their hardware resource comparisons are presented in Section 8. Section 9 concludes the chapter.

Typical Look-Up-Table based Solution

Look-up-table based rectification methods can be distinguished by two different image warping flows: forward mapping and inverse mapping. Forward mapping computes the rectified target pixel locations based on the given pixel locations in the original image. Inverse mapping computes the original source pixel locations based on the given pixel locations in the rectified image. The mapping requires separate tables for X and Y coordinates, and for the right and left images. Therefore, four tables are needed. The formulations for forward and inverse mappings are presented in equations ( 1) and ( 2), respectively. In these equations, ForwT is the forward mapping table, InvT is the inverse mapping table, Ori represents the original image taken from the camera, Rec represents the rectified image. YRec, XRec, Yori and Xori represent the Y and X coordinates.

Fig. 2. Forward mapping with fractional precision coordinates

A typical rectification process utilizes fractional pixel precision which requires the linear interpolation of four pixels. The linear interpolation schemes for forward and inverse mappings are represented in Figure 1 and Figure 2, respectively. The linear interpolation process for forward mapping is more complex than the linear interpolation process of inverse mapping, since it requires additional computations and an intermediate [START_REF] Bouguet | Camera Calibration Toolbox for Matlab[END_REF] (2) memory consumption to find the closest target pixels in the rectified image. The lookup-table based rectification hardware architectures presented in [START_REF] Vancea | LUT-based Image Rectification Module Implemented in FPGA[END_REF][START_REF] Gribbon | A Real-time FPGA Implementation of a Barrel Distortion Correction Algorithm with Bilinear Interpolation[END_REF][START_REF] Park | Real Time Rectification Using Differentially Encoded Lookup Table[END_REF] use inverse mapping due to its simplicity.

The size of the look-up-table depends on the size of the rectified image and the fractional precision. For example, for the rectification of 1024×768 resolution stereo images with 6 bits fractional precision, the rectification map alone requires approximately 6 MB of space in a memory. This amount of data is excessive to fit into the on-chip memory of a mid-range FPGA. Therefore, dumping look-up-tables into an external memory is preferred in the hardware implementations of [START_REF] Vancea | LUT-based Image Rectification Module Implemented in FPGA[END_REF][START_REF] Gribbon | A Real-time FPGA Implementation of a Barrel Distortion Correction Algorithm with Bilinear Interpolation[END_REF].

Proposed Compression Scheme of the CLUTR Algorithm

In contrast to the selection of the hardware implementations of [START_REF] Vancea | LUT-based Image Rectification Module Implemented in FPGA[END_REF][START_REF] Gribbon | A Real-time FPGA Implementation of a Barrel Distortion Correction Algorithm with Bilinear Interpolation[END_REF][START_REF] Park | Real Time Rectification Using Differentially Encoded Lookup Table[END_REF], a forward mapping based rectification scheme is selected for the proposed CLUTR algorithm. In CLUTR, fractional precision is ignored. Ignoring fractional precision allows an efficient compression scheme. The negligible distortion in the rectified images originating from this simplification is analyzed in Section 8.

The compression scheme is presented in the flow graph in Figure 3. The proposed compressed rectification algorithm produces four compressed tables. The compression scheme requires eight steps. The details of steps 1-2 can be found in [START_REF] Bouguet | Camera Calibration Toolbox for Matlab[END_REF]. The details of steps 3-8 are detailed in this section.

In the third step, integer coordinate precision forward mapping is extracted from the fractional precision inverse mapping. The extraction scheme is demonstrated in Figure 4. The example original and rectified pictures have a size of 4×5 pixels. First, inverse mapping is applied to find the fractional source pixel locations of all pixels in the rectified image. Due to the 3D rotation, some of the pixels in the rectified image cannot be related to their source pixels in the 4×5 original image, as shown in Figure 4(a). The nearest integer coordinates of all fractional source coordinates are computed, and they are targeted onto the integer pixel coordinates in the rectified image, as presented in Figure 4(b). Thus one-to-one mapping is provided in the third step.

The integer pixel precision forward mapping extracted for the example picture in Figure 4 yields the look-up-tables of X and Y coordinates shown in Figure 5. The pixels that are not targeted to any location are identified with NT. Ori [START_REF] Chang | Algorithm and Architecture of Disparity Estimation with Mini-Census Adaptive Support Weight[END_REF][START_REF] Chang | Algorithm and Architecture of Disparity Estimation with Mini-Census Adaptive Support Weight[END_REF] and Ori [START_REF] Chang | Algorithm and Architecture of Disparity Estimation with Mini-Census Adaptive Support Weight[END_REF][START_REF] Georgoulas | A Real-time Occlusion Aware Hardware Structure for Disparity Map Computation[END_REF] are adjacent pixels, and both of them target "row no 2" of the rectified image; Ori [START_REF] Chang | Algorithm and Architecture of Disparity Estimation with Mini-Census Adaptive Support Weight[END_REF][START_REF] Chang | Algorithm and Architecture of Disparity Estimation with Mini-Census Adaptive Support Weight[END_REF] and Ori [START_REF] Georgoulas | A Real-time Occlusion Aware Hardware Structure for Disparity Map Computation[END_REF][START_REF] Chang | Algorithm and Architecture of Disparity Estimation with Mini-Census Adaptive Support Weight[END_REF] are adjacent pixels and both target "column no 1" of the rectified image. This regular order is more apparent with higher resolution images. According to our experiments with a 1024×768 image, repetition of a single target coordinate up to 220 times is observed in the integer precision forward mapping table of X coordinates.

The method governing compressed rectification is similar to the run-length encoding technique. In the proposed coding scheme, instead of coding the run-length of the regular order, the locations where the regular order changes are encoded. These locations are called breakpoints. Moreover, the proposed scheme includes additional specific techniques to compress the integer precision forward mapping efficiently.

The regular order of the Y coordinate mapping is encoded following a row-by-row scheme, and the regular order of the X coordinate mapping is encoded following a column-by-column scheme. The resulting look-up-tables after encoding Figure 5(a) and Figure 5(b) are presented in Figure 6(a) and Figure 6(b). In Figure 6(a), the elements of the compressed table are represented as (column number, new value in row). In Figure 6(b), the elements of the compressed table are represented as (row number, new value in column).

The high number of NT pixels dramatically increases the number of breakpoints. This issue becomes more pronounced for high resolution images. Therefore, the fourth step of the compression algorithm fills the NT pixel locations to keep the regular order. In order to fill the NT pixel locations, the same order is repeated vertically and horizontally for Y and X locations, respectively. After the fourth step, Figure 5 is transformed into Figure 7, and Figure 6 into Figure 8. After the first two steps, two or more source fractional coordinates can have the same pixel coordinate in the original image as their nearest neighbor, as presented in Figure 9 (a). However, after step three and four, every integer pixel coordinate of the original image is targeted to a single coordinate in the rectified image. Consequently, some pixels in the rectified image may be void, as presented in Figure 9 (b). These pixels will remain as void, i.e. black pixels, in the rectified image if they are not filled. The fifth step is applied to fill these voids. As shown in Figure 10, the pixels on the original image which target the pixel coordinates that are located on the row above these voids are marked. Marked pixels are used to fill the voids as source pixels which have double targets (DT). Thus, DT pixels are used to concurrently target two vertically neighboring pixels on the rectified image.

Fig. 3. Flow chart for the proposed compressed look-up-table based stereo image rectification process

The sixth step of the algorithm extracts the breakpoint locations and analyzes the behavior of the breakpoints. As shown in Figure 8, the difference between the new and previous target locations equals plus or minus one, which can be encoded consuming less data than encoding the exact integer coordinates. An example of coding the behavior of the cells in Figure 8 is presented in Figure 11 as (location, behavior). The initialization coordinates are provided in the first column of the look-up-table for Y coordinates, and in the first row of the look-up-table for X coordinates. The next breakpoint values are identified with ±1. Moreover, dummy breakpoints are inserted at the edges of the image to simplify the hardware implementation. Dummy insertions are represented by (5,0) and (4,0) in Figure 11 In the seventh step, the locations and behaviors of the breakpoints are concatenated and stored in a data array. Every BRAM in a Virtex-5 FPGA has 1024 addresses and it can be configured to store one array composed of 1024×36bits or two arrays composed of 1024×18bits. The BRAMs of the FPGAs are configured to store 18-bits in each address in the proposed concatenation scheme. As shown in Figure 12, 3-bits are used for coding the behaviors, and the remaining 15-bits encode the locations of the breakpoints. Therefore, the proposed concatenation scheme can be applied to an image that has a resolution lower than 32767×32767 pixels. In Fig. 12 (a), DT and changing the last targeted row by ±1 are independent breakpoint conditions of Y coordinates, which can be applied to source pixels, concurrently or separately. Therefore, the "X" symbol in the -1 and +1 columns of Fig. 12 implies keeping the last targeted row coordinate. The number of breakpoints in every row of the Y table and the number of breakpoints in every column of the X table depend on the distortion of the lens, the resolution of the image sensor and the mechanical misalignment. The experimental setup used in this chapter consists of 1024×768 resolution cameras. In the experiments, cameras are aligned in parallel configuration without using any sensitive mechanical placement tool. At most 21 breakpoints are observed in any given row of Y tables, and at most 17 breakpoints are observed in any given column of X tables. Data arrays of CLUTR are created for 24 possible breakpoint locations for Y tables and 20 possible breakpoint locations for X tables to support more challenging distortion conditions. Therefore, storing the X and Y tables for the right and left images requires 38 BRAMs which can even be supported by low cost FPGAs. The data arrays that are programmed into the BRAMs are converted into coefficient (COE) files using MATLAB.

In the eighth step, 38 BRAMs are instantiated as single port ROMs. The pre-computed compressed rectification maps are programmed into the BRAMs using the Xilinx ISE 12.4 and COE files.

Real-Time Decompression Hardware of CLUTR

The decompression process is simpler than the off-line compression process in terms of computational complexity. The proposed rectification module can be used as a hardware accelerator taking place between the camera interface hardware and the on-chip memory controller, as shown in Figure 13. The rectification module is used for the left and right cameras separately. The rectification module processes source pixel values as Ori(Yori, Xori) and the respective source row and source column coordinates as Yori and Xori. The rectification module computes the target row and target column coordinates as YRec,and XRec, and the 1-bit DT signal to identify double targeted locations. Ori(Yori, Xori) is delayed for 6 clock cycles and Rec(Yrec, Xrec) is given as an output. Due to the pipelined structure of the hardware, inputs can be consecutively received and outputs can be consecutively provided. the compressed table of X coordinates. Half of 1 additional BRAM is used to store the last break point locations and the last target X coordinates of the row which is located above the row currently being processed.

The block diagram of the decompression hardware of Y coordinates is presented in Figure 15. The hardware resets itself every time Xori is equal to zero which implies that the first pixel of a new row is fetched from the camera. The target Y coordinate of the first incoming pixel in a new row is loaded from the ROM and written to the output register of YRec. For every consecutive pixel, Xori is compared to the coordinate of the next breakpoint which is loaded from the ROM. When a breakpoint is reached, the YRec value is changed using a multiplexer depending on the coded behaviors of the breakpoints. Meanwhile, the hardware loads the coordinate of the next breakpoints to compare with the upcoming Xori.

The block diagram of the decompression hardware of X coordinates is presented in Figure 16. Pixels are supplied by the camera row-by-row, whereas the X coordinates are compressed column-by-column. This situation causes one important difference between the decompression hardware architectures of the X and Y tables. When the camera provides pixels of a new row, the decompression hardware needs to keep record of the previous XRec coordinates and the last checked breakpoint address in the ROM for the respective column of the previous row. Two 1×1024 size data arrays are needed to store this information. These arrays are named array_last_break_x and array_last_target_x in Figure 14. These arrays are concatenated for respective column coordinates of the original image, and stored into one half of the 1 BRAM, which is named X_last_data_BRAM in Figure 16. The values in X_last_data_BRAM are replaced with the new ones when a breakpoint is reached for the respective XOri. The decompression hardware of the Y coordinates does not comprise these arrays because Y coordinates are compressed row-by-row. Therefore, the last YRec can be directly used for computing the next YRec of the next pixel in the same row of the original image. The decompression hardware of X coordinates operates in a similar fashion as the decompression hardware of Y coordinates, with the exception of keeping record of the information about the previous row.

Fig. 16. Block diagram of the proposed rectification hardware for decompressing the table for X coordinates

The proposed rectification hardware can be used in any stereo-matching system. The stereo matching process can be started when the required amount of rows is buffered in the BRAMs of the stereo matching hardware. Processed rows in these BRAMs can be overwritten by new rows during the stereo matching process.

The hardware architectures presented in [START_REF] Gribbon | A Real-time FPGA Implementation of a Barrel Distortion Correction Algorithm with Bilinear Interpolation[END_REF][START_REF] Park | Real Time Rectification Using Differentially Encoded Lookup Table[END_REF] require large pixel buffers due to the inverse mapping scheme. The proposed decompression does not need large pixel buffers between the camera interface and the rectification modules. In contrast, the hardware requires these pixels buffers for the rectified image. However, typically DE hardware implementations already include BRAMs to buffer the pixels [START_REF] Chang | Algorithm and Architecture of Disparity Estimation with Mini-Census Adaptive Support Weight[END_REF][START_REF] Georgoulas | A Real-time Occlusion Aware Hardware Structure for Disparity Map Computation[END_REF][START_REF] Akin | A Hardware-Oriented Dynamically Adaptive Disparity Estimation Algorithm and Its Real-Time Hardware[END_REF][START_REF] Greisen | An FPGA-based Processing Pipeline for High-Definition Stereo Video[END_REF]. Therefore, these buffers can be used for the proposed decompression hardware. Thus, using the proposed rectification hardware on a complete DE system may not need additional large pixel buffers.

Limitations of the CLUTR

In an ideal case, i.e. where the cameras are perfectly parallel and lenses do not have distortion, two breakpoints are required for every row of the look-up-table for Y coordinates and two breakpoints are required for every column of the look-up-table for X coordinates. One of these two breakpoints is needed to define the initial breakpoint location to target coordinate 0, and the other one is needed to define the final target location as the horizontal or vertical size of the image. In the classical case of a real-time working environment, the mechanical set-up of the stereo-matching system should be carefully designed to be close to an ideal case. Still the main goal of the rectification consists of solving lens distortions and sensitive mechanical misalignments.

According to tests applied to the CLUTR algorithm and hardware, the pixel location difference of two consecutive breakpoints typically reaches more than 15 pixels and the number of breakpoints is smaller than the pre-defined breakpoint capacity of the ROMs of CLUTR. Therefore, CLUTR successfully rectifies the images when the lens distortion and the mechanical misalignments are not excessive. However, unusual conditions bring limitations on the CLUTR hardware.

Two important limitations of the CLUTR hardware must be considered to maintain its suitability to rectify challenging situations. The first limitation relates to the capacity of the ROMs to store a sufficient number of breakpoints. The second limitation relates to the frequency of the breakpoints.

The limitation caused by the number of breakpoints is mainly due to the mechanical misalignment of the cameras. In order to identify the limit of the breakpoint storage capacity of ROMs, two cameras are manually rotated around 3 degrees around opposite directions of all rotational axis. This test can be considered as an excessive misalignment of a carefully designed mechanical setup of the stereo-matching system. Using the compression scheme of CLUTR, 43 breakpoints are needed in the look-up-table of X coordinates for one column, and 69 breakpoints are needed in the look-up-table of Y coordinates for one row. The pre-defined breakpoint capacity of CLUTR does not support this condition. Overcoming this first limitation is straightforward to achieve by increasing the size of ROMs to store more breakpoints.

CLUTR supports rectification if the two breakpoints of the Y coordinates have at least 4 pixel position difference. When a breakpoint location is reached, the hardware needs to read the next breakpoint from ROM. The address computation and reading the next breakpoint from the ROM consume 4 clock cycles. The camera continues to send pixels and the camera controller increases XOri during the address computation and reading breakpoints from the ROM. Therefore, if there are multiple breakpoints in 4 consecutive pixels, CLUTR is not able to apply a breakpoint condition to those pixels. Hence, the limits of the CLUTR hardware to successfully rectify stereo images is exceeded if breakpoints are frequent, i.e. if two breakpoints of Y coordinates in a row have less than 4 pixel position difference. Since the breakpoints of X coordinates are coded column by column but the camera sends pixels row by row, the time to process consecutive breakpoints in same column and consecutive rows is sufficiently long. Therefore, frequent breakpoints in the same column of the look-up-table of X coordinates does not cause a limitation. The frequency limitation of CLUTR is visually explained in Fig. 17.

The main reason for the occurrence of frequent breakpoints is the high number of adjacent void pixels, which is caused by excessive camera misalignment or lens distortion. An example of one such challenging condition is presented in Fig. 18. As shown in Fig 18 (b), 4 out of 5 consecutive pixels are marked as breakpoints. 3 out of these 4 breakpoints are DT breakpoints which are coded in the look-up-table of Y coordinates, and the other one is located at Ori [START_REF] Bouguet | Camera Calibration Toolbox for Matlab[END_REF][START_REF] Scharstein | A Taxonomy and Evaluation of Dense Two-Frame Stereo Correspondence Algorithms[END_REF] which requires changing the target row number from 2 to 3. This challenging example case exceeds the limits of CLUTR, since CLUTR is not able to apply a breakpoint condition if there are multiple breakpoints in 4 consecutive pixels. Adjacent void pixels may occur not only horizontally but also vertically. Vertically adjacent void pixels may cause void pixels, which cannot be filled by the DT feature of CLUTR. As visualized in Fig. 19, DT pixels can fill the voids located one row below the targeted pixel on the rectified image. If there are two voids which are vertically adjacent, the void below can not be filled by CLUTR since the pixel above is not targeted directly by any source pixel.

Another limitation of CLUTR is related to the usage of ROMs for the hardware implementation which increases the off-line processing duration. Using ROMs is suitable to demonstrate the efficiency of compressed look-up-table based rectification. However, after each adjustment of the camera settings and alignments, creating new compressed tables of the hardware requires re-synthesis and place & route of the implementation. Thus it takes a long time to initialize CLUTR hardware.

6

The

compression scheme of the E-CLUTR algorithm

The E-CLUTR algorithm and its hardware implementation are designed to overcome the limitations of CLUTR while maintaining the low complexity decompression scheme. The limitations of CLUTR are mainly solved by improving the design of the decompression hardware. Moreover, algorithmic enhancements are applied to further improve the efficiency of the compression scheme to handle challenging lens distortions and mechanical misalignments.

Algorithmic enhancements are explained in this section. The flow chart of the compression scheme of E-CLUTR is identical to the flow chart presented in Fig 3 . The steps 5, 6 and 7 that are shown in Fig. 3 are enhanced in E-CLUTR. The algorithmic enhancement for the compression scheme is proposed to reduce the frequency of DT breakpoints. In order to decrease the amount of consecutive breakpoints in a row, the condition type of breakpoints for filling the voids are improved. As explained in Section 4 in step 5, DT breakpoints are used to fill the voids that are located one pixel below the targeted pixel. To avoid any confusion, the DT condition of CLUTR is renamed as below-DT (B-DT) in E-CLUTR. In addition to B-DT, below-backward-DT (BB-DT), below-forward-DT (BF-DT), upper-DT (U-DT), upper-backward-DT (UB-DT) and upper-forward-DT (UF-DT) breakpoint conditions are defined in E-CLUTR. Using extra DT conditions, the source pixel can be targeted not only to one pixel below the target, but additional options are provided to fill any of 6 possible neighbors of the targeted pixel of the rectified image. These additional options are visualized in Fig. 20.

As presented in Fig. 21, the frequency of breakpoints in the same row is reduced compared to Fig. 18, by using BF-DT, BB-DT and U-DT breakpoint conditions. As presented in Fig. 19, vertically adjacent void pixels are problematic for CLUTR. However, these voids can be filled by multiple DT options of E-CLUTR as presented in Fig. 22.

In step 6, the breakpoints are coded considering the new DT breakpoint conditions. In the challenging example, there should be support of 2 consecutive breakpoints at least in 3 consecutive pixel coordinates, as presented for the breakpoints at Ori [START_REF] Bouguet | Camera Calibration Toolbox for Matlab[END_REF][START_REF] Akin | A Hardware-Oriented Dynamically Adaptive Disparity Estimation Algorithm and Its Real-Time Hardware[END_REF] and Ori [START_REF] Bouguet | Camera Calibration Toolbox for Matlab[END_REF][START_REF] Scharstein | A Taxonomy and Evaluation of Dense Two-Frame Stereo Correspondence Algorithms[END_REF]. Therefore, the algorithmic enhancement requires the support of at least 2 breakpoints for 3 consecutive locations as an additional constraint. The hardware based enhancement to provide this support is explained in section 7.

At the edges of the rectified image, black pixels occur. These stem from the 3-D rotation of the original image and the mapping of the rectified image to the original resolution. Consequently, the effective resolution slightly decreases in the rectified image [START_REF] Bouguet | Camera Calibration Toolbox for Matlab[END_REF]. Due to this fact, changing row and column coordinates stay in the range of ±1. However, this range may not be guaranteed for all possible extreme conditions. Thus, a generic solution should cover all possible extreme situations. In order to cover the cases that create a situation beyond the challenging camera misalignment tests, ±2 row and ±2 column coordinate change options are included as a breakpoint condition in step 6 of E-CLUTR.

The concatenation scheme presented in Fig. 12 of CLUTR is modified for E-CLUTR as presented in Fig. 23, Fig. 24 and Fig. 25 using improved breakpoint conditions. The breakpoint conditions for DT conditions and row changing can be applied concurrently or separately to the source pixel. Consequently, 7×5=35 different conditions occur for the concatenation of the conditions for the breakpoints of Y coordinates. The brief concatenation scheme of the Y coordinates is presented in Fig. 23. The concatenation scheme for DT codes and row changing are separately presented in Fig. 24 (a) and Fig. 24 (b). The concatenation scheme of the X coordinates is presented in Fig. 25. As presented in Figure 27, the decompression hardware of E-CLUTR pertaining to X coordinates is similar to the hardware used in CLUTR. The multiplexing stage is adapted to provide the ±2 target pixel column change feature for the X breakpoints.

The decompression hardware of E-CLUTR pertaining to Y coordinates is redesigned to support frequent breakpoints and the six different DT options. As presented in Figure 28, the E-CLUTR hardware reads the first six breakpoints from the RAM as soon as the camera starts to send a new row. The first breakpoint is used to initialize the target (a) (b) row and the next five breakpoints are buffered in a local cache. Whenever a new breakpoint location is reached, the next breakpoint location is read from the RAM and the cache shifts the existing upcoming breakpoint locations. Using this local cache of the breakpoints, the original pixel coordinates can be compared to the pixel locations in the cache. Therefore, the breakpoint conditions can be applied to all passing pixels even if the breakpoints are frequent.

The multiplexing stage for the computation of the next target row is improved to provide the ±2 target pixel row change feature. Moreover, the hardware sends the 3-bit DT condition to the BRAM controller together with the pipelined source pixel and its target locations YRec,and XRec, synchronously.

The BRAM controller that is shown in Fig. 13 writes the pipelined source pixels to the decompressed target row and target column coordinates. E-CLUTR hardware is verified by merging it with the DE hardware presented in [START_REF] Akin | A Hardware-Oriented Dynamically Adaptive Disparity Estimation Algorithm and Its Real-Time Hardware[END_REF], which buffers pixels of rows in its own, separate BRAMs. YRec is used to select and enable the BRAM to write target pixel. XRec is used to determine the write address of the enabled BRAM of the DE hardware. If a DT condition exists, the same source pixel is written to two BRAMs concurrently by enabling two BRAMs that buffer two consecutive rows. If the DT condition is pointing into forward or backward positions, the address port of the BRAM that is targeted by the DT condition receives the computed target BRAM address ±1.

Implementation Results

The proposed rectification hardware architectures of CLUTR and E-CLUTR are implemented using Verilog HDL, and verified using Modelsim 10.1d. The Verilog RTL models are mapped to a Virtex-5 XCUVP-110T FPGA comprising 69k Look-Up-Tables (LUT), 69k DFFs and 148 BRAMs. One rectification module of CLUTR consumes 0.32% of the LUTs, 0.28% of the DFF resources and 14% of the BRAM resources of the Virtex-5 FPGA. One rectification module of E-CLUTR consumes 0.63% of the LUTs, 0.51% of the DFF resources and 38% of the BRAM resources of the Virtex-5 FPGA. The proposed E-CLUTR hardware operates at 212 MHz after place & route.

Therefore, it can process up to 269 fps at a 1024×768 XGA video resolution. In addition, the proposed rectification hardware of CLUTR and E-CLUTR are merged with the DE hardware presented in [START_REF] Akin | A Hardware-Oriented Dynamically Adaptive Disparity Estimation Algorithm and Its Real-Time Hardware[END_REF]. The merged DE systems are also verified using Modelsim 10.1d.

The proposed rectification hardware of CLUTR and E-CLUTR do not need the support of external memory if the cameras are synchronized. The cameras can be perfectly synchronized by driving the cameras with same clock source and using one common I2C module for the initialization the cameras [START_REF] Akin | Hemispherical Multiple Camera System for High Resolution Omni-Directional Light Field Imaging[END_REF].

The proposed compression and decompression algorithms are evaluated using the pictures taken by the stereo camera system presented in [START_REF] Akin | A Hardware-Oriented Dynamically Adaptive Disparity Estimation Algorithm and Its Real-Time Hardware[END_REF]. If the cameras are not extremely misaligned, CLUTR and E-CLUTR provide identical visual and numerical results. The 1024×768 size original left and right pictures that are shown in Figure 29 are taken under a camera misalignment condition that exceed the limits of CLUTR implementation. The example pictures and the test results of CLUTR for typical camera misalignment are presented in [START_REF] Akin | Compressed look-up-table based real-time rectification hardware[END_REF]. The original images in Figure 29 are rectified using the Caltech rectification algorithm [START_REF] Bouguet | Camera Calibration Toolbox for Matlab[END_REF] and the proposed E-CLUTR algorithm. The rectification results of the E-CLUTR is presented in Figure 30. The extreme rotation of the rectified image can be visually observed in Figure 30. The breakpoint locations for the X and Y coordinates of the left image are presented in Figure 31.

The PSNR between the rectification results of E-CLUTR and Caltech rectification algorithm are evaluated in Table I. The PSNR of the left image is 43.10 dB, and the PSNR of the right image is 42.02 dB. Generally, a PSNR larger than 30 dB is considered acceptable to the human eye. Therefore, E-CLUTR provides very high quality rectification results. The PSNR between the original images and the rectification results of Caltech are also provided in Table I for comparison. The hardware implementation of the E-CLUTR is compared with the stereo image rectification hardware implementations in Table II. The hardware architecture of [START_REF] Son | A Rectification Hardware Architecture for an Adaptive Multiple-Baseline Stereo Vision System[END_REF] requires a significant amount of hardware resources to support complex operations for solving the lens distortion models. Hardware architectures of look-up-table based implementations [START_REF] Vancea | LUT-based Image Rectification Module Implemented in FPGA[END_REF] and [START_REF] Gribbon | A Real-time FPGA Implementation of a Barrel Distortion Correction Algorithm with Bilinear Interpolation[END_REF] require an external memory. Combining the E-CLUTR with BRAM controller consumes less LUT and DFF resources than [START_REF] Son | A Rectification Hardware Architecture for an Adaptive Multiple-Baseline Stereo Vision System[END_REF][START_REF] Vancea | LUT-based Image Rectification Module Implemented in FPGA[END_REF][START_REF] Gribbon | A Real-time FPGA Implementation of a Barrel Distortion Correction Algorithm with Bilinear Interpolation[END_REF] and it does not require an external memory. The DFF and LUT consumption of [START_REF] Park | Real Time Rectification Using Differentially Encoded Lookup Table[END_REF] is not available (NA). Nevertheless, the capacity of E-CLUTR to fit the look-up-tables into the on-chip memory of the Virtex-5 FPGA is approximately two times more efficient than [START_REF] Park | Real Time Rectification Using Differentially Encoded Lookup Table[END_REF], as a benefit of its efficient compression scheme. Moreover, high precision hardware of the Caltech rectification algorithm is implemented, and its hardware consumption results are presented in Table II The hardware resource consumption of E-CLUTR is higher than CLUTR hardware. However, if the cameras are extremely misaligned, the limitations of CLUTR can be exceeded. In these extreme conditions, E-CLUTR still supports rectification. Whereas, using CLUTR hardware can be more profitable if the stereo cameras are carefully aligned.

Conclusion

In this chapter, two novel compressed look-up-table based image rectification algorithms and their hardware implementations are presented. The proposed CLUTR and E-CLUTR algorithms are based on off-line compression of the rectification information to fit the tables into the on-chip memory of a Virtex-5 FPGA. The presented decompression hardware implementations of CLUTR and E-CLUTR consume negligible amounts of hardware resources, and they do not require any external memory to store the look-uptables. The proposed hardware implementations are advantageous if using external memory is considered as an additional cost, or if the disparity estimation system has external memory bandwidth limitations. The proposed rectification hardware implementations would be even more profitable if they are adapted for high resolution multiple camera disparity estimation systems.

Fig. 1 .

 1 Fig. 1. Inverse mapping with fractional precision coordinates. Corners indicate integer pixel coordinates.

Fig. 4 .Fig. 5 .

 45 Fig. 4. Third step of the compression flow (a) selection of nearest source pixels from fractional inverse mapping (b) extraction of forward mapping with integer coordinates

Fig. 6 .Fig. 7 .Fig. 8 .Fig. 9 .Fig. 10 .Fig. 11 .

 67891011 Fig. 6. Coded regular orders after the third step (a) coded mapping of Y coordinates (b) coded mapping of X coordinates

  (a) and Figure 11 (b), respectively.

Fig. 12 .

 12 Fig. 12. Concatenation of the locations and behaviors at the seventh step (a) for the mapping of Y coordinates (b) for the mapping of X coordinates

  The top-level block diagram of the rectification module is presented in Figure14. The rectification module involves (768×24×18)/(1024×36) = 9 BRAMs to store the compressed table of Y coordinates and (1024×20×18)/(1024×36) = 10 BRAMs to store (a) (b)

Fig. 13 .

 13 Fig. 13. Example utilization of the proposed rectification hardware

Fig. 17 .Fig. 18 .Fig. 19 .

 171819 Fig. 17. Visualization of the breakpoint frequency capacity of the X and Y coordinate mappings (a) breakpoints for the mapping of Y coordinates (b) breakpoints for the mapping of X coordinates

Fig. 20 .Fig. 21 .

 2021 Fig. 20. Filling the voids on the rectified image in the fifth step (a) DT option of CLUTR (b) DT options of E-CLUTR

Fig. 23 .Fig. 24 .

 2324 Fig. 23. Brief representation for the concatenation of the locations and behaviors at the seventh step of E-CLUTR for the mapping of Y coordinates

Fig. 26 .

 26 Fig. 26. Top-level block diagram of the proposed rectification hardware of E-CLUTR

Fig. 29 .Fig. 30 .

 2930 Fig. 29. Original images have distortions as observed on the lines (a) left image (b) right image. (a) (b)

Fig. 31 .

 31 Fig. 31. Breakpoint locations of the left image (a) breakpoints of the targeted Y coordinates; coded row-by-row. (b) breakpoints of the targeted X coordinates; coded column-by-column

  for comparison. Hardware implementations of CLUTR and E-less hardware resource than the hardware implementation of Caltech rectification while providing almost same rectification results.

Table 1 .

 1 PSNR (𝑑𝐵) with the Rectified Images Produced By[START_REF] Bouguet | Camera Calibration Toolbox for Matlab[END_REF]

	Comparison

with Rectified Left Image [1] Comparison with Rectified Right Image [1]

	Original Image	17.99	19.34
	Proposed (E-CLUTR)	43.10	42.02

Table 2 .

 2 Hardware Resource Comparison of the Rectification Hardware Implementations

		Device	Resolution	LUT	DFF	On-Chip Memory (KB)	External Memory
	[7]	Virtex-4	752×480	3418	5932	0	✔
	[8]	Virtex-E	640×512	2459	2075	99	✔
	[9]	Spartan-2	640×480	≈2396	≈2396	16	✔
	[10]	Virtex-5	1280×720	NA	NA	1300	X
	Caltech Hardware Virtex-5	1024×768	24384	25346	0	✔
	CLUTR	Virtex-5	1024×768	227	197	90	X
	2×(CLUTR + BRAM Contr.)	Virtex-5	1024×768	784	427	176	X
	E-CLUTR	Virtex-5	1024×768	434	350	252	X
	2×(E-CLUTR + BRAM Contr.)	Virtex-5	1024×768	2278	956	500	X


Andreas Minwegen
email: minwegen@ice.rwth-aachen.de

Dominik Auras
email: auras@ice.rwth-aachen.de

Gerd Ascheid
email: ascheid@ice.rwth-aachen.de

A Flexible ASIC for Time-Domain Decision-Directed Channel Estimation in MIMO-OFDM Systems

Keywords: MIMO-OFDM, VLSI, ASIC, channel estimation, expectation maximization, SAGE

Channel estimation is a crucial task for the overall communication performance of a wireless receiver. Compared to traditional approaches the estimation of the wireless channel can be improved by using iterative estimation with feedback from other receiver components, however the VLSI implementation of such iterative channel estimation in multiple-input multiple-output (MIMO) orthogonal frequency division multiplexing (OFDM) systems is challenging due to the high computational complexity. In this chapter we introduce the first ASIC for Decision-Directed MIMO-OFDM channel estimation which tracks channel variations using feedback from a decoder and supports M-QAM. Furthermore, timing and power dissipation trade-offs are analyzed.

Orthogonal frequency division multiplexing (OFDM) and spatial multiplexing over multiple-input multiple-output (MIMO) transmissions schemes are adopted by several recent wireless communication standards such as 3GPP Long Term Evolution (LTE) or IEEE 802.11n and beyond. Due to the concept of coherent detection in MIMO-OFDM receivers the channel estimation (CE) is a crucial and computational intensive part of the overall system and has a significant impact an the communication performance in terms of frame error rate and thus influence directly the maximal achievable throughtput. Traditionally, pilot-aided channel estimation (PACE) is applied where the channel is estimated at predefined pilot positions. The complete channel over all subcarriers is obtained via interpolation. Iterative channel estimation can be used to improve the estimates which delivers promising SNR gains [START_REF] Gao | Low-complexity MAP channel estimation for mobile MIMO-OFDM systems[END_REF]. In the case of iterative channel estimation algorithms, a priori knowledge from the detector or the decoder is used to improve the channel estimates iteratively.

However, the gain in terms of algorithmic performance is payed for in terms of high latencies since each block that participates in the iterations has to complete its processing before the next one can start using the improved input. An alternative that does not increase the latency but still provides significant benefits is the channel tracking approach. This approach ( [START_REF] Li | Channel estimation for OFDM systems with transmitter diversity in mobile wireless channels[END_REF] and [START_REF] Ylioinas | Iterative joint detection, decoding, and channel estimation in turbo-coded MIMO-OFDM[END_REF]) provides channel estimation updates for time instance n + U x based on the detector or decoder decisions for time instance n, as shown in Fig. 1.

Especially, fast fading channels are interesting scenarios for such a solution. An algorithmic investigation for a communication system similar to LTE is performed in [START_REF] Xie | Two EM-type channel estimation algorithms for OFDM with transmitter diversity[END_REF] using the simplified frequency domain (FD) SAGE (spacealternating generalized expectation-maximization) algorithm, which calculates updates of the channel impulse response (CIR) estimates for M-QAM constellations. The authors of [START_REF] Ylioinas | Iterative joint detection, decoding, and channel estimation in turbo-coded MIMO-OFDM[END_REF] present a modification that provides a gain of 2 dB for 64-QAM. This modification requires a non trivial matrix inversion. Fortunately, this matrix inversion can be avoided by iteratively processing each tap of the CIR in the time domain, as proposed by the time-domain (TD) SAGE algorithm presented in [START_REF] Ylioinas | Avoiding matrix inversion in DD SAGE channel estimation in MIMO-OFDM with M-QAM[END_REF].

Furthermore, the authors of [START_REF] Ketonen | Decision directed channel estimation for improving performance in LTE-A[END_REF] present an analysis of the TD-SAGE algorithm in the context of LTE-Advanced. The results show that it has the potential to double the system throughput at high user mobility due to a better channel estimate and therefore, a lower frame error rate. Apart from that the authors present simulation results about the impact on the system performance of usage of variable number of feedback symbols from the decoder to the channel estimation block. These results show interesting trade-offs between the computational complexity and the algorithmic performance for the TD-SAGE algorithm. Therefore, the number of feedback symbols seems to be interesting parameter for a trade-off analysis between energy dissipation and algorithmic performance of a dedicated VLSI architecture.

Contributions: This chapter introduces an extension of the first ASIC implementation of the TD-SAGE algorithm which is presented in [START_REF] Minwegen | A multimode decision-directed channel estimation ASIC for MIMO-OFDM[END_REF]. The TD-SAGE algorithm is transformed to a novel variant, further reducing the computational complexity, to what is termed the tap alternating (TA) SAGE. Apart from that this work discusses two options to further reduce the computational complexity with the penalty of a loss in algorithmic performance. Either reducing the number of feedback symbols as discussed in [START_REF] Ketonen | Decision directed channel estimation for improving performance in LTE-A[END_REF] or reducing the frequency of updating the channel estimate. A suitable VLSI architecture is presented and area, timing and power numbers are provided. The support for a variable number of feedback symbols introduces a modification to the state machine that has an impact on the critical path. Therefore the implementation results for a architectures with and without variable feedback support are presented. This work evaluates the hardware costs of channel tracking in a MIMO-OFDM system.

Outline of the chapter: The remainder of the chapter is organized as follows. In section 2 the system setup is introduced. Section 3 presents the implemented, modified algorithm followed by Section 5 deriving the ASIC architecture and details about the processing units and the memory are given. Finally Section 6 discusses post-layout implementation results.

System Model

π -1 π SISO Detector SISO Decoder PACE DD-CE Mapper r[n] Ĥ r[n + 1] {ĉ}[n] t[n] b Fig. 2. Receiver model
The system considered in this chapter is a MIMO-OFDM system using N T = 2 transmit antennas and N R = 2 receive antennas, with N K = 512 subcarriers and a cyclic prefix length of N L = 32. First, information bits {b} are encoded using a convolutional encoder with the generator polynomial [133 o 171 o ]. These code bits {c} are interleaved by a random interleaver, mapped to complex symbols using a 4-, 16-or 64-QAM modulation and then multiplexed over N T spatial streams, each corresponding to a transmit antenna. Each symbol stream is expressed as a vector

t i [n] = [t i [n, 0], . . . , t i [n, N K -1]]
T ∈ C NK×1 , where t i [n, k] is the complex symbol at time n on sub-carrier k transmitted over the ith antenna. Each of the spatial FD streams is processed by an OFDM modulator which outputs the final TD vector s i transmitted by the ith transmit antenna.

The channel model used in this chapter is a frequency-selective Rayleigh fading channel with a power delay profile according to the typical urban COST259 model. It is time variant with the correlation according to Jake's model with a normalized Doppler frequency of f d = 1.4468 • 10 -5 , a sub-carrier spacing of 15 kHz, a user velocity v = 50 km/h and a carrier frequency f c = 2.4 GHz.

The frame structure of the system setup used throughout this chapter is shown in Fig. 1. The first OFDM symbol of a frame is a preamble following an orthogonal preamble scheme over the transmit antennas as proposed in [START_REF] Li | Simplified channel estimation for OFDM systems with multiple transmit antennas[END_REF]. The subsequent OFDM symbols are only consisting of data symbol vectors across all sub-carriers.

After the DFT processing of the received TD samples the N K × 1 vector r j is obtained at the jth receive antenna and can be written as:

r j [n] = NT-1 i=0 H i,j [n]t i [n] + w j [n] (1)
with

H i,j [n] = [H i,j [n, 0], . . . , H i,j [n, N K -1]]
T ∈ C NK×1 being the channel frequency response between the ith transmit antenna and the jth receive antenna. w j [n] is additive white complex Gaussian noise. For the sake of a shorter notation the time index n will be omitted in the remainder of the chapter.

The receiver model considered throughout this chapter is depicted in Fig. 2. Each received data symbol vector is iteratively processed by a soft-in soft-out (SISO) detector and a SISO decoder. A data symbol vector is defined as the vector over all receive antennas at the kth sub-carrier. Detection is performed by a max-log MAP SISO sphere detector (SD) with QR decomposition [START_REF] Studer | Soft-input soft-output sphere decoding[END_REF], while the channel decoder is a BCJR decoder providing soft information of the coded bits {c}. For the simulation results each processing block is executed twice, corresponding to one complete iteration between the detector and the decoder is performed.

The channel estimation provides the required estimate of the channel frequency response to the detector. As depicted in Fig. 2 the CE is split into two parts. First, the PACE processing block calculates an initial CIR estimate based on the preamble. This initial estimate is used to decode the second OFDM symbol of the frame. Second, the DD-CE block uses the decoder decisions of time n in order to provide an updated estimate at time n + U x for the detection of the next OFDM symbol. Thus, the third OFDM symbol is detected and decoded based on the updated channel estimate. The update frequency (U x ) can vary depending on the considered Doppler frequency.

There is an option to adjust the computational complexity of calculation of an update of the CIR. It was presented in [START_REF] Ketonen | Decision directed channel estimation for improving performance in LTE-A[END_REF]. The idea is to reduce the feedback from the decoder and therefore reduce the number of computations. Either the full feedback is used, meaning that all modulation symbols of the previous OFDM symbol are used to calculate the CIR estimate or only every second, third and so on symbol is used in the calculation.

The TA-SAGE Algorithm

The H i,j from (1) can be expressed as the DFT of the CIR: H i,j = DFT(h i,j ), where h i,j is the CIR between the ith transmit antenna and the jth receive antenna. Only the calculation of the CIR estimate ĥi,j = ( ĥi,j [0], ..., ĥi,j [N L -1]) T will be investigated in the remainder of this chapter.

The description of the TD-SAGE algorithm in [START_REF] Ylioinas | Avoiding matrix inversion in DD SAGE channel estimation in MIMO-OFDM with M-QAM[END_REF] is modified in this chapter to remove redundant calculations for an efficient VLSI implementation.

First, a permutation matrix P l is defined such that P l ti cyclically shifts the vector ti by l elements with ti containing the N K complex remapped decoder hard decisions in the frequency domain for the ith antenna. Second, the TD vectors ŝi and z i,l are defined as follows:

ŝi = F H NK ti ∈ C NK (2)
z i,l = P l ŝi ∈ C NK (3)
where F NK is the N K dimensional DFT matrix, s j is the vector with the timedomain samples of the decisions from the decoder that were sent via the ith transmit antenna. z i,l is the by l elements cyclically-shifted time-domain vector.

With these definitions the signal model from ( 1) can be reformulated in the time domain to

y j = F H NK r j = NT-1 i=0 NL-1 l=0 h j,i,l z i,l + wj ( 4
)
where wj is the transformed noise. The estimate of h j,i,l at SAGE iteration m is denoted as ĥ(m) j,i,l . A SAGE iteration is defined as the calculation of an update for one single tap of the CIR. Therefore, the SAGE iteration range is m = 1, . . . , N L N T N i , where N i specifies how often the algorithm iterates over the complete CIR. The initial estimate provided by the PACE is ĥ(0) j,i,l . The iterations are done for each receive antenna independently. The whole processing of the algorithm is done on the TD samples and can be split into four steps that have to be executed for each receive antenna. For these steps a new variable is introduced. The residual (m) j is the vector of the values that are remaining after subtracting the reconstruction of the observation given the current CIR estimate from the real observation y j , which is basically the current estimation error. Then, the steps of the modified SAGE algorithm are the following:

Step 1: Initialize the residual and calculate the norm for each transmit antenna vector:

(0) j = y j - NT-1 i=0 NL-1 m=0 ĥ(0) j,n,m z n,m (5)
||ŝ i || 2 = NK-1 k=0 Re{ŝ i [k]} 2 + Im{ŝ i [k]} 2 (6)
Step 2: Select the current tap by

l = (m -1) mod N L (7) i = m -1 N L mod N T . (8)
Step 3: Calculate a new δ based on the decoder decisions and the previous .

δ (m) = z H i,l (m-1) ||ŝ i || 2 (9)
Step 4: Update the selected tap and the residual by ĥ(m) j,i,l = ĥ(m-1)

j,i,l + δ (m) ( 10
) (m) j = (m-1) j -δ (m) z i,l . (11)
Step 2 to 4 are repeated N i times for every tap of the estimated CIR.

A variable feedback from the decoder was evaluated in [START_REF] Ketonen | Decision directed channel estimation for improving performance in LTE-A[END_REF]. There the authors explained that they use for a reduced feedback either every second, fourth symbol and so on. In the time domain this means that only the first N K F = N K /F b time-domain samples are used to calculate an update of the CIR estimate. The variable F b is defined to be a power of two, where F b greater 8 is not considered. This can be directly realized using the formulars ( 5) - [START_REF] Meinerzhagen | Towards generic low-power area-efficient standard cell based memory architectures[END_REF] with the first elements of the vectors z, s i and j respectively. Fig. 3 shows the block error rate (BLER) for the investigated modulation schemes 4-,16-and 64-QAM using a floating-point implementation. A block is defined as one code word which is spread over one OFDM symbol. The simulations for 4-, 16-and 64-QAM were performed with N i = 3 iterations. The number of estimated taps is N L = 32 which equals the length of the cyclic prefix. This is a worst case assumption for the presented OFDM system which is used throughout this work. Besides the floating-point simulations the results for a fixed-point implementation are shown in Fig. 3. The degradation due to the fixed-point arithmetic is negligible.

Algorithm Evaluation

For the follwing evaluation the modulation 4-QAM was chosen exemplary and the number of internal iteration is N i = 1. In Fig. 1 the update frequency is given by U x . In the following two different options to reduce the computational complexity by a factor of four are evaluated for two different exemplary operation points. The plot depicted in Fig. 4 shows the algorithmic performance for 4-QAM and a mobile device speed of v = 50 km/h and an U x = 1 and U x = 4. The loss in algorithmic performance is about 1 dB at a BLER of 1 %. Additionally the BLER for a reduced feedback from the decoder is plotted. There F b equals 4 which means that every 4th symbol from decoder is used to calculate the update of the CIR estimate. In this case leads to the same computational complexity than using the full feedback but updating the CIR estimate only every 4th OFDM symbol. From an algorithmic perspective it can be concluded from Fig. 4 that for the given mobile speed and the same computational complexity it is a better choice to update the CIR estimate every 4th OFDM symbol using the full decoder feedback.

Fig. 5 shows a BLER plot with the same parameters of U x and F b but evaluated at a mobile device speed of v = 100 km/h. The loss in terms of algorithmic performance is in this case for an U x = 4 about 7 dB at a BLER of 10 % using the full decoder feedback, compared to updating the CIR estimate every OFDM symbol. Apart from that it also shows a clear error floor at a BLER of 3 • 10 -2 . The second option using only one 4th of the decoder feedback for the calculations shows at v = 100 km/h a loss of 2 dB at a BLER of 10 %. The plot depicted in Fig. 5 also shows that at the given speed it is better to use a fourth of the decoder feedback every OFDM symbol instead, when the computational complexity should be kept constant.

SNR [dB] Block Error Rate

Ux = 1, F b = 1 Ux = 4, F b = 1 Ux = 1, F b = 4

TA-SAGE VLSI Architecture

The architecture is split into three processing units and three different memories. The processing units are the residual update (RU) unit, the scalar product (SP) unit and the tap unit. The memories are the residual memory, the TX memory which stores the re-modulated symbol decisions from the decoder and the tap memory used to store all taps of the channel impulse response h. Fig. 6 shows these blocks and depicts the memory accesses from each block.

Processing Schedule

The processing is split into four different phases: load, pre-computation, iteration and write-back. The load and the write-back phases do not include any computation but are necessary to load the input data in the memories and write back the results. These phases are considered for completeness of the hardware complexity analysis. In the load phase the received data y j , the current CIR estimate ĥ and the decoder feedback ŝi are loaded into the memories depicted in Fig. 6. The first processing phase (pre-computation) corresponds to step 1 of the algorithm description. First, the scalar ||ŝ i || 2 is calculated for all transmit antennas in the SP unit. Second, the residual vector (0) j for all receive antennas is calculated in the SP and RU unit. Both units are running concurrently, processing different receive antennas. In parallel to the (0) j calculation the reciprocal of ||ŝ i || 2 is pre-computed for all transmit antennas, since it does not change over the internal iterations.

The second processing phase is the iteration phase, which corresponds to steps 2, 3 and 4 of the algorithm. Step 2 is reflected in the dedicated address generation of each memory. Step 3 is executed by the SP unit calculating the inner product of ( 9) and the multiplication with the scaling factor 1 ||ŝi|| 2 . The last steps of the algorithm are [START_REF] Simko | Implementation aspects of channel estimation for 3GPP LTE terminals[END_REF], executed on the RU unit and (10) calculated by the tap unit. To achieve full utilization of the processing units and account for the data dependencies between ( 9) and ( 11) the SP and RU unit are separated by a pipeline register and execute the calculations concurrently for different receive antennas. This is possible since there is no data dependency between different receive antennas, which is a property of the SAGE algorithm.

In the write-back phase the new calculated estimate of the CIR is written from the tap memory to the output ports.

Processing Units

The processing units are the parts of the architecture that are executing the calculations of ( 5) to [START_REF] Simko | Implementation aspects of channel estimation for 3GPP LTE terminals[END_REF]. Apart from the algorithmic parameters defined in the previous sections, the main architectural design parameter is the data path parallelism w.

SP Unit Section 5.1 discussed that the SP unit is used in two phases and calculates ( 5), ( 6) and [START_REF] Studer | Soft-input soft-output sphere decoding[END_REF]. It can be seen from ( 5) that all complex multiplications can be executed in parallel. Therefore, it is possible to have a data path parallelism up to N K . In ( 6) and ( 9) it is necessary to accumulate the result of the concurrent calculations. This is implemented via an adder tree. Due to the high data path parallelism (up to w = 32) the maximum achievable frequency is determined by these adder trees. This leads to the design decision to have a dedicated pipeline stage as shown in Fig. 7 (third pipeline stage). The separation into the first and second pipeline stage is done to avoid two real multipliers in chain. This unit includes 6 • w multipliers and 3 • log(w) + 7 • w adders. The multipliers in the first pipeline stage are active in the pre-computation phase and the iteration phase. The dotted part in Fig. 7 is only active during the pre-computation phase, where first ||ŝ i || 2 is computed and written into a register file and then the initial residual vector (0) j is computed while the sequential divider concurrently outputs all 1 ||ŝi|| 2 . The dashed part is active in the iteration phase calculating δ (m) .

{z i,l } { (m) } {z i,l } { (m) } {z i,l } { (m) } {z i,l } { (m) } Seq. Divider ||ŝi|| 2 1 ||ŝ i || 2 {δ} {δ} { (0) r } { (0) r } NT w Fig. 7. SP Unit

RU Unit

The RU unit also computes [START_REF] Ylioinas | Avoiding matrix inversion in DD SAGE channel estimation in MIMO-OFDM with M-QAM[END_REF]. Thus, the structure of this unit differs only slightly from the SP unit. The calculation of (11) allows for parallel complex multiplications up to N K . The separation into three pipeline stages for this unit is done to achieve a balanced design. The output registers are added to ensure the same latency for the RU and SP unit, which eases the scheduling of the memory accesses. The complexity in terms of multipliers and adders of the RU unit is 4 • w multipliers and 4 • w adders.

Tap Update Unit This unit updates the current tap [START_REF] Meinerzhagen | Towards generic low-power area-efficient standard cell based memory architectures[END_REF]. Due to the low requirements in terms of throughput and the low complexity (2 adders) of this unit compared to the RU and SP units it is no longer discussed separately.

Memory Architecture

As shown in Fig. 6 the design has three different memories. Each of these memories has a dedicated controller that includes an address generation unit and multiplexers to realize the different data access schemes. The first memory is the tap memory, which stores the N R N T N L taps of the CIR. This memory has the most relaxed constraints in the architecture. During the initialization phase it is read in every cycle from the RU and the SP unit with a linear addressing scheme. In the iteration phase the tap memory is read and written once per tap update, i.e. every NK w cycles. Therefore, one read/write port is sufficient. The TX memory stores the N T N K TD samples of the complex symbols of the remapped decisions from the decoder. The circular shift in (3) is realized

{z i,l } {δ (m) } {z i,l } {δ (m) } {z i,l } {δ (m) } {z i,l } {δ (m) } { (m-1) } { (m-1) } w { (m) } { (m) } Fig. 8. RU Unit
as part of the address calculations. This memory is read by the RU and SP unit during the pre-computation phase and the iteration phase every cycle in parallel and written in the load phase. Each access reads/writes w elements in parallel. Therefore, two read/write ports with a word width of w elements are implemented.

The third memory is the residual memory. It stores the N R N K -values and needs to be read by the SP unit [START_REF] Studer | Soft-input soft-output sphere decoding[END_REF] and read and written by the RU unit [START_REF] Simko | Implementation aspects of channel estimation for 3GPP LTE terminals[END_REF] independently and concurrently with a word width of w elements. Furthermore, during the pre-computation phase the SP and RU units read and write independently the residual memory [START_REF] Ylioinas | Avoiding matrix inversion in DD SAGE channel estimation in MIMO-OFDM with M-QAM[END_REF]. Thus, the residual memory has two read and two write ports with a data width matching the data parallelism w.

With w and the algorithmic parameters N T , N R , N L , N K and N i , the cycle count of each phase can be calculated using the following equations.

Load phase:

c load = N T N R N L (12)
Pre-computation phase:

c precomp = N T N K (1 + N L )/w + 2 (13)
Iteration phase:

c iter = N K w N T N R N L N i + 3 + 2 (14)
Write-back phase:

c wb = N T N R N L + 2 (15)
The additive constants are needed because the pipeline of the processing units needs to be empty before the next phase can start. In the pre-computation phase the RU and the SP unit are running completely independently, leading to an overhead of two cycles. In the iteration phase the complete pipeline is the concatenation of the SP and RU unit. Thus the latency is 5. The latency of two in the write-back phase is due to the pipelined access to the memory. The total cycle count for one update of the CIR is the sum of the cycle counts of the four phases.

c total = c load + c precomp + c iter + c wb (16)
6 Implementation Results

The architecture was synthesized using a 90 nm, 1.0 V standard-performance standard cell library with Synopsys Design Compiler 2010.12-SP2 and layouted with Cadence SoC Encounter 9.1. In the following first the implementation results for the full feedback architecture are presented assuming Three different configurations of the architecture were implemented. A configuration is defined by its data path parallelism w = {8, 16, 32}. Each configuration was synthesized and layouted for its maximum achievable frequency and additionally for 400 MHz and 200 MHz. There are only two different design points for w = 32 since the maximum achievable frequency is 400 MHz.

N i = 3, N T = 2, N R = 2, N L =
The area-time trade-off diagram for the architecture variants is shown in Fig. 9. In this diagram T exec is defined as the time that the specific architecture requires to calculate a complete update of the CIR.

The best AT exec = 53.55 mm 2 µs product is the configuration with w = 32 and a synthesis and layout constraint of 400 MHz. However, the following discussion will focus on the configurations with an execution time around 70 µs. The configuration with data path parallelism w = 8 @ 400 MHz has the AT exec = 122.81 mm 2 µs product. Doubling the parallelism w = 16 and halving the frequency (AT exec = 138.9 mm 2 µs), leads to the same execution time and only a slight increase in terms of area. This stems from the fact that this architecture is memory dominated while an increase in the data path parallelism does not influence the memory as much as the data path (Table 1).

The memories in the presented architecture are implemented using standard cell based memories (SCM) [START_REF] Meinerzhagen | Towards generic low-power area-efficient standard cell based memory architectures[END_REF]. In this work flip-flop SCMs are used. The TX memory needs to be split into w banks each providing one word to allow for non-aligned vector accesses. This would lead to 32 macro cells for the maximum configuration, rendering floor-planning difficult. Therefore, the SCMs were utilized in this architecture.

Besides area and timing analysis, post-layout simulations were performed to obtain power estimates for the different configurations. The post-layout simulations with timing annotations were executed for independent test vectors for each configuration in order to obtain statistic toggling information. Synopsys Power Compiler uses the post-layout netlist and the annotated toggling information to calculate the average power estimates.

The results of the power analysis are shown in Fig. 10. The execution time T exec is the same as in the AT plot depicted in Fig. 9. The power-based analysis leads to different conclusions than the AT -based analysis. The energy for a certain configuration to compute one CIR estimate is the P T exec product. Comparing the configuration w = 8 @ 400 MHz with w = 16 @ 200 MHz in terms of smallest AT product leads to choose the configuration w = 8 @ 400 MHz. The same comparison achieving the lowest energy dissipation leads to the choice of the configuration w = 16 @ 200 MHz. The main reason that a doubled data path parallelism and a halved clock frequency leads to a better design decision in terms of energy dissipation is the fact that this architecture is memory dominated. This means that halving the frequency saves more energy in the memories than it is added by doubling the data path parallelism, because the memory including the memory controller is not affected as much as the processing units by an increase of the data path.

The area requirements of the two aforementioned design points are split into the different parts of the architecture in Table 1. As expected the tap memory is the smallest one since it only has to store N L N T N R elements. From this breakdown the aforementioned memory domination of this architecture is obvious. This is also supported by the layout for the configuration with w = 8 @ 400 MHz in Fig. 11. The layout shows that the SCM based memory approach makes it possible that the memory is placed where needed. In section 5.3 the discussion of the residual memory revealed that the SP and the RU unit are either writing it in parallel with w parallel accesses in the initialization phase or reading and writing it during the iteration phase. It can be seen that due to this constraints the RU unit is surrounded by the residual memory and the SP unit is placed close to it.

Flexible Feedback Architecture

The extension of the architecture to support the flexible feedback involved adjustements in the state machine and therefore in the schedule. Fig. 12 compares the implementation results for the full and the flexible feedback architectures.

It can be seen that for 400 MHz the area increase becomes visible in the diagram. This comes from the fact, that the modifications have an influence on the critical path and therefore to achieve the given timing constraint more area has to be invested. The average increase in area is up to 10 %. The same can be observed for the power dissipation depicted in Tab. 2. The procentual increase in the power is the same as the area increase for this architecture configuration. Therefore, the additional power only comes from the fact that the extension for the flexible feedback support has an influence on the critical path.

Conclusion

In this chapter we present to the best of our knowledge the first ASIC implementation of a decision directed channel estimation for MIMO-OFDM for M-QAM. The architecture is described and formulas for the calculation of the run-time of the algorithm depending on its parameters on the architecture are presented. The implementation is characterized in terms of area-time trade-offs and power dissipation.

It was shown, that the additional hardware costs of a channel tracking algorithm like the TA-SAGE are high compared to traditional PACE as presented in [START_REF] Simko | Implementation aspects of channel estimation for 3GPP LTE terminals[END_REF] but it is possible and therefore worth further investigations. Future work will include the influence of using latched based SCMs as presented in [START_REF] Meinerzhagen | Towards generic low-power area-efficient standard cell based memory architectures[END_REF] and the evaluation of mixing macro cell memories (e.g. for the residual and the tap memory) with the SCMs approach for the TX memory. Furthermore, this architecture will be compared with simplified algorithms as for example in [START_REF] Qiao | Decision-directed channel estimation for MIMO-OFDM systems[END_REF].

Fig. 1 .

 1 Fig. 1. The frame structure of the system

Fig. 3 .

 3 Fig. 3. Block error rate for the TA-SAGE with 4-,16-and 64-QAM. Using all available time-domain samples (NK )

Fig. 4 .Ux = 1 , 4 Fig. 5 .

 4145 Fig. 4. Block error rate for the TA-SAGE with 4-QAM at a speed of v = 50 km/h

Fig. 6 .

 6 Fig. 6. High level architecture

Fig. 9 .

 9 Fig. 9. Area-time trade-offs for different degrees of parallelism w and different synthesis/layout constraints. The algorithmic parameters are Ni = 3, NT = 2, NR = 2, NL = 32 and NK = 512.

Fig. 10 .

 10 Fig. 10. Power-time trade-offs for different configurations of the TA-SAGE and different synthesis/layout constraints. The algorithmic parameters are Ni = 3, NT = 2, NR = 2, NL = 32 and NK = 512.

Fig. 11 .

 11 Fig. 11. Layout of the TA-SAGE ASIC for a parallelism degree of w = 8.

Fig. 12 .

 12 Fig. 12. Area-time trade-offs for different degrees of parallelism w and different synthesis/layout constraints. The algorithmic parameters are Ni = 1, NT = 2, NR = 2, NL = 32 and NK = 512.

Table 1 .

 1 Area breakdown for the TA-SAGE w = 16 @ 200 MHz and w = 8 @ 400 MHz area w = 16 [µm 2 ] area w = 8 [µm 2 ]

	Tap mem	74284.8 (3.98%)	76822.6 (4.6%)
	TX mem	515499.5 (27.62%) 518432.5 (31.08%)
	Resid mem 747217.1 (40.03%) 767431.7 (46.01%)
	RU unit	154939.6 (8.3%)	90029.9 (5.4%)
	SP unit	374723.8 (20.07%) 215378.1 (12.91%)

Table 2 .

 2 Power comparision of the full and flexible feedback architecture with different frequencies and a w = 8 Frequency TA-SAGE TA-SAGE flex

	100 MHz	77.63 mW	82.84 mW
	200 MHz	152 mW	162.24 mW
	400 MHz	374 mW	400.1 mW

Acknowledgments This work has been supported by the UMIC Research Centre, RWTH Aachen University. The authors would like to thank Ernst Martin Witte, David Kammler, Martin Senst, Filippo Borlenghi and Uwe Deidersen for the valuable discussions and their feedback.