Elian Dib

Mikaël Le Pendu

Christine Guillemot

LIGHT FIELD COMPRESSION USING FOURIER DISPARITY LAYERS

In this paper, we present a compression method for light fields based on the Fourier Disparity Layer representation. This light field representation consists in a set of layers that can be efficiently constructed in the Fourier domain from a sparse set of views, and then used to reconstruct intermediate viewpoints without requiring a disparity map. In the proposed compression scheme, a subset of light field views is encoded first and used to construct a Fourier Disparity Layer model from which a second subset of views is predicted. After encoding and decoding the residual of those predicted views, a larger set of decoded views is available, allowing us to refine the layer model in order to predict the next views with increased accuracy. The procedure is repeated until the complete set of light field views is encoded. Following this principle, we investigate in the paper different scanning orders of the light field views and analyse their respective efficiencies regarding the compression performance.

INTRODUCTION

Light field imaging is gaining in popularity for a variety of vision applications, thanks to the emergence of real light fields capturing devices, and commercially available cameras. However, light fields represent very large volumes of high dimensional data, hence the need for designing efficient compression algorithms. Many solutions proposed so far adapt standardized image and video compression solutions (in particular HEVC) to light field data (e.g. [START_REF] Liu | Pseudo-sequence-based light field image compression[END_REF] [2] [START_REF] Li | Compression of unfocused plenoptic images using a displacement intra prediction[END_REF]). The authors in [START_REF] Jiang | Light field compression with homography-based low-rank approximation[END_REF], [START_REF] Dib | Super-ray based low rank approximation for light field compression[END_REF] investigate the use of homography-based low rank models for reducing the angular dimension, while local Gaussian mixture models in the 4D ray space are considered in [START_REF] Verhack | Steered mixture-of-experts for light field coding, depth estimation, and processing[END_REF]. A depth based segmentation of the light field into 4D spatio-angular blocks is used in [START_REF] Tabus | Lossy compression of lenslet images from plenoptic cameras combining sparse predictive coding and jpeg 2000[END_REF] for prediction, and the prediction residue is encoded using JPEG-2000.

In this paper, we describe a novel compression algorithm based on the Fourier Disparity Layer representation (FDL) introduced in [START_REF] Pendu | A fourier disparity layer representation for light fields[END_REF]. This representation has been shown to be efficient for a variety of processing applications, namely rendering, view synthesis and denoising. The proposed compression approach iteratively reconstructs the light field using the FDL representation. More precisely, an initial subset of views is first encoded with traditional video compression methods, and is used to construct a FDL representation of the light field. The FDL representation is then used to synthesize a second subset of light field views (neighboring ones) according to a pre-defined order. Given that there may remain correlations between the prediction residue of these synthesized views, the residual signal is best encoded using a video encoder with inter coding (in the experiments we used HEVC-inter, and in particular the HM 16.10). Then, in order to predict and encode the next subset of views, a more accurate FDL representation is constructed from the previously encoded and decoded subsets. The FDL representation is thus iteratively refined after the encoding of each view subset until all the views are encoded. Two types of scanning orders are considered for synthesis and coding, one following a hierarchical scheme and the other one following a circular order. In both cases, prediction schemes with either two or four view subsets were tested. For the best performing schemes, experimental results with real and synthetic light fields show average rate savings of more than 50%, using the Bjontegaard measure, with respect to the JPEG-Pleno verification model (VM 1.1) [START_REF]Jpeg pleno light field coding vm 1.1[END_REF].

NOTATIONS AND SCHEME OVERVIEW

Let us consider the 4D representation of light fields proposed in [START_REF] Levoy | Light field rendering[END_REF] and [START_REF] Gortler | The lumigraph[END_REF] describing the radiance along rays by a function Lpx, y, u, vq. This representation is based on a parameterization of orientations of light rays with two parallel planes with pairs px, yq and pu, vq respectively representing spatial and angular coordinates of light rays.

The coding algorithm is depicted in Fig. 1. First, the Fourier Disparity Layer calibration described in [START_REF] Pendu | A fourier disparity layer representation for light fields[END_REF] determines a set of disparity values as well as the angular coordinates of each view. These parameters are needed for the FDL construction and view prediction steps, and are transmitted as metadata to the decoder. Our coding scheme partitions the light field into subsets of views, the first subset being directly encoded as a group of pictures using a video encoder. This first set of views is used for an initial construction of the FDL representation which allows us to synthesize (or predict) the views of the second subset. The prediction residue is coded, decoded and added to the prediction to reconstruct the corresponding views. The views of the two first subsets are then used to re-compute the FDL representation that will then be used for synthesizing the views of the third subset. The algorithm continues iterating until all the light field views have been coded.

FOURIER DISPARITY LAYERS (FDL)

For simplicity of notation, let us consider only one 2D slice of the light field with only one spatial and one angular dimension. A view Lu 0 of the light field at angular coordinate u 0 is then defined by Lu 0 pxq " Lpx, u 0 q. It has been shown in [START_REF] Pendu | A fourier disparity layer representation for light fields[END_REF] that, given a set of n disparity values td k u kPv1,nw , the Fourier Transform Lu0 of Lu 0 can be decomposed as:

Lu0 pωxq " ÿ k e `2iπu0d k ωx Lk pωxq. (1)
where ωx denotes the spatial frequency, and Lk is defined by

Lk pωxq " ż Ω k e ´2iπxωx Lpx, 0qdx. (2)
Each function Lk can be interpreted as the Fourier transform of the central view only considering a spatial region Ω k of disparity d k , hence the name Fourier Disparity Layers (FDL). Given m input views Lu j (j P v1, mw) and by computing their Fourier Transforms Luj , the FDL representation can be learned by solving a linear regression problem for each frequency coefficient ωx. Constructing the FDL amounts to solving the equation Ax " b, with a Tikhonov regularization, where A, x, b are matrices and vectors of dimensions m ˆn, n ˆ1 and m ˆ1 respectively. The matrix A is defined as A jk " e `2iπuj d k ωx , while the vector x containing the Fourier coefficients of the disparity layers (for the frequency ωx) is defined as x k " Lk pωxq, and b containing the Fourier coefficients of the input image j is defined as b j " Luj pωxq.

VIEW SYNTHESIS USING FDL

Knowing the layers and their disparities d k , any view Lu 0 , at angular position u 0 , can be synthesized in the Fourier domain using Eq.( 1), and by computing the inverse Fourier transform. One important issue that can have a strong impact on the quality of the synthesized views used for prediction, is the view synthesis (or prediction) order. Typically, views that are close to the input views of the FDL construction step are better synthesized than more distant views. Furthermore, it was observed in [START_REF] Pendu | A fourier disparity layer representation for light fields[END_REF] that constructing a FDL from all the outer views at the periphery of the light field is an ideal configuration for synthesizing the inner views.

The circular prediction order illustrated in Fig. 2(a) was designed following those insights. To give an example with a light field of dimension 9 ˆ9, the circular scheme starts by encoding the views in the middle of each segment forming the periphery (in dark blue in Fig. 2(a)). The four adjacent views (light blue) are then synthesized using the FDL constructed from the base views. The prediction residue for these 16 views is then coded with HEVC inter, decoded and added to the synthesised views. The reconstructed views at the corresponding positions will then be used together with the base views to refine the FDL representation. The algorithm iterates considering next the 12 views "closing" the circle (green in Fig. 2(a)). Note that we preferred forming a circle rather than a square in order to avoid using the corner views for predicting the inner views. This choice is motivated by the fact that, in real light field datasets captured by plenoptic cameras, the corner views are generally of lower quality. A simplified circular prediction scheme was also designed as shown in Fig. 2(c), where only two subsets are used, and all the views forming the circle are contained in the initial subset.

For the comparison, we have also designed and tested more traditional hierarchical configurations as shown in Fig. 2(b) and (d) respectively with four and two subsets. The prediction schemes in Fig. 2(a),(b),(c), and(d) are referred to as Circular-4, Hierarchical-4, Circular-2, and Hierarchical-2 respectively.

CODING SCHEME

FDL Calibration

To construct the FDL, one first needs to estimate precise angular coordinates u j of input views as well as the layers' disparity values d k . The joint estimation of these two sets of parameters allows the method to be robust to common issues with real light field data, e.g. to the fact that the angular coordinates of all the views may not form precisely a square grid as it is generally assumed. These parameters u j and d k are found by minimizing over the Q frequency components ω q

x (Q being the number of pixels in each input image) [START_REF] Pendu | A fourier disparity layer representation for light fields[END_REF] min

x q ,u,d Q ÿ q"1 ´ Apω q x , u, dqx q ´bq 2 2 `λ Γx q 2 2 ¯, (3)
where the input view positions u j and the disparity values d k are arranged in the vectors u and d respectively, and where rApω q x , u, dqs j,k " e `2iπuj d k ωx . The vectors x q and b q contain the Fourier coefficients of, respectively, the disparity layers and the input images at the frequency ω q x (i.e. x q k " Lk pω q x q and b q j " Bj pω q x q). The regularization matrix Γ is defined as a discrete approximation of the second order differential operator:

Γ " » - - - - - - ´2 1 1 ´2 1 . . . . . . . . . 1 ´2 1 1 ´2fi ffi ffi ffi ffi ffi fl . ( 4
)
Note that the calibration results depend on the regularization parameter λ. We optimize this parameter by performing several calibrations with different values of λ. Given a prediction order in our coding scheme (see Fig. 2), we keep the calibration results yielding the best predictions in the least square sense. For this step, the prediction of a subset of views does not take into account the compression loss of the previous subsets, in order to avoid performing the complete compression algorithm only to determine λ. Therefore, the disparity values d k and view positions u j are learnt for optimal predictions within our coding scheme.

Prediction Scaling

Additionally, in order to cope with possible variations of illumination between views, we perform a scaling of each predicted view. Let the image I be one of the views, and Ipr its prediction obtained using the FDL view synthesis. The final prediction Ă Ipr is then computed as Ă Ipr " s ¨Ipr, where the scaling factor s is derived by the encoder to minimize the sum of squared error of the scaled prediction as follows:

arg min s I ´s ¨Ipr 2 2 " xI, Ipry Ipr 2 2 , (5)
where x¨, ¨y is the inner product. The value of s is transmitted as metadata in order to perform the same prediction scaling on the decoder side.

Encoding of metadata

All the transmitted parameters (i.e. scaling factors, set of disparity values and angular coordinates), are encoded in the double precision floating point format with a fixed length encoding, resulting in 64 bits for each parameter. The additional cost is accounted for in the results presented in Section 6.

Encoding of initial views and view residuals

The initial set of views is directly encoded as a video sequence using HEVC 8 bits with inter coding. For the other views, the prediction residual is also encoded in HEVC. Since the prediction residual requires one additional bit of precision compared to the original signal, we avoid the risk of precision loss by using HEVC 10 bits. Furthermore, due to imperfections either in the original signal or in the FDL view synthesis, correlations may remain between the prediction residual of the different predicted views. Therefore, at each iteration, the current subset of views to encode is arranged in a sequence and encoded using HEVC-inter.

For both the circular and the hierarchical prediction schemes, the views are arranged in a spiral order starting from the center of the light field, and only considering the views of the current subset.

In order to optimize the bitrate allocation, we use different QP parameters in the HEVC encoding of the different subsets of views. Let QP 1 be the base parameter used for the initial subset of views, then we use QP t " QP t´1 `1 as the parameter for the subset of index t in the coding order.

EXPERIMENTAL RESULTS

Experimental Setting

Experiments were performed using the luminance component of light fields coming from the HCI [START_REF] Wanner | Datasets and benchmarks for densely sampled 4D light fields[END_REF], INRIA [START_REF]INRIA Lytro image dataset[END_REF] and ICME 2016 Grand Challenge [START_REF]ICME 2016 Grand Challenge dataset[END_REF] datasets. The HCI dataset contains synthetic light fields with 9x9 views of 768x768 pixels (Buddha, Butterfly, StillLife) and with 9x9 views of 512x512 pixels (Greek, Sideboard), the INRIA and ICME datasets originally contain light fields captured by a Lytro Illum camera with 15x15 views of 625x434 pixels from which we keep the 9x9 central views (to avoid strong vignetting issues) of 616x424 pixels (to remove the black borders and ensure the size is a multiple of 8 for the HEVC coding). The Lytro LFs have been extracted using the Matlab Light Field Toolbox v0.4 [START_REF] Dansereau | Light Field Toolbox for Matlab[END_REF] with gamma correction.

The performances of the proposed FDL coding schemes have been compared to JPEG Pleno VM 1.1 [START_REF]Jpeg pleno light field coding vm 1.1[END_REF] and to the HEVC coding of all the views arranged in a spiral order. The latter method can be seen as a particular case of our approach where all the views are in the initial set. It is therefore referred to as the Circular-1 scheme. For Circular-1 and our FDL coding schemes, the HEVC encoding is performed using the HM 16.10 with the Main-RExt profile and the random access configuration. We have used a range of quality parameters QP P t5, 10, 15, 20, 30, 40u (which corresponds to the QP parameter of the initial view subset in our approach). For our FDL coding schemes, we have fixed the number of layers in the FDL method to n " 30. Although the random access configuration was used for HEVC coding (intra period of 32), the random access capability of our approach is limited by the fact that the views of the last subset are dependant on those of the previous subsets. In practice, for the four variants in Fig. 2, the percentage of the bitstream required to decode one view in the worst case is about 80% at high bitrates (QP " 40) and about 95% at low bitrates (QP " 5).

For JPEG Pleno, the required input central disparity maps have been generated using [START_REF] Jiang | Depth estimation with occlusion handling from a sparse set of light field views[END_REF] and the configuration files provided with JPEG Pleno VM 1.1 have been adapted for use on 9x9 light fields.

Results

The rate-distortion results are presented in Fig. 3 for two real and two synthetic light fields. For all the tested light fields, the bitrate savings of the different variants of our coding scheme, as well as the direct HEVC encoding (i.e. 'C1') are presented in Table 1. They have been computed with the Bjontegaard metric [START_REF] Bjontegaard | Calculation of average PSNR differences between RD curves[END_REF] using JPEG-Pleno as a reference.

We observe in Table 1 that the Circular-2 (respectively Circular-4) scheme systematically outperform Hierarchical-2 (respectively Hierarchical-4). This result supports the idea that using a set of connected peripheral views to construct the FDL rather than a more scattered set substantially improves the view prediction. Additionally, the variants of our coding scheme with only two subsets (i.e. Circular-2, Hierarchical-2) have significantly better performances than the ones with 4 subsets (i.e. Circular-4, Hierarchical-4). This may be due to high frequency artifacts in the predicted views caused by the FDL construction from a too sparse set of known views. Compared to JPEG-Pleno, more than 50% rate gains are obtained on average for our variants with two subsets, This is a substantial improvement over the average 35% gain of the Circular-1 method which directly encodes all the views in HEVC with the same configuration. It can also be observed in Figs. 3 and4 that particularly large PSNR gains are obtained for very low bitrates, especially compared to JPEG-Pleno. This can be explained by the very limited amount of additional data needed to perform our predictions. On the other hand, JPEG-Pleno requires the transmission of a complete disparity map, which results in a significant overhead for low bitrate coding.

Note finally that another advantage of our approach is to allow scalable light field coding. Since the FDL model can be used to synthesize any view of the light field, a complete light field can be reconstructed by the decoder at each iteration. By decoding additional view subsets, a more accurate light field is obtained. In that sense, despite their lower performance, the schemes with 4 iterations (i.e. Circular-4, Hierarchical-4) may be preferable in some scenarios as they provide additional levels of scalability. Furthermore, Circular-4 still clearly outperforms the JPEG-Pleno anchor in most cases (20% bitrate savings on average).

CONCLUSION

We have presented a novel light field compression method based on the Fourier Disparity Layer representation. In the proposed scheme, the light field is partitioned into several subsets of views, where the first subset is encoded as a video sequence using HEVC. The next subsets are iteratively predicted from the previously encoded and decoded ones thanks to the Fourier Disparity Layer view synthesis, and only the prediction residual is encoded with HEVC. We have shown that defining view subsets forming a circular pattern is better suited to this prediction scheme than the more conventional hierarchical coding order. In the tested configurations with only two subsets of views, significant bitrate savings were obtained in comparison to both JPEG-Pleno and a direct encoding of all the views with HEVC. Using more subsets with our circular coding order still outperforms JPEG-Pleno, but reduces the performance compared to the case of two subsets. However, it provides additional levels of scalability since the FDL model refined with each decoded subset can be used to synthesize any view of the light field.

  This work was supported in part by the EU H2020 Research and Innovation Programme under grant agreement No 694122 (ERC advanced grant CLIM) and in part by the Science Founda-tion Ireland (SFI) under the Grant Number 15/RP/2776.

Fig. 1 .

 1 Fig. 1. Overview of the encoding (a) and decoding (b) scheme.

Fig. 2 .

 2 Fig. 2. View subsets for different view prediction orders: (a) Circular-4, (b) Hierarchical-4, (c) Circular-2, (d) Hierarchical-2.

Fig. 3 .

 3 Fig. 3. PSNR-Rate performance on real (top) and synthetic (bottom) light fields from (a) INRIA, (b) ICME and (c,d) HCI datasets.

Fig. 4 .

 4 Fig. 4. Visual comparison of the reconstructed top left view image for JPEG Pleno VM 1.1, Circular-1 (C1) and Circular-2 (C2).

Table 1 .

 1 Bjontegaard percentage rate savings for the Hierarchical-2

	Method	H2	C2	H4	C4	C1
	Bench	-25.53	-34.3	2.76	-15.05 -10.29
	Fruits	-35.54 -46.02 -4.05 -22.12 -11.67
	Toys	-48.17	-51.2	-4.56 -16.24 -52.67
	Bikes	-41.18 -47.31	1	-15.82 -20.05
	D de M	-29.78 -40.42 -3.62 -21.61	1.54
	F & V 2	-45.89 -48.64 -8.18 -27.44 -28.37
	Friends 1	-46.41 -48.75 -21.5 -37.39 -21.41
	S P I	-28.25 -33.22 41.08	7.19	-41.06
	Vespa	-43.87 -47.55 21.07	-7.49	-23.45
	Greek	-58.44 -67.31 20.47	-2.14	-52.24
	Sideboard -56.31 -64.14 54.27 27.09 -13.63
	Buddha	-62.69 -68.21 -8.19 -26.97 -31.73
	Butterfly	-86.03 -86.96 -72.6 -75.16	-75.6
	StillLife	-49.8	-60.99 44.85 16.51 -60.11
	Average	-50.11	-56.5	-1.54 -20.42 -35.42
	(H2), Circular-2 (C2), Hierarchical-4 (H4), Circular-4 (C4), and
	Circular-1 (C1) schemes with respect to JPEG Pleno VM 1.1 (nega-
	tive values represent gains).


Mikael Le Pendu

Christine Guillemot

Aljosa Smolic

A Fourier Disparity Layer representation for Light Fields

Keywords: Light Fields, Fourier domain, rendering, refocusing, view interpolation, denoising

In this paper, we present a new Light Field representation for efficient Light Field processing and rendering called Fourier Disparity Layers (FDL). The proposed FDL representation samples the Light Field in the depth (or equivalently the disparity) dimension by decomposing the scene as a discrete sum of layers. The layers can be constructed from various types of Light Field inputs including a set of sub-aperture images, a focal stack, or even a combination of both. From our derivations in the Fourier domain, the layers are simply obtained by a regularized least square regression performed independently at each spatial frequency, which is efficiently parallelized in a GPU implementation. Our model is also used to derive a gradient descent based calibration step that estimates the input view positions and an optimal set of disparity values required for the layer construction. Once the layers are known, they can be simply shifted and filtered to produce different viewpoints of the scene while controlling the focus and simulating a camera aperture of arbitrary shape and size. Our implementation in the Fourier domain allows real time Light Field rendering. Finally, direct applications such as view interpolation or extrapolation and denoising are presented and evaluated.

I. INTRODUCTION

Light Fields are commonly represented as 4 dimensional functions with 2 spatial and 2 angular dimensions [START_REF] Levoy | Light field rendering[END_REF], [START_REF] Gortler | The lumigraph[END_REF]. They can be seen as 2D arrays of images (called sub-aperture images), each having an unlimited depth of field, and differing from their neighbour images only by a slight shift of the view angle. The sampling in the angular dimensions is key in Light Field imaging [START_REF] Chai | Plenoptic sampling[END_REF]. In particular, densely sampled Light Fields make it possible to directly render images with shallow depth of field while controlling the focus depth. Such rendering, often referred to as Light Field refocusing, does not require knowledge of the scene's geometry. It is usually performed either by shifting and averaging the sub-aperture images [START_REF] Ng | Light Field Photography with a Hand-Held Plenoptic Camera[END_REF] or by selecting a 2D slice in the 4D Fourier domain [START_REF] Ng | Fourier slice photography[END_REF]. However, a dense angular sampling comes at the expense of very high requirements in terms of capture, storage, processing power and memory. A too sparse angular sampling, on the other hand, does not allow for a smooth transition between viewpoints and causes angular aliasing in the refocused images, characterized by sharp structures in the out of focus regions. The importance of a dense angular sampling is clearly shown by the vast literature on viewpoint interpolation. Several approaches exist including depth image Fig. 1. Fourier Disparity Layer representation. For the visualization, the layers are shown in the spatial domain (i.e. after inverse Fourier Transform). The magnitude spectrum of each layer is also shown in the red boxes. Note that the FDL representation differs visually from a focal stack: the out-of-focus regions tend to disappear while the regions in focus have enhanced contrasts. based rendering techniques [START_REF] Georgiev | Spatio-angular resolution tradeoffs in integral photography[END_REF]- [START_REF] Penner | Soft 3d reconstruction for view synthesis[END_REF], deep learning methods either exploiting a depth map estimation [START_REF] Kalantari | Learning-based view synthesis for light field cameras[END_REF], [START_REF] Flynn | Deepstereo: Learning to predict new views from the worlds imagery[END_REF] or not [START_REF] Wu | Light field reconstruction using deep convolutional network on epi[END_REF]- [START_REF] Yeung | Fast light field reconstruction with deep coarse-to-fine modeling of spatial-angular clues[END_REF], and approaches leveraging sparsity priors of the Light Field data in a transformed domain [START_REF] Shi | Light field reconstruction using sparsity in the continuous fourier domain[END_REF]- [START_REF]Light field reconstruction using shearlet transform[END_REF]. However, although viewpoint interpolation greatly simplifies the capture of dense Light Fields, it also increases the amount of data to store and process for the final rendering application.

Alternatively, Light Fields can be represented as a focal stack, that is, a set of shallow depth of field images (e.g. photos taken with a wide aperture) with different focusing depths. This representation has the advantage of allowing an unlimited angular density with few images because the sampling is performed on the depth dimension instead of 2 angular dimensions. However, for rendering tasks such as simulating a different camera aperture size or shape, or a change of viewpoint, the common approach is to first convert the focal stack into the 4D representation. For instance, Levin and Durand [START_REF] Levin | Linear view synthesis using a dimensionality gap light field prior[END_REF] retrieve sub-aperture images by the deconvolution of shifted and averaged focal stack images. A similar deconvolution technique is used in [START_REF] Kodama | Efficient reconstruction of all-in-focus images through shifted pinholes from multi-focus images for dense light field synthesis and rendering[END_REF] to first synthesize the 4D Light Field from a focal stack in order to render images with arbitrary aperture shapes. More recent methods have also been proposed to reconstruct the 4D Light Field from a focal stack either in the spatial domain using depth from focus [START_REF] Mousnier | Partial light field tomographic reconstruction from a fixed-camera focal stack[END_REF], or via optimization in the Fourier domain [START_REF] Pérez | Lightfield recovery from its focal stack[END_REF], [START_REF] Alonso | Reconstruction of perspective shifts and refocusing of a three-dimensional scene from a multi-focus image stack[END_REF].

The main motivation of this work is to simplify the handling of Light Fields through the definition of a compact representation which, unlike focal stacks, can be directly used for any Light Field rendering task. In this aim, we propose the Fourier Disparity Layer (FDL) representation illustrated in Fig. 1. It can be easily constructed either from a 4D Light Field, a focal stack, or even a hybrid Light Field combining subaperture and wide aperture images with varying parameters (i.e. focusing depth, aperture, point of view). Given the 2D Fourier transform of each input image with their parameters (i.e. angular coordinates, aperture, focus), the FDL model is constructed using linear optimization. For each frequency component, a linear least squares problem is solved to determine the corresponding Fourier coefficients of the different layers, each layer being associated to a given disparity value. The layers can then be directly used for real time rendering. For instance, sub-aperture images are obtained by shifting the layers proportionally to their associated disparity value and by averaging them. This is directly implemented in the Fourier domain as a simple linear combination with frequencydependent coefficients. More general rendering with arbitrary point of view, aperture shape and size, and focusing depth is performed with the same computational complexity without the need to first reconstruct the 4D Light Field.

In the case where the input is a set of sub-aperture images, we propose a gradient descent based calibration method to determine their angular coordinates as well as the optimal set of disparity values. The formulation of the optimization problems for the calibration and the layer construction are closely related. Nevertheless, we define two regularization schemes with different properties to better suit each situation.

Additionally, we demonstrate the effectiveness of our approach for several direct applications. First, when the input is a sparse set of sub-aperture images, view interpolation and extrapolation is obtained by constructing the FDL representation and by rendering views at intermediate angular coordinates with an infinitely small aperture. In a second application, the same viewpoints as the input are rendered to produce a denoised result. For this use case, we present a possible extension of the model where the shift applied to each layer is not constrained to be proportional to the angular coordinates of the view to reconstruct. The relaxed model allows a more accurate representation of occlusions and non-Lambertian effects in the scene.

Since the computational complexity is a key aspect motivating the need for a new Light Field representation, our implementation makes efficient use of the GPU at every step of the processing chain (i.e. calibration, layer construction, rendering). The proposed algorithms are built upon simple linear algebra operations performed independently at each frequency component, which makes our approach particularly suitable for GPU parallelization.

In summary, the contributions are:

' Definition of the Fourier Disparity Layer representation and its construction from other Light Field representations (e.g. sub-aperture images, focal stack, combination of focal stack images and sub-aperture images). ' Calibration method jointly determining the input view positions and disparity values of the layers. ' Fast and advanced Light Field rendering from the FDL representation with simultaneous control over the viewpoint, aperture size, aperture shape and focusing depth. ' Analysis of other application scenarios: view interpolation and denoising.

II. RELATED WORK

Related Light Field representations have been used in the design of several Light Field displays [START_REF] Wetzstein | Layered 3D: Tomographic image synthesis for attenuation-based light field and high dynamic range displays[END_REF]- [START_REF] Takahashi | From focal stack to tensor light-field display[END_REF]. These displays reproduce the Light Field using a stack of light attenuating LCD layers placed in front of a backlight. Thanks to the distance separating the LCD panels in the display, the image perceived depends on the observer's position and is proportional to the product of the layers. This layer representation has similarities with the one presented in this paper, and it can be constructed either from the Light Field views [START_REF] Wetzstein | Layered 3D: Tomographic image synthesis for attenuation-based light field and high dynamic range displays[END_REF], [START_REF] Wetzstein | Tensor Displays: Compressive Light Field Synthesis using Multilayer Displays with Directional Backlighting[END_REF] or from a focal stack [START_REF] Takahashi | From focal stack to tensor light-field display[END_REF]. The main difference however, is that, because of the physics of the light attenuating LCD layers, the sub-aperture images of the Light Field are reconstructed as a product of the layers' pixels instead of a sum. Hardware limitations also impose constraints on the layer representation. For instance, the number of layers is generally small (e.g. 3 to 5), which is often insufficient to accurately represent the whole Light Field. Furthermore, the layers must only have positive pixel values in order to be displayed on the LCD panels. This constraint is not required in our model, which allows us to efficiently construct the layers in the Fourier domain.

Similarly to the FDL method proposed in this paper, Alonso et al. [START_REF] Alonso | Reconstruction of perspective shifts and refocusing of a three-dimensional scene from a multi-focus image stack[END_REF] construct layers by an optimization in the Fourier domain. However, their method is limited to a focal stack input. In this configuration, the problem is well conditioned because the input images already contain dense angular information and each constructed layer is associated to the focusing depth of one of the focal stack images. Hence, no regularization scheme was considered for this application. The method we propose is more generic as it can also construct the layers from sub-aperture images. Therefore, specific regularization strategies are studied, which allows us to address a much larger range of applications including calibration, view interpolation, denoising, etc.

Finally, the proposed FDL representation directly relates to the dimensionality gap Light Field prior described by Levin and Durand [START_REF] Levin | Linear view synthesis using a dimensionality gap light field prior[END_REF]. It states that the support of the Light Field data in the 4D Fourier domain is a 3D manifold which was later characterized as a hypercone in [START_REF] Dansereau | Light field image denoising using a linear 4d frequency-hyperfan allin-focus filter[END_REF]. By additionally considering the limited depth range of a scene, Dansereau et al. [START_REF] Dansereau | Light field image denoising using a linear 4d frequency-hyperfan allin-focus filter[END_REF] determined that the frequency-domain support of the Light Field forms a hyperfan. They define this shape as the intersection of the hypercone with a dual fan previously described in [START_REF] Dansereau | A 4-d dual-fan filter bank for depth filtering in light fields[END_REF]. In this paper, we derive the FDL representation from the dimensionality gap prior assuming a discrete set of depths instead of a continuous range. For the discrete depth case, we show formally in Section IV-A that this prior is itself derived from the assumption of a non-occluded Lambertian scene. This is a limitation for any method directly enforcing the dimensionality gap prior. For instance, as observed in [START_REF] Levin | Linear view synthesis using a dimensionality gap light field prior[END_REF], [START_REF] Dansereau | Light field image denoising using a linear 4d frequency-hyperfan allin-focus filter[END_REF], in the reconstructed sub-aperture images, occluding objects may appear transparent near the occlusion boundaries. However, semi-transparent objects and reflections on flat surfaces are accurately reproduced, which is particularly challenging for depth image based rendering methods. We also present in Section V-C a possible generalization of our layer model to allow a better representation of other non-Lambertian effects and occlusions.

III. LIGHT FIELD NOTATIONS

Let us first consider the 4D representation of Light Fields in [START_REF] Levoy | Light field rendering[END_REF] and Lumigraph in [START_REF] Gortler | The lumigraph[END_REF] parameterized with two parallel planes, as illustrated in Fig. 2. The 4D representation describes the radiance along rays by a function Lpx, y, u, vq where the pairs px, yq and pu, vq respectively represent spatial and angular coordinates. For simplicity of notation, we consider a 2D Light Field Lpx, uq with one spatial dimension x and one angular dimension u, but the generalization to a 4D Light Field Lpx, y, u, vq is straightforward.

In this paper, we use the notion of disparity instead of depth. Given the depth z 0 of the focal plane in Fig. 2, a depth z can be directly converted into a disparity d with d " z0´z z (i.e. objects at depth z 0 from the camera plane have zero disparity). Refocusing the Light Field then consists in defining a new Light Field L 1 px, uq " Lpx ´us, uq. The refocus parameter s is defined such that the regions of disparity d " s in the Light Field L, have a disparity equal to zero in the Light Field L 1 . A refocused image, noted B s , is then formed by refocusing the Light Field with parameter s and integrating the light rays over the angular dimension as described in [START_REF] Ng | Light Field Photography with a Hand-Held Plenoptic Camera[END_REF]:

B s pxq " ż R Lpx ´us, uqψpuqdu. (1)
where ψpuq represents the aperture of the imaging system. In the case where B s is captured by a camera with an aperture area A (typically a disk), then ψ is defined by: ψpuq "

" 1 if u P A 0 otherwise. (2)
In this paper, we use a more general version of Eq. ( 1) where the image can be observed at any position u 0 on the camera plane with any given aperture ψ :

B s u0 pxq " ż R Lpx ´us, u 0 `uqψpuqdu. (3)
Note that when the aperture area is infinitely small, the function ψ is equal to the dirac delta function δ. In this case, the image formed by Eq. ( 3) is the sub-aperture image noted L u0 and defined by L u0 pxq " Lpx, u 0 q.

The different notations used in the article are summarized in Table I. Transpose of the matrix A. A Complex conjugate of A (without transpose). A ˚Conjugate transpose of A (A ˚" A J ).

A j,k

Element of A on the j th row and k th column. A k k th column of matrix A.

x k k th element of vector x.

IV. FOURIER DISPARITY LAYER REPRESENTATION

A. Light Field Prior and FDL Representation

For the derivations, we assume that the scene is Lambertian, without occlusion, and can be divided into n spatial regions Ω k with constant disparity d k . Formally, this can be written: @k P v1, nw, @px, uq P Ω k ˆR, Lpx ´ud k , uq " Lpx, 0q. (4)

Here, the spatial regions Ω k are defined for the central view at u " 0. From this assumption, we prove in Appendix A that the Fourier transform of the Light Field can be decomposed as follows:

Lpω x , ω u q " ÿ k δpω u ´dk ω x q Lk pω x q, ( 5
)
where Lk is defined by

Lk pω x q " ż Ω k e ´2iπxωx Lpx, 0qdx. (6)
Each function Lk can be interpreted as the Fourier transform of the central view only considering the region Ω k of disparity d k . Hence, we call these functions Fourier Disparity Layers (FDL). One can note from Eq. ( 5) that the Light Field information is entirely contained in the Fourier Disparity Layers and their associated disparity values d k . We will show that the layers Lk can be constructed either from sub-aperture or wide aperture images without knowing the regions Ω k . In other words, our method does not require a disparity map. However, the disparity values d k are necessary. Hence, a method for estimating these values is also presented in Section V-B. Another interpretation of Eq. ( 5) is that the Light Field is sparse in the Fourier domain and the non-zero values are located at frequencies such that ω u " d k ω x for each disparity value d k of the Light Field. This is similar to the Light Field prior referred to as dimensionality gap prior in [START_REF] Levin | Linear view synthesis using a dimensionality gap light field prior[END_REF], and resulting in the definition of the hyperfan filter in [START_REF] Dansereau | Light field image denoising using a linear 4d frequency-hyperfan allin-focus filter[END_REF]. However, these previous works consider a continuous disparity range, which does not allow for a practical layer representation, since it would result in an infinite number of layers.

B. FDL Decomposition of Images

Sub-aperture images: Here, we derive the relationship between the FDL representation of the Light Field and the Fourier Transform Lu0 of a sub-aperture image L u0 . Note that in Lu0 , the Fourier Transform only applies to the spatial dimension. Hence, Lu0 is obtained from the Light Field's spectrum L by applying the inverse Fourier transform in the angular dimension as follows:

Lu0 pω x q " ż `8 ´8 e `2iπu0ωu Lpω x , ω u qdω u . (7)
Using the FDL decomposition of Lpω x , ω u q from Eq. ( 5), the transformed sub-aperture image can be directly derived:

Lu0 pω x q (8) " ÿ k Lk pω x q ż `8 ´8 e `2iπu0ωu δpω u ´dk ω x qdω u (9)
"

ÿ k e `2iπu0d k ωx Lk pω x q. ( 10
)
General image model: Now, considering the general model in Eq. ( 3) the Fourier Transform of an image B s u0 is:

Bs u0 pω x q " ż `8 ´8 e ´2iπxωx "ż `8 ´8 Lpx ´us, u 0 `uqψpuqdu  dx (11)
"

ż `8 ´8 ψpuq ż `8 ´8 e ´2iπpx`usqωx Lpx, u 0 `uqdxdu (12)
"

ż `8 ´8 e ´2iπusωx ψpuq Lu0`u pω x qdu. (13)
Using the FDL decomposition of Lu0`u from Eq. ( 10), we obtain:

Bs u0 pω x q " ż `8 ´8 e ´2iπusωx ψpuq ÿ k e `2iπpu0`uqd k ωx Lk pω x qdu (14)
"

ÿ k e `2iπu0d k ωx Lk pω x q ż `8 ´8 e `2iπupd k ´sqωx ψpuqdu (15)
"

ÿ k e `2iπu0d k ωx ψpω x ps ´dk qq ¨L k pω x q. ( 16
)
One can easily verify that in the particular case of an infinitely small aperture such that ψ " δ, we have ψ " 1, thus Eq. ( 16) becomes equivalent to the sub-aperture case in Eq. ( 10).

V. CONSTRUCTION OF THE LAYERS

The layer construction is performed from a set of m images noted tB j u jPv1,mw that follow the model of Eq. ( 3) (e.g. images captured by cameras on the same plane and oriented perpendicularly to this plane). Each image B j is associated with its angular coordinate u j , aperture function ψ j , and refocus parameter s j . Thanks to the genericity of this model, different forms of input data may be used. For example, in the case of a set of sub-aperture images, the aperture function ψ j of each image is equal to the dirac delta function δ. Note that, in this case, the refocus parameter s j has no influence in the equations (e.g. see Eq. ( 10)) and is not required. Another notable example is that of a focal stack input (i.e. images with the same angular coordinates u j and aperture functions ψ j , but with different refocus parameters s j ).

The layer construction problem is formulated in the general case in the next sub-section, assuming all the image parameters are known. Then, in sub-section V-B, we show how to determine the positions u j of the views and the layer's disparity values d k when the input data is a set of sub-aperture images.

A. Problem Formulation

The results of Eqs. ( 10) and [START_REF] Yeung | Fast light field reconstruction with deep coarse-to-fine modeling of spatial-angular clues[END_REF] show that for a fixed spatial frequency ω x in the Fourier domain, the images B j can be simply decomposed as linear combinations of the Fourier Disparity Layers. Therefore, the FDL representation of the Light Field can be learned by linear regression for every coefficient in the discrete Fourier domain.

We first compute the discrete Fourier Transform Bj of each image B j and we construct, for each frequency component ω x , a vector b such that b j " Bj pω x q. Given the disparity value d k of each layer (k P v1, nw), we also construct the matrix A P C mˆn with:

A jk " e `2iπuj d k ωx ψj pω x ps j ´dk qq. (17)
By defining the vector x with x k " Lk pω x q, the FDL decomposition in Eq. ( 16) is reformulated as Ax " b. A simple layer construction method then consists in solving an ordinary linear least squares problem independently for each frequency component ω x . In practice, however, this may be an ill-posed problem. It is typically the case when the number of available images is lower than the number of layers required to represent the scene (i.e. m ă n). More generally, depending on the input configuration, and for some frequency components, the matrix A may be ill-conditioned. This results in over-fitting and extreme noise amplification when the layer model is used to render new images of the scene (e.g. view interpolation or extrapolation). In order to avoid this situation, we include a Tikhonov regularization term in the problem formulation:

x " arg min x Ax ´b 2 2 `λ Γx 2 2 , ( 18
)
where Γ is the Tikhonov matrix, and λ is a parameter controlling the amount of regularization. The problem [START_REF] Takahashi | Good group sparsity prior for light field interpolation[END_REF], has a well-known closed form solution:

x " pA ˚A `λΓ ˚Γq ´1A ˚b, [START_REF] Vagharshakyan | Accelerated shearletdomain light field reconstruction[END_REF] where ˚is the Hermitian transpose operator (i.e. complex conjugate of the transposed matrix).

In our implementation of the layer construction, the Tikhonov matrix is defined according to the 2 nd order view regularization scheme presented in subsection V-D1 (see Eq. ( 29)). It encourages smooth variations between close viewpoints generated from the FDL model, which is intuitively a desirable property for a Light Field.

B. FDL Calibration

We additionally propose a calibration method that jointly estimates the layers' disparity values d k and the angular coordinates u j of the input images B j . For simplicity, the calibration is restricted to the case of sub-aperture images with infinitely small aperture such that ψ j " δ. In this case, Eq. ( 17) has a simpler expression A jk " e `2iπuj d k ωx . Note that in the general configuration where the aperture is unknown, both the aperture functions ψ j and the refocus parameters s j would also need to be estimated. However, that generalization is out of the scope of this paper and is left for future work.

In what follows we express A as a matrix function A : R mˆn Ñ C mˆn such that, for a matrix M P R mˆn , ApMq j,k " e `2iπM j,k . The calibration problem can then be stated in a similar way as the layer construction problem [START_REF] Takahashi | Good group sparsity prior for light field interpolation[END_REF] by treating the calibration parameters as unknowns in the minimization. However, these parameters do not depend on the frequency. Hence, the function to minimize is expressed as a sum over the Q frequency components ω q

x (Q is equal to the number of pixels in each input image):

min x q ,u,d Q ÿ q"1 ´ Apω q x ud J qx q ´bq 2 2 `λ Γx q 2 2 ¯, (20)
where the input view positions u j and the disparity values d k are arranged in the column vectors u and d respectively (ud J P R mˆn ). The vectors x q and b q contain the Fourier coefficients of, respectively, the disparity layers and the input images at the frequency ω q x (i.e. x q k " Lk pω q x q and b q j " Bj pω q x q). Unlike the layer construction problem, the matrix Γ is defined according to the 2 nd order layer regularization approach detailed in subsection V-D2 (see Eq. ( 30)).

In order to solve this problem, we perform a gradient descent along the vectors of parameters u and d. At each iteration, the current estimate of u and d is first used to update the layers values in each vector x q using Eq. [START_REF] Vagharshakyan | Accelerated shearletdomain light field reconstruction[END_REF]. The layers values are then used to compute the gradients ∇u and ∇d of the objective function in Eq. ( 20) along u and d respectively.

By differentiation with respect to each element of u and d, one can show that the gradients are expressed as:

∇u j " ÿ k d k ∇P j,k (21)
∇d k " ÿ j u j ∇P j,k , (22)
with P " ud J , and each column k of the corresponding gradient matrix ∇P is:

∇P k " 4π ÿ q Im ´ωq x x q k ¨Apω q x Pq k ˝pApω q x Pqx q ´bq q ¯, (23)
where ˝is the Hadamard product (i.e. element-wise multiplication), and Im is the imaginary part. Note that the matrix Apω q

x Pq was previously constructed to determine x q , and can thus be re-used in Eq. ( 23) to efficiently compute the gradients.

The computations are further accelerated by performing a stochastic gradient descent, where only a small subset of the Q frequency components ω q

x is selected randomly at each iteration for the gradients computation. In our implementation, subsets of 4096 frequency components were selected. Therefore, the computational cost per iteration does not depend on the image resolution.

Finally, given the gradients, the updated vectors u 1 and d 1 are then computed as

u 1 " u ´α ∇u ` ∇u 2
, and

d 1 " d ´α ∇d ` ∇d 2 , (24)
where is a small value encouraging the stability of the algorithm when the gradients become small (i.e. when u and d are close to an optimum). In our experiments a fixed value α " 0.2 was used.

C. Calibration of a Relaxed FDL Model

Due to the assumptions of non-occluded Lambertian Light Field in the problem definition, the FDL representation may introduce distortions in occluded regions or for strong non-Lambertian effects. In order to better cope with this limitation, we propose a relaxation of the FDL model. In the original FDL calibration (see previous subsection), the parameters in u and d form a parameter matrix P " ud J , where each element P j,k corresponds to the shift applied to the k th layer to approximate the j th input sub-aperture image. Since u and d are column vectors, the rank of the shift parameter matrix P is equal to 1. Hence, the shift applied to each layer is proportional to the angular coordinates of the view to reconstruct. However, this property is not suitable for non-Lambertian surfaces and occlusions. Therefore, in the relaxed FDL model, we remove the rank 1 constraint by searching directly for the matrix P instead of u and d in the calibration algorithm. For that purpose, the gradient descent is applied using the gradients ∇P defined in Eq. [START_REF] Mousnier | Partial light field tomographic reconstruction from a fixed-camera focal stack[END_REF]. In order to obtain satisfying convergence, the matrix P is initialized by the solution of the original rank 1 constrained problem.

The advantage of the relaxed model for handling occlusions or strongly non-Lambertian effects is illustrated in Section VIII-B for the Light Field denoising application. Note that, for this application, we only need to reconstruct the input views for which the shift parameters were obtained in the matrix P. Hence, at each frequency ω q

x , the Fourier coefficients of the reconstructed views are given by Apω q

x Pqx q . However, for more general rendering tasks (e.g. view interpolation, refocusing), the rank 1 constraint is necessary in order to apply the rendering equation [START_REF] Yeung | Fast light field reconstruction with deep coarse-to-fine modeling of spatial-angular clues[END_REF]. Therefore, in the rest of the paper, we consider the original FDL model unless otherwise stated.

D. Regularization Schemes

In order to keep the layer construction and calibration problems in Eqs. [START_REF] Takahashi | Good group sparsity prior for light field interpolation[END_REF] and ( 20) easy to solve, we have used a Tikhonov regularization. The definition of the Tikhonov matrix Γ depends on the intended objective. In the most basic form, classical l 2 regularization is obtained by simply taking Γ equal to the identity matrix I n . This prevents the values in x (i.e. Fourier coefficients of the disparity layers) from taking too high values, which reduces the noise, but also results in a loss of details. Furthermore, for the calibration problem, using the l 2 regularization may not accurately estimate the disparity distribution of the Light Field as illustrated by the calibration results of Fig. 3. The sorted disparity values shown in Fig. 3(b) can be interpreted as an estimation of the cumulative disparity distribution.

In this section, we present two schemes referred to as 2 nd order view regularization and 2 nd order layer regularization which are better suited to the layer construction and the calibration respectively.

1) 2 nd Order View Regularization: For the layer construction, we want to encourage smooth variations between close viewpoints generated from the layers. For that purpose, at each frequency, we penalize the second derivative of sub-aperture images generated from the model with respect to the angular coordinate. From the expression of sub-aperture images in Eq. ( 10), the second derivative at a coordinate u 0 is given by:

B 2 pω x q Bu 2 ˇˇu "u0 " ´4π 2 ÿ k w x 2 d k 2 e `2iπu0d k ωx Lk pω x q.
(25) Let us now consider a set R of angular coordinates to regularize. From Eq. ( 25), the 2 nd order view regularization then consists in constructing a matrix Γ where each row l is associated to an angular coordinate u l P R and is defined by Γ l,k " w x 2 d k 2 e `2iπu l d k ωx . For simplicity, we ignore the constant factor ´4π 2 in the definition of Γ since it can be accounted for in the regularization parameter λ. Intuitively, in order to apply the regularization at every coordinate in the camera plane, one would need to define a matrix Γ with an infinite number of rows (i.e. infinite set R). Although this might seem impractical, we show in what follows that it can be done by observing that the solution to the regularized problem in Eq. ( 19) only requires the knowledge of Γ ˚Γ of finite size n ˆn.

First, let us take a continuous interval R " r´r{2, r{2s of size r, instead of a discrete set. Then, each column Γ k must be interpreted instead as a function Γ k : u Þ Ñ w x 2 d k 2 e `2iπud k ωx . Thus, each element of Γ ˚Γ is given by the following inner product:

rΓ ˚Γs k1,k2 " xΓ k1 , Γ k2 y " ż r{2 ´r{2 Γ k1 puq ¨Ě Γ k2 puqdu, (26)
" ω x 4 d k1 2 d k2 2 ż r{2 ´r{2 e 2iπupd k 1 ´dk 2 qωx du, (27)
" For convenience, we replace Γ ˚Γ in Eq. ( 19) by the scaled version Ć Γ ˚Γ " Γ ˚Γ r , so that the amount of regularization remains independent of the size of the integration domain in Eq. ( 26) (i.e. area of the portion of the camera plane covered by the regularization). In our implementation of the 2 nd order view regularization, we consider the full camera plane (i.e. R " p´8, 8q) by making r tend towards the infinity. In this case, Ć Γ ˚Γ simply becomes a diagonal matrix such that:

ω x 4 d k1 2 d k2 2 r sincprpd k1 ´dk2 qω x q. ( 28
)
r Ć Γ ˚Γs k,k " ω x 4 d k 4 . (29)
Note that for ω x " 0, Eq. ( 29) gives Ć Γ ˚Γ " 0, which does not produce any regularization. Therefore, in practice, a small value is added to the diagonal elements of Ć Γ ˚Γ. This is equivalent to adding a l 2 regularization term ¨ x 2 2 to the problem [START_REF] Takahashi | Good group sparsity prior for light field interpolation[END_REF].

The advantage of the proposed 2 nd order view regularization over the simple l 2 regularization is demonstrated in Fig. 4 for the view interpolation problem.

2) 2 nd Order Layer Regularization: For the calibration step, however, it is preferable to use a regularization that does not depend on the parameters d k . Otherwise, the expression of the gradients ∇d would need to take the regularization term into account, which can result in more complex gradient computations and a slower convergence of the algorithm. Instead, we encourage smooth variations between successive layers by defining Γ as a discrete approximation of the secondorder differential operator as follows:

Γ " » - - - - - - ´2 1 1 ´2 1 . . . . . . . . . 1 ´2 1 1 ´2fi ffi ffi ffi ffi ffi fl . ( 30
)
The 2 nd order layer regularization thus penalizes large differences between neighboring layers, which results in a more uniform distribution of the disparity values. As shown in Fig. 3(b), the calibration using the l 2 regularization (i.e. Γ " I n ) tends to find too many disparity values close to the dominant disparity (between 0 and 0.5 in the figure) and may underestimate other parts of the light field (e.g. with disparities close to -1). The 2 nd order layer regularization attenuates that effect by encouraging a more uniform distribution. Note, however, that for the layer construction, smoothness along the layers is not desirable and produces artifacts along the edges in the images rendered from the FDL model as shown in Fig. 4(c).

VI. LIGHT FIELD RENDERING A. Implementation

Knowing the layers and their disparities, any view with arbitrary aperture and focus can be rendered in the Fourier domain by applying the FDL decomposition equation ( 16) and by computing the inverse Fourier transform. The interpretation in the spatial domain is that each layer is shifted (i.e. multiplied by e `2iπu0d k ωx in the Fourier domain) and filtered (i.e. multiplied by ψpω x ps ´dk qq in the Fourier domain).

However, except for specific aperture shapes (e.g. square or dirac), the Fourier transform ψ of the aperture function has no analytical expression. Hence, in our implementation, an approximation is obtained by drawing the aperture shape as an image, and by taking its discrete Fourier transform. A linear interpolation is used to determine ψpω x ps ´dk qq from the discrete frequency samples. In order to increase the accuracy of this approximation, the aperture image is zero-padded before computing its Fourier transform, which increases the resolution in the spectral domain. The effect of the padding is illustrated in Fig. 5 showing that the resulting increased spectral domain resolution allows a better preservation of the contrasts in the rendered images. Note also that the spatial resolution of the aperture image can be taken arbitrarily large without affecting the aperture size in the rendered image. In practice, we control the aperture size by replacing ψpω x ps ´dk qq in Eq. ( 16) by ψpf ω x ps ´dk qq with a scaling factor f . For example, taking f " 0 simulates a sub-aperture image. On the other hand, taking a large factor f produces a shallow depth of field effect without affecting the complexity of the method.

B. Visual Results and Discussion

In Fig. 6, we show a comparison of our approach with the conventional Shift and Sum refocusing method [START_REF] Ng | Light Field Photography with a Hand-Held Plenoptic Camera[END_REF]. Thanks to the regularization used in the layer construction (see Section V-D), the light field is extended angularly, which allows rendering images with reduced angular aliasing as in Fig. 6(b). As shown in Fig. 6(c), the aperture can also be taken larger than the baseline of the input Light Field from which the FDL model is constructed. The possibility to change the shape of the aperture also gives better control over the bokeh (i.e. the quality of the out-of-focus regions).

Another important and often overlooked factor influencing the bokeh is the color space of the input image data. The refocusing Eqs. ( 1) and ( 3) assume that the Light Field rays Lpx, y, u, vq are proportional to luminance data. However, images are traditionally represented in a non-linear color space (e.g. gamma corrected) in order to account for the non-linearity of the human visual system as well as typical displays. This non-linearity affects the appearance of the refocused images. More realistic bokeh can be achieved by applying inverse gamma correction to the input images before the layer construction, and by applying gamma correction back after the rendering step, as shown in Fig. 7. One limitation of this approach, however, is that computing the layers in a linear color space also produces more visible artifacts. This is due to the fact that solving the layer construction problem [START_REF] Takahashi | Good group sparsity prior for light field interpolation[END_REF] in a linear color space does not take the non-linearity of human perception of luminance into account.

A demonstration of our rendering application is available in the supplementary materials. It shows that our approach can be used for controlling simultaneously the viewpoint, the aperture shape and size, and the focus in real-time.

VII. COMPLEXITY

In this section, we analyse the complexity of the different processing steps proposed in our approach.

First, regarding the spatial resolution, our implementation takes advantage of the symmetry property of the Fourier transform for real signals. Given a real valued function g, its Fourier transform ĝ is such that ĝp´ωq " ĝpωq. Hence, for the layer construction and the rendering algorithms, only half of the discrete spectrum must be computed. The remaining frequencies are directly obtained by copying and taking the complex conjugate of the first half. This symmetry property is also used in the calibration step by selecting random frequencies only in one half of the spectrum for the stochastic gradient descent. Regarding memory requirements, using only half of the frequencies in the Fourier representation compensates the fact that complex numbers require twice as much memory as real numbers. In addition to the spatial resolution, another important factor to consider for the complexity is the number of layers used. In order to analyse the number of layers needed for our model, we have performed an experiment where Light Field sub-aperture images are used for FDL calibration and construction with varying numbers of layers. The input images are then reconstructed via FDL rendering. Several Light Fields from various datasets were tested including natural Light Fields captured with a first generation Lytro camera ('Totoro Waterfall' [START_REF] Jiang | Light field compression with homography-based low-rank approximation[END_REF]), a Lytro Illum camera ('Fruits' and 'Toys' [START_REF] Pendu | Light field inpainting propagation via low rank matrix completion[END_REF]), a Raytrix R8 camera (first frame of the Light Field 'ChessPieces' [START_REF]Light field video dataset captured by a r8 raytrix camera[END_REF]), a traditional camera moving on a gantry ('Lego Knights' [START_REF]The stanford light field archive[END_REF]), and synthetic Light Fields ('papillon' and 'buddha' [START_REF] Wanner | Datasets and benchmarks for densely sampled 4d light fields[END_REF]). Fig. 8 presents the average Peak Signal to Noise Ratio (PSNR) of all the reconstructed sub-aperture images with respect to the number of layers. The saturation of the PSNR curves show that 15 to 30 layers are generally sufficient, and using more layers does not significantly change the results. Light Fields with a large baseline, and thus large differences between extreme viewpoints (e.g. Lego 5x5, Lego 17x17), typically require more layers than smaller baseline Light Fields. However a very large range of viewpoints also implies large occlusions, which results in a significant loss when reconstructing the sub-aperture images, as observed with Lego 5x5 and Lego 17x17. This makes the analysis difficult for very large baseline Light Fields which may still benefit from more than 30 layers in the ideal non-occluded case. Note that for the small baseline Light Fields captured with Lytro cameras, the relatively low PSNR is due to inaccuracies in the input data (e.g. noise, color differences between views) which are reduced in our reconstruction. More details and results on the reduction of noise and color inconsistencies are presented in section VIII-B.

The computing times for the different processing steps are presented in Table II for all the Light Fields in Fig. 8. The table also provides the number of iterations needed for the calibration since it has a large impact on the computing time. For this experiment, 30 layers were used although less layers would be sufficient for some of the Light Fields (e.g. see 'ChessPieces', 'Toys', 'Totoro Waterfall' in Fig. 8). The processing times were measured using an Intel Core i7- 7700 CPU and a Nvidia GeForce GTX 1080 GPU. All the proposed steps were implemented in Matlab using the parallel computing toolbox in order to process the different frequencies in parallel on the GPU. Although the processing time for the calibration and layer construction steps is affected by the number of input views (i.e. angular resolution), the rendering time only depends on the spatial resolution and the number of layers. For example, the same rendering time is obtained for 'Lego Knights' using either 5x5 or 17x17 input views. This is particularly advantageous when rendering images with a large aperture size that are traditionally obtained by averaging a large number of sub-aperture images with the Shift and Sum algorithm [START_REF] Ng | Light Field Photography with a Hand-Held Plenoptic Camera[END_REF]. Note that similarly to our method, the Fourier Slice refocusing algorithm [START_REF] Ng | Fourier slice photography[END_REF] first transforms the Light Field in an intermediate representation in order to perform fast refocusing. In that method the refocusing complexity does not depend on the number of input views. However, since the intermediate representation is the 4D Fourier transform of the Light Field, the memory requirement remains proportional to the number of views. Furthermore, it does not allow changing the aperture or the viewpoint.

VIII. DIRECT APPLICATIONS

A. View Interpolation and Extrapolation

As shown in Fig. 6, our approach reduces the angular aliasing and allows extending the aperture size for an input Light Field consisting of a sparse set of sub-aperture images. Equivalently, the layer construction can be seen as a view interpolation and extrapolation method since sub-aperture images with arbitrary angular coordinates can be rendered from the FDL model.

1) Comparison and Combination with Sparsity Based Approach: For evaluating the view interpolation capability of the method, we have compared our results with the recent method of Vagharshakyan et al. [START_REF]Light field reconstruction using shearlet transform[END_REF] that enforces the sparsity of the epipolar images of the Light Field in the shearlet domain.

For the comparison, the synthetic Light Field 'papillon' was used with several sampling configurations of the input views as illustrated in Fig. 9. The average PSNR of the interpolated views and the computation times are shown in Table III. The table also includes the results for two schemes combining the FDL and the Shearlet approaches. In the 'Shearlet(full)+FDL', the full Light Field is reconstructed with [START_REF]Light field reconstruction using shearlet transform[END_REF]. Then, a FDL model is computed from the input views and the previously reconstructed intermediate views. The final reconstruction is obtained with FDL rendering. The 'Shearlet(border)+FDL' scheme is similar, but only the views at the periphery of the Light Field are interpolated using [START_REF]Light field reconstruction using shearlet transform[END_REF]. Our experiments with the Shearlet method were performed using the author's implementation of the epipolar image reconstruction with 100 iterations. The reconstruction of the full Light Field was performed by scanning the epipolar images in the 'direct' order as detailed in [START_REF]Light field reconstruction using shearlet transform[END_REF].

The results in Table III show that our FDL view interpolation performs particularly well in the 'border' configuration. The Shearlet method, on the other hand, does not fully take advantage of this configuration as it processes each epipolar image independently using only the input views located on the corresponding row or column in the Light Field. In the other configurations, a higher PSNR is also observed with the FDL method. However, when the input Light Field is too sparse (e.g. 2x2), it produces ringing artifacts as shown in Fig. 10. In this situation, the Shearlet method better preserves the edges than our FDL method, but it tends to blur the fine details and introduce color distortions, which explains the lower PSNR. The best strategy for such sparse inputs is then to combine the two methods with the 'Shearlet(border)+FDL' scheme which does not produce ringing artifacts. In comparison with the Shearlet method alone, this scheme better keeps the details and color consistency as shown in Fig. 10. It is also significantly faster since the reconstruction of the interior views is performed using our FDL approach which requires a negligible computing time compared to the Shearlet method.

2) Comparison with Deep Learning Approaches: We have additionally performed comparisons with three state-of-the art deep learning approaches [START_REF] Kalantari | Learning-based view synthesis for light field cameras[END_REF], [START_REF] Wu | Light field reconstruction using deep convolutional network on epi[END_REF], [START_REF] Yeung | Fast light field reconstruction with deep coarse-to-fine modeling of spatial-angular clues[END_REF]. The method of Kalantari et al. [START_REF] Kalantari | Learning-based view synthesis for light field cameras[END_REF] uses two Convolutional Neural Networks (CNN) where the first one determines a disparity map, while the second one refines the disparity compensated views. Yeung et al. [START_REF] Yeung | Fast light field reconstruction with deep coarse-to-fine modeling of spatial-angular clues[END_REF] have proposed a CNN defining an end-to-end mapping between sparse input views and the dense Light Field. Finally, the method of Wu et al. performs superresolution of the epipolar images. For the comparison, the natural Light Fields 'Fruits' and 'Toys' [START_REF] Pendu | Light field inpainting propagation via low rank matrix completion[END_REF] as well as 'Friends 1' and 'Vespa' [START_REF] Rerabek | New light field image dataset[END_REF] were used with gamma correction. They were captured with a Lytro Illum camera, similarly to the Light Fields used for the training in the tested deep learning approaches. The PSNR results and average runtimes are presented in Table IV for the configurations with 2x2 input and 8x8 output views (2x2Ñ8x8) and with 3x3 input and 7x7 output views (3x3Ñ7x7). Note that for the method of Yeung et al. [START_REF] Yeung | Fast light field reconstruction with deep coarse-to-fine modeling of spatial-angular clues[END_REF], the available CNN model only applies to the 2x2Ñ8x8 configuration. Although the method of Kalantari et al. [START_REF] Kalantari | Learning-based view synthesis for light field cameras[END_REF] was designed for 2x2 input views, we have used it in the 3x3Ñ7x7 configuration by treating the quadrants of the Light Field separately. Furthermore, in [START_REF] Kalantari | Learning-based view synthesis for light field cameras[END_REF] and [START_REF] Yeung | Fast light field reconstruction with deep coarse-to-fine modeling of spatial-angular clues[END_REF], the output Light Field is cropped in the spatial dimensions. For the comparison, the same crop (12 pixels on each side) was applied to compute all the PSNR values in Table IV.

The CNNs used in [START_REF] Kalantari | Learning-based view synthesis for light field cameras[END_REF] and [START_REF] Yeung | Fast light field reconstruction with deep coarse-to-fine modeling of spatial-angular clues[END_REF] are particularly efficient for the 2x2Ñ8x8 configuration for which they were designed and trained. In this case, they obtain better results than our FDL approach which is penalized by ringing artifacts, similarly to the results shown in Fig. 10 (a). However, the performances of Kalantari et al's method drop in the 3x3Ñ7x7 configuration, indicating that the model learned for 2x2 input views does not generalize well to more dense inputs. Our FDL model, on the other hand, provides significantly better interpolation in this case thanks to the higher density of input views, leading to a better conditioning of the FDL construction problem. While Yeung et al. have also reported in [START_REF] Yeung | Fast light field reconstruction with deep coarse-to-fine modeling of spatial-angular clues[END_REF] high quality results for the 3x3Ñ7x7 scenario, an adapted version of their model must be designed and trained for each configuration. Our PSNR results are more comparable to those of Wu et al. [START_REF] Wu | Light field reconstruction using deep convolutional network on epi[END_REF]. Visual comparisons in Fig. 11 additionally reveals that our method better preserves the details. However, highly textured occluded areas become apparent on the fore- ground, making the borders of occluding objects appear semitransparent. Note also that for the tested Light Fields and the two configurations 2x2Ñ8x8 and 3x3Ñ7x7, the average computing time of [START_REF] Wu | Light field reconstruction using deep convolutional network on epi[END_REF] was about 16 minutes against only 6 seconds for our method (including the calibration step).

B. Denoising and Other Filtering Effects

Another direct application of our FDL construction algorithm is the Light Field denoising problem. Denoising is naturally obtained by constructing a FDL model from a Light Field, and by rendering images using the same angular coordinates as the input sub-aperture images. The noise is filtered in this approach thanks to the model definition that enforces the Light Field dimensionality gap prior [START_REF] Levin | Linear view synthesis using a dimensionality gap light field prior[END_REF] (see Section IV-A). Further noise reduction is also allowed by the proposed regularization that better ensures consistency between views. Note that for the same reasons, the method also ensures color and illumination consistency along the angular dimensions. Hence it can serve simultaneously as a denoising and color correction method for Light Fields with variations of color and illumination between the sub-aperture images. In practice, this is particularly useful for Light Fields captured with plenoptic cameras (e.g. Lytro), as shown in figure 12. Similarly, by ensuring consistency between multiple views, the FDL model and its regularization act as a spatial antialiasing filter in the case of spatially aliased Light Fields (see Fig. 13). This situation typically arises with Light Fields captured with lenslet-based plenoptic cameras because of the low lenslet resolution. It has been shown in [START_REF] Liang | A light transport framework for lenslet light field cameras[END_REF] that such Light Fields are spatially non-bandlimited and that despite the low spatial sampling, higher frequencies than the Nyquist frequency can be recovered in theory. However, we construct the layers at the same spatial resolution as the input views. Therefore, while an anti-aliasing effect is produced, these higher frequencies are not recovered in our method. Constructing higher resolution layers to fully exploit the 4D Light Field sampling would allow Light Field super-resolution and is left for future work.

For Light Fields with a wide baseline, our processing may also introduce distortions due to the assumptions of nonoccluded Lambertian Light Field. In this case, better results can be obtained using the relaxed version of the FDL model presented in Section V-C. Denoising results are presented in Fig. 14 for the Light Field 'Tarot' that exhibits strong non-Lambertian effects. The reflections in the crystal ball are better preserved by the relaxed FDL model than the original one. For the comparison, the results obtained with the LFBM5D [START_REF] Alain | Light field denoising by sparse 5d transform domain collaborative filtering[END_REF] and the Hyperfan 4D filter (HF4D) [START_REF] Dansereau | Light field image denoising using a linear 4d frequency-hyperfan allin-focus filter[END_REF] are also shown. The corresponding PSNR and SSIM results are given in Table V.

The best results were obtained with the LFBM5D method that completely removes the noise while preserving most of the details in the image. However, this method performs heavy processing and typically requires several hours for a Light Field. The HF4D method is faster since it simply consists in multiplying the 4D Fourier Transform of the Light Field by a 4D filter. This hyperfan filter is designed to attenuate the frequencies outside of the theoretical region of support of the Light Field in the 4D Fourier domain, under the assumption of a non-occluded Lambertian scene. With comparable processing times, our FDL approach removes more noise thanks to the linear optimization used instead of the direct filtering in the HF4D method. Furthermore, our model can be more easily generalized for the case of scenes with occlusions or non-Lambertian surfaces, as shown with the proposed relaxation.

Note however that for a high level of noise (i.e. σ " 50), the relaxation of the FDL model results in lower objective scores despite the better preservation of the non-Lambertian surface observed in Fig. 14. The reason is that the relaxed model also preserves more noise from the input data. Hence the choice of the original or the relaxed version of the model is a tradeoff between the noise level and the amount of occlusions or non-Lambertian effects.

IX. CONCLUSION AND PERPECTIVES

We have presented a new representation for Light Fields called Fourier Disparity Layers. The light information in the scene is decomposed into several layers, each corresponding to a depth plane parameterized by a disparity value. The decomposition step is formulated in the Fourier domain using a simple linear least square problem per-frequency, hence allowing fast processing with GPU parallelization. We have demonstrated the advantages of the FDL representation for several Light Field processing tasks. Those include real time rendering, view interpolation, denoising, as well as a calibration step that determines the angular coordinates of the input sub-aperture images along with an optimal set of disparity values for each layer.

The computational efficiency of our layer construction method coupled with its flexibility regarding the type of input images opens perspectives for an even larger range of applications. For example, our supplementary video presents results of FDL model construction from a focal stack, where the camera aperture and the focus parameter of each image are known. More generally, the method could even take advantage of a hybrid capture system taking images with different apertures, focusing depth, and varying positions on the camera plane. Furthermore, since our approach does not require a specific pattern of the input view positions, it can greatly simplify the capture of Light Fields by removing the need for the different viewpoints to be regularly sampled on a grid.

However, despite this flexibility, we have shown that the quality of the FDL model produced with our method depends on the input configuration. For example, while very accurate view interpolation is obtained from the sub-aperture images located on the borders of the Light Field, a too sparse viewpoint sampling may lead to ringing artifacts. Furthermore, for Light Fields with a large baseline, the FDL model is more likely to produce visible errors due to large occlusion areas or non-Lambertian effects. In the paper, we have partially addressed these limitations. In the case of a too sparse viewpoint sampling, a combination of our view interpolation with a state of the art approach was shown to outperform either of the two methods taken individually. A better handling of occlusions and non-Lambertian effects was also obtained, for the denoising application, using a relaxed version of the model.

In the aim of further extending the applicability of the FDL approach, future work may focus on generalizing the relaxed model to view interpolation, or including an additional prior directly in the layer construction (e.g. sparsity in the Shearlet domain similarly to [START_REF]Light field reconstruction using shearlet transform[END_REF]) to better cope with very sparse or very noisy Light Fields. A generalization of the calibration to wide aperture images would also be a valuable tool to facilitate the creation of Light Fields from less conventional input data such as focal stacks.

APPENDIX A PROOF OF SPARSITY PRIOR

The spatial regions Ω k are defined for the central view at u " 0. The corresponding regions Ω u k can also be defined for any view u by Ω u k " tx P R | x `ud k P Ω k u. The expression in Eq. ( 4) is then formulated as: @k P v1, nw, @px, uq P Ω u k ˆR, Lpx, uq " Lpx`ud k , 0q. [START_REF] Wanner | Datasets and benchmarks for densely sampled 4d light fields[END_REF] In the assumption of a non-occluded Light Field, the regions Ω u k are such that @u P R, x P Rz Ť k Ω u k , Lpx, uq " 0, and for any fixed u P R, the sets Ω u k are pairwise disjoint. Hence, the Fourier Transform of the Light Field is given by: Lpω x , ω u q " ij R 2 e ´2iπpxωx`uωuq Lpx, uqdxdu (

) 34
where we define Lk pω x q " ş Ω k e ´2iπxωx Lpx, 0qdx.

Fig. 2 .

 2 Fig. 2. Two-plane parameterization. The focal plane px, yq is parallel to the camera plane pu, vq and located at a distance z 0 .

Fig. 3 .Fig. 4 .

 34 Fig. 3. Example of calibration results for the synthetic Light Field 'papillon' [31] with 30 layers: (a) Estimated angular coordinates of the views (the ground truth is the 9x9 regular grid). (b) Sorted disparity values. The vertical axis represents the normalized index of each of the sorted values. The blue curve was obtained by sorting all the pixels of the ground truth disparity map in order to visualize the cumulative disparity distribution.

Fig. 5 .

 5 Fig. 5. Top row: no padding. Bottom row: with zero-padding. From left to right: Aperture image; Real part of the Fourier transform of the aperture image (up-sampled with linear interpolation for the top row for the comparison); Rendered result with magnified details. Zero-padding removes the aliasing in the Fourier transform of the aperture, which better preserves contrasts in the final rendered image (best viewed zoomed in).

Fig. 6 .

 6 Fig.6. Refocused image using 5x5 views of the Light Field 'Lego Knights' from the Stanford dataset[START_REF]The stanford light field archive[END_REF]: (a) Shift and Sum method[START_REF] Ng | Light Field Photography with a Hand-Held Plenoptic Camera[END_REF], (b) Our method using a square aperture covering the angular coordinates of the input 5x5 views. (c) Our method using a larger circular aperture. Angular aliasing is reduced in our method by considering all the angular coordinates within the aperture.

Fig. 7 .

 7 Fig.7. Rendered result when the processing is performed in (a) gamma corrected space, (b) linear space. Linear space processing results in brighter out-of-focus regions where the bokeh shape (i.e. aperture shape shown in the top left corner) appears more clearly. However, it also increases the artifacts (see details in the green rectangles).

Fig. 8 .

 8 Fig. 8. Reconstruction PSNR (average over views) for different Light Fields with respect to the number of layers used in the FDL model. For the Light Field 'Lego Knights', two different results are shown using either the full 17x17 views or only the 5x5 central views.

Fig. 9 .

 9 Fig. 9. Tested sampling configurations of the input views (shown in red) for the view interpolation: (a) '2x2', (b) '3x3', (c) '5x5', (d) 'border'.

Fig. 10 .

 10 Fig. 10. Interpolation of the central view from the 4 corner views of the Light Field 'papillon' (i.e. 2x2 configuration) using: (a) only FDL, (b) only Shearlet [20], (c) 'Shearlet(full)+FDL', (d) 'Shearlet(border)+FDL'. The bottom images show the residual computed with the ground truth central view (magnified 5x). A video of the results is available in the supplementary materials.

Fig. 11 .

 11 Fig. 11. Interpolated view (5,5) of the Light Field 'Vespa' from the configuration 3x3Ñ7x7: (a) Wu et al. (b) FDL. Our FDL method produces more detailed results (see red boxes). However, textures of occluded regions become apparent on the occluding objects (see green boxes).

Fig. 12 .

 12 Fig. 12. Simultaneous denoising and color correction of the Light Field 'TotoroWaterfall' (captured with a Lytro camera): (a) Original (b) Processed with the FDL method. The top image corresponds to a view on the side of the Light Field. Epipolar images corresponding to the horizontal red line are shown on the bottom to visualize the color consistency between views.

Fig. 13 .

 13 Fig. 13. Spatial anti-aliasing effect of the FDL model: (a) Original central view (detail of the Light Field 'papillon'), (b) Reconstructed view from the FDL model. The effect is particularly visible on the butterfly's antennas.

Fig. 14 .

 14 Fig. 14. Denoising results for the Light Field 'Tarot'. A view on the corner of the Light Field is shown. For the test, the Light Field was corrupted by gaussian noise with standard deviation σ " 50.

"ee( 33 )"

 33 [START_REF] Wanner | Datasets and benchmarks for densely sampled 4d light fields[END_REF], and by change of variable we obtain:Lpω x , ω u q " ´2iπxωx Lpx `ud k , 0qdx ff du " ´2iπpx´ud k qωx Lpx, 0qdx  du.By re-arranging the terms the result is:Lpω x , ω u q " δpω u ´dk ω x q Lk pω x q ı .

TABLE I TABLE OF NOTATION

 IOF Sub-aperture image at position u 0 (Lu 0 pxq " Lpx, u 0 q).

	Symbols	Description
	Lpx, uq	Light Field with spatial coordinate x and angular coordinate u.
	ωx, ωu	Respectively spatial and angular frequencies.
	Lu 0
	B s u 0	Image with refocus parameter s and position u 0 .
	ψ	Aperture function.
	u 0	Angular coordinate of the view to reconstruct.
	u l	Angular coordinate of a known input view.
	d k	Disparity value in the Light Field.
	Ω u k	Region of disparity d k in an image of view position u.
	δ	Dirac delta function
	f	Fourier transform of a function f .
	A J

TABLE II COMPUTING

 II TIMES FOR OUR MATLAB IMPLEMENTATION USING AN INTEL CORE I7-7700 CPU AND AN NVIDIA GEFORCE GTX 1080 GPU. FOR EVERY LIGHT FIELD 30 LAYERS WERE USED.

	Input Light Field (resolution)	Calibration time (#iter)	Layer construction	Render
	Lego Knights (1024x1024x5x5)	9.2 s (184)	2.2 s	35 ms
	Lego Knights (1024x1024x17x17) 56 s (403)	7.2 s	35 ms
	Chess Pieces (1920x1080x5x5)	5.8 s (90)	4.3 s	53 ms
	Fruits (625x434x9x9)	5.1 s (70)	0.9 s	12 ms
	Toys (625x434x9x9)	6.2 s (91)	0.9 s	11 ms
	Totoro WaterFall (379x379x7x7)	3.8 s (72)	0.4 s	7 ms
	papillon (768x768x9x9)	7.8 s (109)	1.9 s	20 ms
	buddha (768x768x9x9)	8.1 s (113)	1.9 s	19 ms

TABLE III VIEW

 III INTERPOLATION RESULTS FOR THE LIGHT FIELD 'PAPILLON'. THE RESULTS GIVEN ARE THE AVERAGE PSNR OF THE INTERPOLATED VIEWS AND THE TOTAL RECONSTRUCTION TIME (INCLUDING CALIBRATION FORTHE FDL METHOD).

	Method \ Input views	2x2	3x3	5x5	border
	Shearlet [20]	34.8 dB 143 mn	40.8 dB 157 mn	41.8 dB 185 mn	37.2 dB 105 mn
	FDL	36.8 dB 5.4 s	40.9 dB 7.8 s	42.9 dB 8.9 s	43.1 dB 8.8 s
	Shearlet(full)+FDL	35.5 dB 143 mn	41.4 dB 157 mn	42.6 dB 185 mn	40.1 dB 105 mn
	Shearlet(border)+FDL	38.4 dB 53 mn	42.0 dB 53 mn	42.7 dB 53 mn	x

TABLE IV VIEW

 IV INTERPOLATION RESULTS FOR NATURAL LIGHT FIELDS. THE RESULTS GIVEN ARE THE AVERAGE PSNR OF THE INTERPOLATED VIEWS.

			Fruits	Toys Friends 1 Vespa	average runtime
	Yeung et al. [16]	2x2Ñ8x8	33.1	32.9	37.5	34.4	40s
	Kalantari	2x2Ñ8x8	34.2	33.2	37.1	34.4	10mn30s
	et al. [12]	3x3Ñ7x7	31.7	32.1	36.4	34.0	7mn
	Wu	2x2Ñ8x8	28.3	32.0	36.6	32.4	19mn
	et al. [14]	3x3Ñ7x7	35.3	35.7	39.4	37.7	13mn30s
	FDL	2x2Ñ8x8 3x3Ñ7x7	32.4 35.9	32.5 35.3	36.8 40.0	32.1 36.9	6s 6s

TABLE V DENOISING

 V RESULTS FOR THE LIGHT FIELD 'TAROT' USING 11X11 INPUT VIEWS AND NOISE WITH STANDARD DEVIATION σ " 10 AND σ " 50.

		σ " 10	σ " 50	average
		PSNR	SSIM	PSNR	SSIM	runtime
	noisy input	28.1	0.787	14.2	0.283
	LFBM5D [38]	39.3	0.983	32.7	0.937	"5h30mn
	HF4D [29]	31.2	0.888	22.3	0.603	"40s
	FDL (original)	32.0	0.945	26.6	0.785	"1mn30s
	FDL (relaxed)	33.4	0.946	26.1	0.760	"2mn15s

This project has been supported in part by the Science Foundation Ireland (SFI) under the Grant Number 15/RP/2776 and in part by the EU H2020 Research and Innovation Programme under grant agreement No 694122 (ERC advanced grant CLIM).