Navaneeth Kamballur Kottayil

Senior Member, IEEE Giuseppe Valenzise

Fellow, IEEE Frederic Dufaux

Senior Member, IEEE Irene Cheng

Blind Quality Estimation by Disentangling Perceptual and Noisy Features in High Dynamic Range Images

High Dynamic Range (HDR) image visual quality assessment in the absence of a reference image is challenging. This research topic has not been adequately studied largely due to the high cost of HDR display devices. Nevertheless, HDR imaging technology has attracted increasing attention because it provides more realistic content, consistent to what the Human Visual System perceives. We propose a new No-Reference Image Quality Assessment (NR-IQA) model for HDR data based on convolutional neural networks. The proposed model is able to detect visual artifacts, taking into consideration perceptual masking effects, in a distorted HDR image without any reference. The error and perceptual masking values are measured separately, yet sequentially, and then processed by a Mixing function to predict the perceived quality of the distorted image. Instead of using simple stimuli and psychovisual experiments, perceptual masking effects are computed from a set of annotated HDR images during our training process. Experimental results demonstrate that our proposed NR-IQA model can predict HDR image quality as accurately as state-of-the-art full-reference IQA methods.

I. INTRODUCTION

High dynamic range (HDR) images can present a much larger range of luminance compared to conventional images. This larger range of luminance is achieved by using 16-32 bit floating point values, instead of a conventional 8 bits per pixel integer representation. Viewers are able to perceive more vivid colors and scene content compared to viewing on a conventional Standard Dynamic Range (SDR) or Low Dynamic Range (LDR) display. This generates a better quality of viewing experience [START_REF] Dufaux | High dynamic range video : from acquisition to display and applications[END_REF]. The advancement of HDR imaging technology has changed the landscape of the entire multimedia communication pipeline from capturing, processing and transmission, to the visualization of HDR content [START_REF] Chalmers | {HDR} video past, present and future: A perspective[END_REF]. This technology has become an important development in the consumer market, e.g., TV and photography, with the support of industrial investments.

In this work, we consider evaluating HDR image quality on a HDR compatible display. Image quality assessment (IQA) can be broadly categorized into Full-Reference (FR) and No-Reference (NR). In FR-IQA, the quality of a given image is evaluated by comparing a distorted image with an undistorted version of the same image. In NR-IQA, the quality is evaluated by judging the distorted image only. Since the target audience of the HDR content is the Human Visual System (HVS), one IQA solution is to conduct user subjective tests. However, subjective tests are tedious and time-consuming. Even with the help of massive crowdsourcing projects, e.g., mturk, HDR IQA is difficult due to the high cost of acquiring systems capable of displaying HDR content. A better solution is to develop an efficient NR-IQA model, which can automatically assess HDR content, matching the human perception.

Our proposed model is capable of predicting the perceived HDR image quality and localizing the distortions. We use a convolutional neural network (CNN) based architecture to achieve this goal. Our model addresses low-level image distortions such as artifacts caused by image compression. We do not consider changes in image quality due to highlevel effects, e.g., artistic intent, where complex aesthetical considerations need to be taken into account.

Our contributions in HDR image quality assessment lie in: 1) Proposing an NR-IQA model based on a convolutional neural network architecture, which can separate pixelwise errors from their impact on perception in a distorted image. Our model outperforms other NR-IQA models and is competitive with state-of-the-art HDR full-reference IQA algorithms. 2) Providing an accurate error prediction in a distorted image without a reference image. 3) Predicting the visual masking effects without the need of explicit psychovisual subjective tests. The rest of this paper is organized as follows: Section II discusses the previous work in related fields. Section III explains our motivation for the new approach. In Section IV we discuss the conceptual idea and implementation of the proposed method. Section V deals with performance comparisons and test of our method against other algorithms. Followed by this, in Section VI, we discuss the functionality of each of the subcomponents of our system and show how they work. Section VII concludes the paper.

II. PREVIOUS WORK

In pace with the rapid development of HDR imaging technology, Full-Reference Image Quality Assessment (FR-IQA) of HDR images has been gaining attention in recent years, and a number of studies have addressed the evaluation of quality for HDR images and video [START_REF] Narwaria | HDR-VDP-2.2: a calibrated method for objective quality prediction of highdynamic range and standard images[END_REF], [START_REF] Narwaria | HDR-VQM: An objective quality measure for high dynamic range video[END_REF], [START_REF] Zerman | An extensive performance evaluation of full-reference HDR image quality metrics[END_REF]. However, much less has been done on No-Reference Image Quality Assessment (NR-IQA) of HDR images. Relevant work in the NR-IQA literature can be classified into two general categories:

• LDR NR-IQA, where the quality of images is estimated when it is visualized on LDR displays. • Tone-Mapped NR-IQA, in which color and contrast values of a HDR image are mapped onto a smaller range of color and contrast values, using Tone Mapping Operators (TMO), and the output is evaluated on LDR screens. TMOs are often used to compress the dynamic range of HDR images. Note that in our work, we assume instead that HDR images are directly evaluated on HDR displays.

A. LDR NR-IQA

Machine learning approaches are often used in LDR NR-IQA. These approaches start by creating a feature image and fit a relevant distribution. The parameters of this distribution are used as the feature vector of the distorted image. An early method based on machine learning in LDR NR-IQA is the Blind Image Quality Index (BIQI) [START_REF] Moorthy | A two-step framework for constructing blind image quality indices[END_REF]. It is a two-step process: (1) from a set of features, a Support Vector Machine (SVM) predicts the type of distortion, and (2) a set of Support Vector Regressors (SVRs) predict the score for each distortion type. The final quality score is computed as

score = m i=1 p i • q i , (1)
where m is the distortion type, p i represents the probability of each distortion obtained from the SVM and q i represents the quality score given by each of the SVR's. BIQI [START_REF] Moorthy | A two-step framework for constructing blind image quality indices[END_REF] used Daubechies 9/7 wavelet as the feature image. Many methods follow a similar approach and show good performances on assessing LDR content on LDR screens without a reference image. Examples are Blind/Referenceless Image Spatial QUality Evaluator (BRISQUE) [START_REF] Mittal | No-reference image quality assessment in the spatial domain[END_REF], Distortion Identification-based Image Verity and INtegrity Evalutation (DIIVINE) [START_REF] Moorthy | Blind image quality assessment: From natural scene statistics to perceptual quality[END_REF] and Spatial-Spectral Entropy based Quality metric (SSEQ) [START_REF] Liu | No-reference image quality assessment based on spatial and spectral entropies[END_REF]. These algorithms can have complex features, but just need a single SVR to predict the final quality of an image. BRISQUE [START_REF] Mittal | No-reference image quality assessment in the spatial domain[END_REF] computes the Mean Subtracted Contrast Normalized (MSCN) image as feature using M SCN (i, j) = I(i,j)-µ I,N,i,j σ I,N,i,j +1 , where µ I,N,i,j and σ I,N,i,j represent the mean and variance computed over a local Gaussian window of size N around the point (i, j). I(i, j) is the image intensity at (i, j). DIIVINE [START_REF] Moorthy | Blind image quality assessment: From natural scene statistics to perceptual quality[END_REF] uses divisive normalized steerable pyramid decomposition coefficients to create the feature image.

SSEQ [START_REF] Liu | No-reference image quality assessment based on spatial and spectral entropies[END_REF] generates features using entropy. A scale space decomposition is used to generate three scales of an image, and entropy is calculated for each image block in the spatial and DCT domain. The entropies are then pooled by percentile pooling. The mean and variance of the spatial and frequency components are used in the feature vector.

An alternative approach in LDR NR-IQA is using Convolutional Neural Networks (CNN) as in [START_REF] Kang | Convolutional neural networks for no-reference image quality assessment[END_REF]. We refer to this method as kCNN throughout the paper. This is an enhanced version of [START_REF] Ye | Real-time no-reference image quality assessment based on filter learning[END_REF]. The basic idea is to learn discriminative features that can perform IQA rather than using a handcrafted method; [START_REF] Ye | Real-time no-reference image quality assessment based on filter learning[END_REF] uses dictionary learning to form discriminative filters and [START_REF] Kang | Convolutional neural networks for no-reference image quality assessment[END_REF] improves the process by redesigning it as a CNN. The CNN has four layers that act on MSCN image blocks of size 32×32. The first layer is a convolution layer with 50 filters (kernels), followed by a pooling layer that reduces the dimensionality of the data, and finally there are two fully connected layers. The network is trained with the Mean Opinion Scores (MOS). This method can produce a "perceptual distortion map", which shows location of the errors of the distorted image.

In general, these LDR NR-IQA algorithms rely largely on statistical characteristics of the distorted images, i.e., they assume image distortions alter the statistical properties exhibited by "natural" undistorted images. However, they do not take in to consideration the human visual resistance to errors, e.g., due to masking phenomena. The natural image statistics are to be captured in the internal representation of the SVM or CNN. Thus, these LDR NR-IQA algorithms need an explicit training stage in order to learn what "natural" is and how the noise change the "naturalness" property. A recent work that tries to alleviate this problem is [START_REF] Zhu | A no-reference objective image quality metric based on perceptually weighted local noise[END_REF], which assesses perceived image quality corrupted by uniform and high-frequency noise. The method uses a combination of features, with appropriate feature weights to scale the errors. However, it cannot estimate compression errors. In contrast, we formulate the naturalness property in our model and our technique can also assess compression errors.

B. Tone-Mapped IQA

As explained earlier, Tone-Mapped IQA is a related research area, where HDR images are tone-mapped to LDR images, which are then evaluated on LDR displays. The Tone-Mapped Quality Index (TMQI) metric [START_REF] Yeganeh | Objective quality assessment of tone-mapped images[END_REF] follows the structural fidelity criterion [START_REF] Wang | Multiscale structural similarity for image quality assessment[END_REF], to compare an HDR image with its tone-mapped version, by embedding the knowledge of the Contrast Sensitivity Function (CSF) at different values of luminance [START_REF] Reddy | Perceptually optimized 3d graphics[END_REF], [START_REF] Luebke | Perceptually driven simplification for interactive rendering[END_REF]. In addition, similar to NR-IQA, a naturalness measure is also included to compare the statistics of the tone-mapped image to those of natural images. This idea is further explored in [START_REF] Kundu | Visual attention guided quality assessment of tone-mapped images using scene statistics[END_REF], where the performance is improved with better error pooling and naturalness measure. Phase congruency is added as a feature in [START_REF] Nafchi | FSITM: A feature similarity index for tone-mapped images[END_REF] for the same purpose to compare two images.

A NR-IQA approach for Tone-Mapped HDR [START_REF] Kundu | No-reference image quality assessment for high dynamic range images[END_REF] employ MSCN images as spatial domain features. It also obtains gradient computations on different neighborhoods of every pixel. This is followed by Gaussian parameter extraction and a Support Vector Regressor (SVR) process like the other techniques described in Section II-A. The idea is to generate a HDR image by fusing images captured using multiple exposures. The HDR image is then tone-mapped to a LDR image. The groundtruth Mean Opinion Score (MOS) is obtained based on subjective quality evaluation using LDR displays. The method is statistics-based and does not incorporate perceptual modeling. In contrast, our evaluations are based on displaying HDR images on compatible HDR displays and the human perception component is an integral part of our model. It should be understood that a tone-mapped HDR image has a reduced gamut of colors and luminance compared to the original HDR image. Experimental evidence of psychophysical differences in viewing HDR and LDR image content is provided in [START_REF] Akyüz | Do HDR displays support LDR content?: a psychophysical evaluation[END_REF]. Their study demonstrates how HDR and LDR (or reduced dynamic range) contents are perceived differently when displayed on a HDR screen. The authors collected opinions, including users' ratings of naturalness, visual appeal, spaciousness, and visibility. Here, visibility refers to the details in the image. The study found statistically significant difference in how users rated visibility for HDR and LDR images when these are displayed on a HDR screen.

III. MOTIVATION

Conventional LDR displays have a maximum luminance of about 300 cd/m 2 . High Dynamic Range (HDR) displays have a luminance of 4,000 cd/m 2 and above, which delivers more realistic scenes and vivid content to the Human Visual System (HVS). In addition to advanced HDR acquisition devices, HDR images can be generated using multi-exposure fusion algorithms [START_REF] Shen | Generalized random walks for fusion of multi-exposure images[END_REF] and tone-mapped onto lower dynamic range images [START_REF] Krawczyk | Lightness perception in tone reproduction for high dynamic range images[END_REF], [START_REF] Reinhard | Photographic tone reproduction for digital images[END_REF] for evaluation on LDR displays. Current NR-IQA methods focus on quality assessment on LDR displays. If these methods are applied to predict the quality of HDR images displayed on a HDR display, the result is not accurate, as shown in Section V. This is because they rely on statistical modeling of noise and fail to take into account how the Human Visual System (HVS) responses to HDR displays [START_REF] Aydın | Extending quality metrics to full luminance range images[END_REF]. Therefore, to design a robust HDR NR-IQA system, which is consistent with how the HVS perceives real-world content, we need a new NR-IQA model that incorporates perceptual factors to predict the visibility of error on HDR screens. To date, we are not aware of any related work that is designed for HDR viewing conditions.

A. HVS response to HDR displays

How the HVS responses differently to HDR and LDR displays has been extensively studied by Aydin et al. [START_REF] Aydın | Extending quality metrics to full luminance range images[END_REF].

In their experiments, the authors used a LDR display with luminance range 1-100 cd/m 2 and another HDR display with luminance range 10-1,000 cd/m 2 to evaluate the difference in perceived image qualities. They found that a distorted image, when viewed on the brighter display, was perceived as worse compared to the same image displayed on the display with lower brightness. This shows that viewing content on LDR displays and HDR displays has different perceptual effects.

Hence, we argue that, in addition to statistical comparison, HDR NR-IQA should incorporate the psychophysical phenomena that can determine the perception of distortion in HDR conditions. Although such perceptual approach is not common in NR-IQA, it is often used in FR-IQA models. HDR-VDP-2.2 [START_REF] Narwaria | HDR-VDP-2.2: a calibrated method for objective quality prediction of highdynamic range and standard images[END_REF] and HDR-VQM [START_REF] Narwaria | HDR-VQM: An objective quality measure for high dynamic range video[END_REF] are two examples. HDR-VDP 2.2, simulates the early processing stages of the HVS, based on psychophysical measurements. HDR-VQM uses sub-band decomposition and spatial-temporal error pooling to simulate visual recognition. Both are FR-IQA methods, which generate a local error visibility map. An alternative approach is to apply SSIM [START_REF] Wang | Image quality assessment: From error visibility to structural similarity[END_REF], which is a FR-IQA method for PU-encoded HDR data. These are the best performance algorithms in FR-IQA for HDR data based on the survey presented in [START_REF] Zerman | An extensive performance evaluation of full-reference HDR image quality metrics[END_REF]. We will compare the performance of our proposed NR-IQA system with these FR-IQA.

IV. PROPOSED SYSTEM

In order to design a perceptual model consistent with how the HVS perceives HDR real-world content, it is necessary to understand the contrast sensitivity associated with a complex image. The traditional approach considers various visual features like contrast, frequency and background luminance. Very often, a handcrafted function is used to compute a quality score based on a combination of these features. The quality score function is obtained by fitting opinions collected from subjects in psychophysical experiments using sample datasets. The research question is how to generalize user study results for all real-world images.

We approach this problem of designing a perceptual HDR NR-IQA model by dividing the visual quality analytic process into sub-components. We represent visual quality perception as the result of two functional units. The first unit takes a distorted image and detects error, and the second unit performs a perceptual scaling of this error to compute a quality score. By using a supervised learning approach, the mathematical behavior of these two units can be modeled. The data required for this training is obtained from an IQA dataset, which contains images and the corresponding quality scores.

A. Model Overview

To formulate the above idea, we design a Convolutional Neural Network (CNN) (Fig. 1) that processes HDR image blocks composed of linear luminance values. We use a block size of 32x32 pixels. This is the same block size that was suggested in [START_REF] Kang | Simultaneous estimation of image quality and distortion via multi-task convolutional neural networks[END_REF]. Our CNN model has three major parts: Enet, P-net and a Mixing function. E-net estimates the Error δ(i, j) of an image block centered at (i, j). P-net computes the Perceptual Resistance T (i, j) of the block. The output of these two systems are then input to a Mixing function, to produce the local block quality. We obtain Differential Mean Opinion Scores (DMOS) for each image block. The block scores are then combined to generate the final image quality score. In our model, DMOS is a number directly proportional to the level of distortion in a HDR image.

B. E-net Error Estimation

The Error δ(i, j) quantifies the change in statistics in a distorted image block. For an image block centered at (i, j), we define the error as,

δ(i, j) = mean( Y R (i, j) -Y D (i, j) ) (2)
where Y R and Y D are, respectively, the original and distorted linear HDR luminance values of the image block centered at (i, j). This does not indicate a Full Reference computation, as the original version is only used during training (a pair of HDR image and its distorted version is used here). The objective is to train E-net with the distortion characteristics, like blocky artifacts, blurring effects, jagged edges, etc. We use 1 norm for error computation (Eq. 2) instead of alternative measures such as 2 norm, to avoid overemphasizing large errors. This is particularly important in HDR images where the histogram of Y is generally very skewed and some pixels take very high luminance values. We use our own CNN approach to design E-net to obtain and estimation δ(i, j) of the error in Eq. ( 2).

C. P-net Perceptual Resistance

For each image block centered at (i, j), we compute the Perceptual Resistance T (i, j). This value represents the difficulty for a viewer to perceive the error δ(i, j) of the block. A high T value implies that it is less likely to see the error, and hence the quality of the block is less affected (high perceptual resistance). Conversely, a low value implies that the image block will be perceptually degraded by error.

Perceptual Resistance T (i, j) aims to represent a combination of all perceptual effects exhibited by an image block centered at (i, j). Though it is functionally similar to the pixel-wise just noticeable error measure used in conventional IQA systems like [START_REF] Zhu | A no-reference objective image quality metric based on perceptually weighted local noise[END_REF], [START_REF] Daly | Visible differences predictor: an algorithm for the assessment of image fidelity[END_REF] and [START_REF] Narwaria | HDR-VDP-2.2: a calibrated method for objective quality prediction of highdynamic range and standard images[END_REF], we introduce Perceptual Resistance as a new term because our model generates local quality scores (DMOS), as opposed to a local probability of error detection.

Instead of following the traditional perceptual modeling method of deriving perceptual thresholds from psychophysical experiments, we solve this problem by a data driven method. We use a convolutional neural network (CNN) based architecture, P-net, to derive the Perceptual Resistance of a block.

The CNN automatically computes the features required to do this task by a training process using real-world images. A detailed analysis of how a generic CNN generates its features is explained in [START_REF] Zeiler | Visualizing and understanding convolutional networks[END_REF]. In our system, P-net approximates a function that maps the image block values onto a perceptual resistance value. Differing from a conventional psychophysical perceptual function, the 'function' that is captured by P-net is represented as weights in the CNN.

D. Mixing function

We use a Mixing function f ( δ, T ), which combines the estimated error and Perceptual Resistance to generate a quality score. This is a critical part of the system because it is this value that is optimized by the training process to match human quality scores. The output of P-net would change based on how the Mixing function is designed.

The Mixing function is designed as follows, with error expressed in multiples of Perceptual Resistance:

DM OS = f ( δ, T ) = G δ T , ( 3
)
where G is a monotonically increasing function. By using this, we express error in JND like measure ( error JN D ), so that the error quantity is mapped onto a more perceptually relevant scale. Such interpretation is common in IQA literature, e.g., [START_REF] Zhu | A no-reference objective image quality metric based on perceptually weighted local noise[END_REF], [START_REF] Aydın | Extending quality metrics to full luminance range images[END_REF] and [START_REF] Narwaria | HDR-VDP-2.2: a calibrated method for objective quality prediction of highdynamic range and standard images[END_REF].

Mapping δ T to quality scores is achieved by the function G. Since increase in visible error always leads to decrease in quality and increase in DMOS, the latter must monotonically increase with error JN D , implying that G also has to be monotonically increasing with δ T . Thus, any monotonically increasing function is sufficient for G. However, choosing a G that is too complex can lead to optimization problems because of unstable data points along the function, or low values for gradients, leading to slow or zero learning. We do not go into the mathematics of CNN convergence and optimization functions as it is beyond the scope of this work.

Based on the above considerations, we use G(x) = 1exp(-|kx|) and DMOS is computed as:

DM OS(i, j) = 1 -exp - k * δ(i, j) T (i, j) (4)
This choice is inspired by the error model proposed in [START_REF] Zhu | A no-reference objective image quality metric based on perceptually weighted local noise[END_REF], but we introduce a scaling factor k. Here the added parameter (weight) k can be tuned during the training process, so that the predicted values of DMOS are as close to the ground truth DMOS as possible.

As seen from the plot of Eq. 4 in Fig. 2, G can characterize different rates of increase for DMOS, depending on the values of δ(i, j), T (i, j) and k. Note how k serves to control how slowly DMOS increases with δ T (Fig. 2 B). The block-wise DMOS scores obtained from Eq. ( 4) is converted to a global score by averaging. A weighted scheme is not required here since the perceptual scaling of errors based on content is already handled by Eq. ( 4) (the T computed by the CNN changes based on image content and handles contentdependent scaling). For P-Net, the ideal training data is a numeric quantity, encapsulating all perceptual effects on the HVS, generated from an image block. Although we cannot get such a final value directly, our system can produce a quality score after the Mixing function process. We use this score for training. This two-stage training forces the P-net to extract a set of perceptual features from the image blocks and to derive a single final Perceptual Resistance value. The optimal value of k in the mixing function is also computed.

We therefore define our two-stage training process as follows:

Stage 1: E-net is trained with distorted image blocks as input and errors δ as target. The error is computed with Eq. (2).

Stage 2: all the training weights of E-net are frozen by setting their learning rate to zero. The whole network is then trained with image blocks as input and ground truth image quality of the whole image, DM OS gt , as target. We use J as the cost function for any image block centered at (i, j), where

J(i, j) = |DM OS(i, j) -DM OS gt | . ( 5
)
DM OS(i, j) is the output of the Mixing function. The P-net and the mixing function (optimal value of k) get trained during this stage. The overall process is illustrated in Fig 3.

Notice that in Eq. ( 5) we assume that the local quality of an image block is the same as the global image quality score, similarly to the setting in [START_REF] Kang | Simultaneous estimation of image quality and distortion via multi-task convolutional neural networks[END_REF]. While this assumption is somehow inaccurate (as distortion can be unevenly spread across a picture), it has been proven to be accurate enough to predict image quality without reference [START_REF] Kang | Simultaneous estimation of image quality and distortion via multi-task convolutional neural networks[END_REF]. We leave to future work how to produce more accurate groundtruth quality labels for training, e.g., as in the very recent method [START_REF] Heng | From Image Quality to Patch Quality: An Image-Patch Model for No-Reference Image Quality[END_REF].

F. Network Architecture

Our network architecture has two major advantages. First, a separation of the perceptual component from the physical error gives more intuitive results that can be used in applied fields of IQA like image and video compression. With adequate calibration, Perceptual Resistance values can be used to optimize compression or transmission. Secondly, it simplifies the learning process leading to improved results and better generalization of those results to real-world cases.

E-net is a typical CNN architecture consisting of five layers. The E-net structure is shown in Fig. 4. Spatial pooling is used after each filtering stage. The final layer consists of one node corresponding to the output. Spatial dropout layers [START_REF] Tompson | Efficient object localization using convolutional networks[END_REF] are added to prevent over-fitting of data.

For P-net, we define a customized CNN layer, namely Augmented Input Layer. In this layer, in addition to the original luminance values of the block, we compute the mean, variance and MSCN images. For the MSCN image, we use the formulation in [START_REF] Mittal | No-reference image quality assessment in the spatial domain[END_REF], i.e., M SCN (y N (i, j)) = y N (i,j)-µ y N (i,j) σ y N (i,j) +0.01 , where µ y N (i,j) is the mean and σ y N (i,j) is the variance. They are computed by replacing every pixel (i, j) with the mean and variance, respectively, over a local Gaussian window of size N around (i, j). We use a smaller value 0.01 as the stabilizing constant in order not to significantly impact the MSCN values. Since neural network training requires the input values in a specific range, we scale the input, variance map and the MSCN map with trainable scalar weights, w i , i = 1, . . . 4, whose values are determined as part of the overall optimization process. Hence the output of the Augmented Input Layer are four scaled feature maps:

w 1 • y N (i, j), w 2 • µ y N (i,j) , w 3 • σ y N (i,j) , w 4 • M SCN (y N (i, j)).
The subsequent layers are convolutional and fully connected layers. The final layer has one node corresponding to the output. The P-net structure is shown in Fig. 5.

The results of E-net and P-net are integrated and analyzed by the Mixing function, whose behavior is modeled by Eq. 4, with trainable parameter k tuned during the training process.

V. EXPERIMENTS AND RESULTS

Fig. 4: E-net structure for error estimation. In this section, we compare the performance of our algorithm with existing methods and show a clear improvement in performance. We conduct two tests: 1) a test of overall performance of the proposed method on a large dataset of subjectively annotated HDR images; 2) a cross-dataset test to assess the generalization capabilities of the proposed approach.

A. Datasets

To obtain a large enough number of HDR images, we combine subjective scores from five different datasets [START_REF] Narwaria | Tone mapping-based high-dynamic-range image compression: study of optimization criterion and perceptual quality[END_REF], [START_REF] Narwaria | Impact of tone mapping in high dynamic range image compression[END_REF], [START_REF] Korshunov | Subjective quality assessment database of HDR images compressed with JPEG XT[END_REF], [START_REF] Valenzise | Performance evaluation of objective quality metrics for HDR image compression[END_REF] and [START_REF] Zerman | An extensive performance evaluation of full-reference HDR image quality metrics[END_REF]. The authors of the respective datasets conducted subjective evaluations using different HDR displays, of different maximum luminance and viewing distances. The characteristics of the datasets are described in Table I. In order to align the subjective scores of the different datasets, we follow the same setup as the evaluation of HDR FR-IQA in [START_REF] Zerman | An extensive performance evaluation of full-reference HDR image quality metrics[END_REF], i.e., we employ the iterated nested least square algorithm (INLSA) proposed by Pinson and Wolf [START_REF] Pinson | An objective method for combining multiple subjective data sets[END_REF]. Details about the application of this method for the considered datasets are given in [START_REF] Zerman | An extensive performance evaluation of full-reference HDR image quality metrics[END_REF]. All the experiments in the following are done using the aligned scores obtained after INLSA. The datasets provide only MOS values of the images. Since our system requires the difference of mean opinion scores (DMOS), we convert MOS to DMOS as follows:

DM OS gt (i) = M OS M AX -M OS(i) M OS M AX , ( 6
)
where DM OS gt (i) is the ground truth DMOS score for image i, M OS M AX represents the maximum MOS in the IQA training dataset and M OS(i) is the MOS of the i th image of combined database after INLSA alignment.

B. Experimental setup

The proposed system was implemented on a computer with an Intel core i7 processor, 16GB RAM, and a Nvidia GTX1070 graphics processor. Python was used as the programming language with Keras on Theano backend, Imageio and Open CV as supporting libraries. We used Adam optimizer ( [START_REF] Kingma | Adam: A method for stochastic optimization[END_REF]) to optimize the weights of the CNN. The parameter values of Adam was learning rate = 0.001, β 1 = 0.9, β 2 = 0.999, = 10 -8 and decay = 0.0. The batch size used was 200. The training was done for 10 epochs.

Performance comparison measures used were Spearman Rank order Correlation Coefficient (SRCC), Kendall Rankorder Correlation Coefficient (KRCC), Pearson Linear Correlation Coefficient (PLCC) and Root Mean Square Error (RMSE). A good NR-IQA is characterized by a higher value for SRCC, KRCC, and PLCC, and a lower value for RMSE. The metric values and subjective DMOS were scaled to [0,100] for evaluation of RMSE.

C. IQA Reference Schemes

Since the research on HDR NR-IQA is still at its preliminary stage and there is no generally accepted benchmark metric, we compared our approach with a number of state-of-theart LDR NR-IQA methods: BRISQUE [START_REF] Mittal | No-reference image quality assessment in the spatial domain[END_REF], SSEQ [START_REF] Liu | No-reference image quality assessment based on spatial and spectral entropies[END_REF], BIQI [START_REF] Moorthy | A two-step framework for constructing blind image quality indices[END_REF], DIIVINE [START_REF] Moorthy | Blind image quality assessment: From natural scene statistics to perceptual quality[END_REF], and kCNN [START_REF] Kang | Convolutional neural networks for no-reference image quality assessment[END_REF], with and without preprocessing operators. The results were obtained after retraining the algorithms on the respective processed HDR datasets. PU encoding was used as a pre-processing operator because it has been shown to perform well in a similar context for HDR FR-IQA [START_REF] Valenzise | Performance evaluation of objective quality metrics for HDR image compression[END_REF]. We also used a number of tone-mapping operators, which include [START_REF] Reinhard | Photographic tone reproduction for digital images[END_REF], [START_REF] Reinhard | Dynamic range reduction inspired by photoreceptor physiology[END_REF], [START_REF] Drago | Adaptive logarithmic mapping for displaying high contrast scenes[END_REF] and [START_REF] Mantiuk | Display adaptive tone mapping[END_REF], in pre-processing. Features were extracted using the implementation provided by the authors. In the case of SSEQ [START_REF] Liu | No-reference image quality assessment based on spatial and spectral entropies[END_REF], we normalized the images by the maximum intensity under the respective schemes (4000 for linear HDR and 455 for PU encoded data). For training the SVM, a grid search on the cost and kernel parameters of the SVM was conducted for a range of 10 -15 to 10 15 before training. The algorithm in [START_REF] Hou | Blind image quality assessment via deep learning[END_REF] was re-implemented using Python.

D. Overall performance

In order to test the overall performance of the proposed method, we (re)train each algorithm on our combined image dataset, by splitting it in training/testing subsets (80% for training and 20% for testing). We repeat this procedure 100 times, similar to [START_REF] Hou | Blind image quality assessment via deep learning[END_REF]. We assure that the training and testing sets do not contain the same contents. Note that the results can vary slightly since the CNN weights initialization is random. The computed median SRCC, KLCC, PLCC and RMSE are shown in Table II.

Considering the NR-IQA originally designed for LDR content, we see an acceptable performance in SRCC (about 0.7 after retraining) for many of the algorithms. Best performances are obtained by using BRISQUE [START_REF] Mittal | No-reference image quality assessment in the spatial domain[END_REF] and kCNN [START_REF] Kang | Simultaneous estimation of image quality and distortion via multi-task convolutional neural networks[END_REF]. The high performances of BRISQUE and kCNN can be attributed to the features they use, i.e., the MSCN coefficients. It is likely that the normalization by variance cancels the effects of the increased dynamic range and yields a similar distortion pattern as LDR images. Practically, kCNN is more useful because it produces a perceptual distortion map in addition to the quality score. The perceptual distortion map indicates what errors a human viewer can notice on the noisy image. Furthermore, we observe a clear performance improvement in LDR NR-IQA algorithms if the data is pre-processed and the dynamic range of the data is reduced to LDR levels. PU encoding improves the performance in most of the cases. The best performance among LDR NR-IQA is obtained while using PU encoding in conjunction with kCNN. The performance of the proposed system is significantly better than the other algorithms in all cases both with or without preprocessing using PU encoding.

E. Cross-dataset testing

In order to demonstrate the generalization capabilities of the proposed NR-IQA technique to different conditions and contents, we train the algorithms using datasets #1, #2 and #3, and test them on datasets #4 and #5. This represents a realworld scenario, where the validating conditions are different from that of the training data. In addition, this testing method allows us to compare performance with FR-IQA algorithms. From a machine learning point of view, we have a sufficient number of examples of each type of distortion in datasets #1, #2 and #3. There is also a combination of all distortion types in dataset #4 and #5. The test set contains DMOS scores uniformly distributed in the range [START_REF] Akyüz | Do HDR displays support LDR content?: a psychophysical evaluation[END_REF]90].

Since the CNN is initialized with a set of randomly generated weights, the results of training can vary. We report the median score after 10 train-test cycles. Our results are given in Table III. By itself, BRISQUE, BIQI, SSEQ and DIIVIINE seem to be unable to adapt to the different image sizes and luminance ranges in the testing set, when these features are different from the training set. This can be explained by the fact that the features used by these algorithms are computed over a global histogram from the entire image.

The kCNN method performs well and shows good adaptability to a different test image size. This can be attributed to the fact that an image block is used to train the kCNN and hence the overall image size becomes less important. The method, however, does not take into account perceptual factors and we can see an improvement if PU encoding is used in preprocessing.

Our proposed algorithm outperforms related methods in all test cases when considering generalization to a real-world scenario. The superior performance demonstrates the strength of our two-stage training design, which successfully adapts to different image and luminance range. Perceptual Resistance values are able to scale the errors based on the luminance and image content. Our result achieves performance close to FR-IQA methods without the need of a reference image. For the full reference methods, the metric values are scaled to the range of [0,100] as suggested in [START_REF] Zerman | An extensive performance evaluation of full-reference HDR image quality metrics[END_REF]. A scatter plot comparing the scores generated by the proposed method and the human judged DMOS is shown in Fig. 6. A high correlation is observed.

Fig. 6: Scatter plot between objective scores from the proposed method and subjective human judged DMOS.

VI. ANALYSIS AND DISCUSSION

In this section, we analyze each sub-component of the proposed system and show the functions performed by them in detail. We highlight some interesting results provided by our system in addition to the quality score of an image. We also discuss the weakness and failure cases. In the following, we color-code the results with red points showing high values, green intermediate and blue low values.

A. Error Estimation

Per block error estimation is performed by E-net (Figure 1). The output δ of E-net on a few example images from Datasets #4 and #5 is reported in Fig. 7, which shows the distorted images, the actual error δ in the image (equation 2) and the error estimations generated by E-net. As a reference, we also include the corresponding Root Mean Square Error (RMSE) between δ and δ in Fig 7 (e).

A more detailed look at the error estimations can be seen in Figure 8, where we show enlarged regions of the image containing actual error and the corresponding error estimations. It is apparent from the error maps that E-net is able to successfully identify major errors in the images.

B. Perceptual Resistance

The Perceptual Resistance value T produced by P-Net shows how unlikely it is that an error in a region can reduce local visual quality. We use the image in Fig. 9 to illustrate how P-net handles contrast sensitivity and masking effects. We consider a JPEG compressed image from the CSIQ dataset [START_REF] Larson | Categorical image quality (CSIQ) database[END_REF]. The image is scaled to a maximum luminance of 100 cd/m 2 . This is done for ease of explanation as this is the luminance range commonly seen in everyday computer and television screens. The image is input to the proposed algorithm as blocks of size 32 × 32. The output from E-net, P-net and the Mixing function are shown in Fig. 9.

a) Analysis of spatial masking: Spatial masking effects can be observed in the output of P-net ( T in Fig. 9). In the color map, the Perceptual Resistance T is lower for pixels in the sky, which is brighter and smoother, compared to other image blocks. Thus, errors in the sky region are easy to notice. On the contrary, in regions with high spatial frequency, e.g., bushes and architecture, T is higher indicating that errors are less noticeable. The error δ estimated by E-net and the actual error δ is shown in Fig. 9. The perceptual distortion map obtained after the mixing function reflects high values in the sky and low values in the more densely textured regions.

b) Analysis of sensitivity to luminance: As explained in Section III.A, [START_REF] Aydın | Extending quality metrics to full luminance range images[END_REF] performed a user study by displaying a distorted image on two screens with maximum luminance of 100 cd/m 2 and 1,000 cd/m 2 , respectively. The study revealed that users rated the perceived quality on the high luminance monitor worse compared to the the same image displayed on the lower luminance monitor. We reproduced this finding by applying our algorithm on the LDR images in the CSIQ database [START_REF] Larson | Categorical image quality (CSIQ) database[END_REF]. We used all the JPEG distorted images in the database and linearly scaled the images to luminance ranges of 100 cd/m 2 and 1,000 cd/m 2 , respectively. We then computed the quality of these images using the proposed method and compared the scores of the two luminance ranges. We found that the mean value of DMOS for the images of low brightness was 68.76 and those of high brightness was 76.72. This implies that the brighter images appear more distorted than their corresponding low brightness counterparts. We then employed ANalysis Of VAriance (ANOVA) to check for statistical significance of the difference between DMOS for the low and high brightness images. The p value was found to be < 0.05 and F < F crit for the hypothesis that the difference in means was zero; proving that the difference is statistically significant. Hence, we can reproduce the results obtained in [START_REF] Aydın | Extending quality metrics to full luminance range images[END_REF] using our data driven method without the need of conducting low level visibility user studies.

We performed further analysis on the image discussed in Fig. 9. Fig. 10 shows the output of our proposed method at different noise levels under 2 different luminance ranges. Notice that:

1) In both luminance ranges, the estimated δ is the same, because the measure is based on content statistics and does not take into account the change in luminance.

2) The estimated T is instead affected by the luminance range, as discussed in Section III. For a given distortion level, the perceptual resistance is lower for higher luminance values, consistently with the the findings in [START_REF] Aydın | Extending quality metrics to full luminance range images[END_REF]. However, the mapping of luminance to perceptual resistance is not a simple global linear scaling, but takes into account the complex content-dependent characteristics that determine contrast and luminance masking. c) Comparison between perceptual resistance and contrast sensitivity: HDR-VDP 2.2 can provide a contrast threshold, i.e., the per pixel contrast (difference) such that an error is visible with a certain probability. This map is related to our Perceptual Resistance values. However, the contrast threshold is the result of psychophysical experiments to determine contrast sensitivity at different luminance levels, while Perceptual Resistance is indirectly obtained through quality scores.

A visual comparison, Fig. 13, can give and idea of the relative values of the two measures. The color maps shown are generated using the same scale and value range, i.e., log scale normalized to the range [0,1]. The comparison shows

C. Comparison using alternative Mixing functions

As explained in Section 4.D, an alternative Mixing function can be used as long as it is monotonically increasing and the network can converge.

We tested the following cases to study the effects of different mixing function formulations:

1) Proposed: DM OS(i, j) = 1 -exp(-k * δ(i,j)

T (i,j) ). 2) Linear mixing: DM OS = δ T . 3) Hyperbolic tangent: DM OS = tanh( δ T ). 4) Logistic function:

DM OS(i, j) = 1 1 + exp(-k(x -x 0 )
where x = δ T . We found that during the train-test cycles in cases 2 and 4 above, the network failed to converge and the training error kept increasing. This happens when the random weight initialization causes these functions go out of bound, interfering the gradient propagation. Cases 1 and 3 do not have convergence problem. The results are shown in Table IV. We observed a similar performance whenever the network converged, implying the possibility of having multiple choices for G.

To investigate further, we selected a distorted image and generate the results using different mixing functions and P-net as shown in Fig. 11. As explained earlier, a mixing function outputs the DMOS and P-net outputs the Perceptual Resistance of an image block. A comparison of HDR-VDP 2.2 error probability and contrast threshold is shown for reference.

The perceptual distortion maps do not change significantly as the mixing function changes. This is because the CNN optimization process tries to minimize the difference between the mixing function output and the actual human judged DMOS. Upon convergence, the results will be similar.

Considering P-net, we see that the results show a similar trend in terms of relative values. For example, visibility values in the sky are generally higher compared to the texturerich forest area. However, the level of change in Perceptual Resistance varies depending on the type of mixing function used. The proposed Mixing function defined in Eq. 4 is used in our system because it generates the closest results to that of the FR-IQA metric HDR-VDP2.2.

Finally, the Mixing function could automatically be learned from data using another CNN. However, this process would involve more weights and difficulties in optimization, with no guarantee that the overall model converge to a good solution. We confirmed this experimentally by using a Deep Belief Network in place of the Mixing function. Even if the system did converge with a DBN, the gain in performance with respect to the proposed mixing function is not substantial, and this function would be rather a "black box" with no intuitive interpretation.

D. Perceptual distortion map

One of the advantages of the proposed CNN based NR-IQA scheme is that it gives the approximate location of the perceived errors in a perceptual distortion map. A comparison of the perceptual distortion maps produced by FR-IQA (HDR-VDP 2.2, PU-SSIM) and NR-IQA (PU-kCNN) algorithms and the proposed method for some example images is shown in Fig. 12. In order to highlight the errors, the probability of error detection is shown for HDR-VDP 2.2, and an inverted PU-SSIM map is shown for PU-SSIM. Notice that the map obtained using HDR-VDP 2.2 is conceptually different from those of PU-SSIM and PU-kCNN: the former represents per pixel probability of detecting distortion, while the latter two convey information about the magnitude of local distortion. Instead, our map summarizes the two types of information into a local perceptual distortion, which can be seen as the inverse of local perceptual quality.

Although they express different properties of distortion, the perceptual distortion maps produced by the proposed scheme, PU-kCNN and HDR-VDP2.2 are consistent in general. To show the performance improvement with respect to the current state-of-the-art PU-kCNN to produce a perceptual distortion map, we compute the MSE between perceptual distortion maps of proposed and PU-kCNN with respect to HDR-VDP2.2, reported in Fig. 12. From both MSE and visual inspection, it is clear that our perceptual distortion map is closer to HDR-VDP visibility map, compared to PU-KCNN, i.e., it discounts errors that are less likely to be perceived due to masking and reduced sensitivity to contrast. This explains the overall better performance of the proposed approach over kCNN in Table III.

E. Effect of preprocessing

In HDR FR-IQA it is common to preprocess images through a perceptual transfer function [START_REF] Aydın | Extending quality metrics to full luminance range images[END_REF], [START_REF] Zerman | An extensive performance evaluation of full-reference HDR image quality metrics[END_REF], [START_REF] Miller | Perceptual signal coding for more efficient usage of bit codes[END_REF] in order to compute perceptually meaningful differences between pixels. However, in our approach we decouple the error computation from its perceptual scaling, and this preprocessing turns out to be unnecessary, as it is implicitly carried out by the P-net.

To support this claim, we test the performance of our system by preprocessing images with the PU encoding [START_REF] Aydın | Extending quality metrics to full luminance range images[END_REF] before training the network. The performance with this setting is reported at the bottom of Table II. Although PU encoding improves the performance of all other NR-IQA (since it provides perceptual scaling), it does not improve the performance of the proposed system -in fact, correlations are slightly lower if pixels are PU encoded. This could be attributed to the loss of information due to perceptual encoding, where information is somehow quantized and some slight variations in the data are compressed with a fixed scaling function.

F. Failure cases

We isolated some cases where the predicted DMOS quality has a large error compared to the groundtruth. In some of those cases, the perceptual distortion maps produced by our method are poorly correlated with the visibility maps produced by HDR-VDP or with the local PU-SSIM estimated distortion. Some examples are reported in Fig. 14. We observe that the perceptual distortion maps produced by our method are not consistent with the ones produced by HDR-VDP and PU-SSIM (e.g., in the sky for the first image, or in the lake for the second or the dress of the woman in Fig. 12): in one case, the FR algorithms estimate minimal distortion, while our method predicts higher errors, while in the other two cases, our method underestimates perceptual distortion. MSE between the estimated distortion map and HDR-VDP error visibility, as well as groundtruth DMOS values, are given for reference.

To analyze this further, we note that the discrepancy appears mainly in large smooth areas, e.g., sky and the water of the lake, where HDR-VDP shows high error detection, and in cluttered ground area where HDR-VDP shows low error detection in smooth region. PU-SSIM estimates a similar spatial distribution of the error into both regions. Another instance of this is in Fig. 12, second row, where the highly textured dress of the woman is shown to have high errors in HDR-VDP and PU-SSIM, but the proposed method suggests a low error.

We can see from intermediate results of the proposed method inside blue box of Fig. 14 that the errors are caused by different components of the system. Over-estimation of errors in the smooth sky is the effect of over-estimation of errors by E-net in smooth areas, whereas under-estimation of errors in highly textured regions is the result of perceptual resistance being too high. These errors are likely due to some bias in our dataset, where a large fraction of examples contain smooth or textured regions have low and high quality respectively.

We notice, however, that those results could indeed be meaningful from a perception point of view, and that the visibility map produced by HDR-VDP alone is not a good indicator of quality [START_REF] Valenzise | Performance evaluation of objective quality metrics for HDR image compression[END_REF]. In order to further assess the prediction of local perceptual error and its impact on overall quality, we would need per block groundtruth of quality scores, which not only is unavailable nowadays, but is very difficult to produce as assessing quality is by definition a global task. Furthermore, the proposed method estimates quality per block independently, which is a major simplification in the model. Extensions to how to embed higher order dependencies between regions of the image, and possibly semantic considerations (which become realistic using deep CNN architectures), are left to future work.

G. Computational complexity

Assuming a fully trained system, the computational complexity of our method is high because of the large amount of convolutions involved. The asymptoptic time complexity is mainly due to E-net because of its convolutional layers. Each convolutional layer uses the results of the previous layer. This makes our method slower compared to kCNN where there is only one layer of convolutional filters. We found that in our GPU implementation, computation of score for an image with 1920 × 1080 resolution takes 0.33 seconds (assuming the theano graph is in memory). For KCNN, a similar coding style gave us execution time of 0.15 seconds per image. Note that the execution speed can be improved with better implementation. This is left as a future work; our focus in this work is to improve NR-IQA performance.

VII. CONCLUSION

We propose a HDR NR-IQA scheme that uses a CNN based architecture, composed of E-net, P-net and Mixing function, to generate values corresponding to Error Estimation and Perceptual Resistance in an image, which are then combined to generate a perceptual distortion map and DMOS. Taking into account perceptual phenomena is directly derived from real-world data driven optimization, and does not involve psychophysical user studies. The derived Perceptual Resistance shows similar characteristics compared with other perceptual models. Experimental results demonstrate that our algorithm accurately predicts visual distortions such as compression artifacts. Our algorithm scores correlate well with human scores. It outperforms other NR-IQA methods and the performance is competitive compared to HDR FR-IQA methods, with the advantage that no reference image is needed.

  N.K. Kottayil and I. Cheng are with Dept of Computing Science, University of Alberta, Edmonton Canada G. Valenzise and F. Dufaux are with the Laboratoire des Signaux et Systèmes (L2S), CNRS, CentraleSupélec, Université Paris-Sud.

Fig. 1 :

 1 Fig. 1: Proposed strategy for a robust HDR NR-IQA. Enet detects the error, P-net detects the perceptual resistance, and Mixing function consolidates the results and computes a DMOS.

Fig. 2 :

 2 Fig. 2: Behavior of the Mixing function. (A) Varying T with k fixed at 1 illustrates the impacts of changing δ and T on DMOS. (B) Varying k with T fixed at 20 illustrates the impacts of changing k and δ values on DMOS. An optimal k value is determined during the training process.

Fig. 3 :

 3 Fig. 3: Two-stage training process.

Fig. 5 :

 5 Fig. 5: P-net structure for perceptual resistance computation.

Fig. 7 :

 7 Fig. 7: Error estimation. (a) reference images; (b) distorted images; (c) ground truth error map provided by the database; (d) error maps estimated by E-net; (e) RMSE between the estimated and ground truth errors on the HDR luminance images.

Fig. 8 :

 8 Fig. 8: Error estimation. Detailed comparison of actual and estimated error. Columns (a) to (d) have the same meaning as Figure 7.

Fig. 9 :

 9 Fig. 9: Perceptual resistance. Example output of E-net and Pnet. For the content of 'Image', Actual error δ, Output of E-net δ -Estimated error, Output of P-net T -Perceptual Resistance and perceptual distortion map DM OS(i, j) from the Mixing function is shown.

Fig. 10 :

 10 Fig. 10: Perceptual resistance: Behavior of proposed system at different luminance levels. Distortion Level of JPEG is shown in dB, δ is the Acutal error, δ is the Estimated error by E-net, T is the Perceptual Resistance by P-net. Within each column, the left and right images are output when the input image is linearly scaled to 1,000 and 100 cd/m 2 respectively.

Fig. 11 :

 11 Fig. 11: Effect of different mixing functions. Comparison of perceptual distortion maps and perceptual resistance generated from different mixing functions.HDR-VDP error probability and contrast threshold are shown for reference.

Fig. 12 :

 12 Fig. 12: Perceptual distortion maps. Comparison of perceptual distortion maps estimated by various IQA schemes. MSE shown is computed between the distortion maps of proposed method and PU-kCNN with HDRVDP2.2.

Fig. 13 :

 13 Fig. 13: Perceptual distortion map. Comparison between P-net Perceptual Resistance and HDR-VDP2.2 contrast threshold. Image pixel values are in log scale and normalized to [0,1].

Fig. 14 :

 14 Fig. 14: Failure cases. There are cases where the perceptual distortion maps from different perceptual IQA show inconsistency, but there is no existing benchmark metric to evaluate distortion maps. Images inside the blue box are intermediate results.

TABLE I :

 I Datasets characteristics

TABLE III :

 III Cross-dataset results for different IQA methods.

TABLE IV :

 IV Performance with various mixing functions similar perceptual results of the two algorithms. The color maps show that the location and relative intensity of the visual errors are similar. While the contrast threshold map is perceptually more accurate, it requires the knowledge of the pristine image. The proposed two-stage network architecture enables to approximate it without any reference.

VIII. ACKNOWLEDGMENTS

We thank Dr Anup Basu for his insightful comments and valuable feedback and Emin Zerman for his help with the HDR datasets.


Nada Chamakhi

Ines Bouzidi

Azza Ouled Zaid
email: azza.ouledzaid@isi.rnu.tn

Frédéric Dufaux
email: frederic.dufaux@l2s.centralesupelec.fr

JPEG based Compression of Digital Holograms

Keywords: Holography, lossy compression, JPEG, angular spectrum reconstruction

Modern holography for 3D imaging allows to reconstruct all the parallaxes that are needed for a truly immersive visualisation. Nevertheless, it represents a huge amount of data which induces higher transmission and storage requirements. To gain more popularity and acceptance, digital holography demands development of efficient coding schemes that provide significant data compression at low computation cost. Another issue that needs to be tackled when designing holography coding algorithms is interoperability with commonly used formats. The upcoming JPEG Pleno standard aims to develop a standard framework for the representation and exchange of new imaging modalities, such as holographic imaging, while maintaining backward compatibility with legacy JPEG decoders. This paper presents a lossy compression method for holographic images that exhibits good coding performance while considering the computation cost and backward compatibility with legacy JPEG standard. To validate our findings, the results of our tests are shown and interpreted.

I. INTRODUCTION

Nowadays, static planar-like images are giving place to more informative volume-like images, thereby providing viewers an immersive experience. Among these volume-like images, the concept of stereoscopic representation is widely used in movie theatres and high end televisions. But this 3D visualization technology only exploits limited depth cues and causes visual fatigue due to the inherent accommodationvergence conflict. These limitations are surpassed when using holographic systems. Indeed, holographic visualisation guarantees that the convergence and focus points are the same as in real life, which makes it possible to have the most realistic depth perception. Thus, holography can be considered as a unique imaging technique which promises to change the way 3D visual information will be generated and consumed in the future, defeating all limitations and constraints of current 3D technologies.

The optical holography has been initiated in 1947 by Denis Gabor to eliminate aberrations in electron microscopy by recording and reconstructing amplitude and phase of a wave field [START_REF] Gabor | A new microscopic principle[END_REF]. Since then, many techniques have been developed to reconstruct and record holograms with the same aim of providing simpler solutions and higher quality [START_REF] Yaras | State of the art in holographic displays: A survey[END_REF]. The physical concept of light propagation can be modelled and simulated on a computer, enabling the production of computer generated holograms (CGHs) [START_REF] Tricoles | Computer generated holograms: an historical review[END_REF], [START_REF] Symeonidou | Computergenerated holograms by multiple wavefront recording plane method with occlusion culling[END_REF] which have the advantage to avoid the need for a physical holographic recording set-up.

The visual richness of digital holograms is naturally associated to large amounts of data, which means huge storage and transmission bandwidth requirements. Consequently, the issue of efficiently compressing this type of data is a major concern. Hologram compression differs to image compression principally because holograms encode 3D information in complex-valued pixels and secondly because of their speckle nature. Due to the fact that holographic data represents the diffraction patterns that correspond to 3D scene information as intensity and/or phase over a bidirectional plane, the majority of research work on holographic data compression have focused on the extension of conventional image/video coding standards and algorithms. Fundamentally there are two modes of image compression: lossless and lossy compression. Lossless compression techniques have been already performed on digital holographic data [START_REF] Naughton | Compression of digital holograms for three-dimensional object reconstruction and recognition[END_REF]. Even if perfect reconstruction is to be attained, reversible compression methods are usually inefficient for holographic data [START_REF] Darakis | Measurement of compression defects in phase-shifting digital holographic data[END_REF]. This is principally due to the low spatial redundancies of hologram's speckle nature. Therefore, lossy compression of holographic data still play a critical role in reducing the storage and transmission costs [START_REF] Xing | Comparative study of scalar and vector quantization on different phase-shifting digital holographic data representations[END_REF].

To date, there is no coding standard appropriate for holographic images. Various reviews, summarizing the pioneering works on holographic data compression were published [START_REF] Xing | Digital holographic data representation and compression[END_REF], [START_REF] Dufaux | Compression of digital holographic data: an overview[END_REF], [START_REF] Blinder | Open access database for experimental validations of holographic compression engines[END_REF], [START_REF] Arrifano | Multiple description coding of digital holograms using Maximum-a-Posteriori[END_REF], [START_REF] Peixeiro | Digital holography: Benchmarking coding standards and representation formats[END_REF]. On the other hand, the upcoming JPEG Pleno standard is actually working on the issue of high efficiency coding of plenoptic images and digital holographic data while maintaining backward compatibility with legacy JPEG decoders [START_REF] Pleno | Executive Summary[END_REF]. It is evident that a common representation framework, which is backward compatible with JPEG would provide interoperability that might enable a wider use of this new type of image modalities without breaking the current JPEG ecosystem.

The aim of this paper is to study the effectiveness of JPEG based solution for DH compression. The main idea consists in using a rate-control compression algorithm, which uses an optimized bit allocation mechanism while guaranteeing a full compliance with the JPEG bitstream. The performance evaluation study includes the reconstruction quality and the computational complexity.

The paper is structured as follows. Section II focuses on the description of the proposed compression algorithm, namely JPEG-OPT [START_REF] Ratnakar | RD-OPT : An efficient algorithm for optimizing DCT quantization tables[END_REF]. Section III assesses the performance of the aforementioned coding algorithm when using holographic data. Finally, Section IV presents the conclusions.

II. DESCRIPTION OF JPEG-OPT ALGORITHM

As cited earlier, the objective of this work is to study the effectiveness of JPEG-based coding algorithm, named JPEG-OPT, for holographic data while considering the compression performance and computation cost. JPEG-OPT has the advantage of ensuring near-optimal rate control in DCT-based compression while guaranteeing backward compatibility with JPEG legacy format. The subsequent algorithmic parts of the JPEG-OPT compression system are illustrated in a block diagram form as shown in Fig. 1.

Let us designate an H × W scaled holographic data by a set I of H×W 64 blocks. For simplicity purpose, the 8 × 8 set of spatial frequencies, in the DCT-domain is ordered into a 1-D array of 64 coefficients using the zig-zag scan. A DCT coefficient with coordinates (i, j), (i, j = 1, ..., 8), is then indexed by n = 8 × i + j. The n th component belonging to the DCT block C b is then denoted by C[n] b , where n is the coefficient position and b is the block index.

Before passing through the Huffman entropy coder, the DCT indices undergo run-length coding, which combines zero DCT indices from different frequency bins into one symbol.

The JPEG syntax leaves the selection of the quantization step sizes to the encoder. Unlike baseline JPEG encoder, JPEG-OPT uses a rate-distortion optimization scheme to obtain an efficient image-adaptive quantization table. The ratedistortion optimization problem is to find a quantization table

Q, such that the rate R(Q) ≤ R * , and D(Q) is minimized, for a given rate budget R * .
Because of the zero run-length coding, it is difficult to obtain a solution to the rate-distortion optimisation problem with classical bit allocation methods. Ratnakar and Miron Livny [START_REF] Ratnakar | RD-OPT : An efficient algorithm for optimizing DCT quantization tables[END_REF] proposed a sub-optimal solution to solve the quantization table design problem using histograms of the probability distribution of the DCT coefficients over all the blocks.

Based on Parseval theorem, the total distortion that affects an input block I b can be replaced by the sum of the coefficient distortions in its corresponding DCT coefficient block C b . On the other hand, the bit rate can be estimated as the sum (over the 64 frequency localizations) of the entropies of the quantized DCT coefficients over all the blocks. Given a quantization table Q, let D(Q) and R(Q) represent the Mean square error (MSE) and estimated bit rate (in bits per sample), respectively, resulting from compression of a scaled holographic data using the quantization table Q. For a given quantizer step size q, D[n][q], the contribution to total distortion at the n th frequency position is given by:

D[n][q] = 1 64 M ean((C[n] -Round( C[n] q )) 2 ), (1)
where the "Mean" is taken over all the 8 × 8 blocks in the quantized and transformed domain. Similarly, the rate contribution, R[n][q n ], at the n th frequency position is estimated by:

R[n][q] = Entropy(Round( C[n] q )), (2)

R[n]

[q] is calculated using the histogram of DCT coefficient distribution in the n th frequency position.

For a specific quantization table Q, the modeled bit-rate is given by:

R(Q) = 63 ∑ n=0 R[n][Q[n]], (3)
whereas the distortion is estimated by:

D(Q) = 63 ∑ n=0 D[n][Q[n]]. (4)
Giving Q m and Q M the quantization tables that, respectively, represent the minimum and maximum values possible for each quantization table entry, such that,

Q m [n] ≺ Q[n] ≤ Q M [n].
The rate-distortion optimisation scheme builds tables of D[n][q] and R[n][q], for each quantizer step size q = Q[n] in the predefined operating range.

The rate-distortion optimisation problem is then formulated by Lagrangian minimization, min(D(Q) + λR(Q)), in such manner that solutions of the latter for any non-negative λ are solutions to the former for a target bit rate R * . Using both Eqs 3 and 4, the Lagrangian minimisation problem can be expressed by:

63 ∑ n=0 D[n][Q[n]] + λR[n][Q[n]]. ( 5
)
For each frequency localization n, a subset of the operating points The slopes of the rate-distortion curve for the n th coefficients at these h n points are given by:

(R[n][.], D[n][.])
λ n (k) = R n (q n (k)) -R n (q n (k + 1)) D n (q n (k + 1)) -D n (q n (k)) (6)
for k = 1, 2, ..., h n -1 and λ n (h n ) = 0. Hence, λ n (k) > λ n (k), for k = 1, 2, ..., h n -1. Notably, for any given λ > 0, the Lagrangian, D[n][q] + λR[n][q] is minimized for q = q n (k), with k being the least index for which λ > λ n (k). More generally, for any given λ > 0, the 64 binary searches can be used to find the quantizer values, belonging to the operating range, such that

R(Q) + λD(Q) = ∑ 63 n=0 D[n][Q[n]] + λR[n][Q[n]] is minimized.
If a target bit rate R * is considered, the bisection method is used to easily search for a λ that would satisfy the target with respect to the fixed tolerance.

To have a better understanding of the effectiveness of JPEG-OPT compression method in the context of digital holograms, experiments have been conducted to evaluate its efficiency.

III. PERFORMANCE EVALUATION

Experiments have been conducted to evaluate JPEG-OPT coding scheme in the context of lossy compression of digital holograms. To complete our study, comparisons have been made with reference methods, on the basis of bit rate, reconstruction error and computational cost.

A. Experimental setup

Due to space limitation, our simulations have been carried out using three CGHs (Girl, 3D Venus and Horse) selected from three databases: Télécom ParisTech, Interfere-I and EmergIMG. The first holographic content dataset was developed in the framework of the PhD research work of Y. Xing at Télécom ParisTech (see [START_REF] Xing | Digital holographic data representation and compression[END_REF]). The Interfere-I dataset was made available by Prof. Peter Schelkens from ETRO-VUB and is available at http://www.erc-interfere.eu./index.html. The EmergIMG database was made available by Prof. Manuela Pereira from University of Beira Interior, Covilh and is available at http://emergimg.di.ubi.pt/downloads.html. The basic setups, for the considered holograms are summarized in Table I. As mentioned in [START_REF] Dufaux | Compression of digital holographic data: an overview[END_REF], the interferograms can be represented using various formatting options. Among these we use shifted distance data (D1/D2) and Real/Imaginary (Rea/Ima) representations. Note that because of the chaotic nature of Amplitude/Phase format, the use of such data is unbeneficial when targeting compression. Also, the three interferograms intensities format requires more components to represent a hologram and thus it is evidently the least efficient format for the compression purpose.

For the assessment purpose, objective comparisons between JPEG-OPT, JPEG [START_REF] Wallace | The JPEG still picture compression standard[END_REF], QT-L [START_REF] Munteanu | Wavelet image compression -the quadtree coding approachs[END_REF], JPEG 2000 (using Kakadu implementation) [START_REF] Skodras | The JPEG 2000 still image compression standard[END_REF], SPIHT [START_REF] Said | A new fast and efficient image codec based on set partitioning in hierarchical trees[END_REF] and HEVC [START_REF] Sullivan | Overview of the High Efficiency Video Coding (HEVC) Standard[END_REF] have been conducted to evaluate the quality of the decoded holograms in the reconstruction domain. JPEG 2000, QT-L and SPIHT are wavelet-based compression techniques. The first two codecs exploit the intra-band dependencies between wavelet coefficients, whereas the third exploits the inter-band dependencies between the wavelet sub-bands and employs a zero-tree structure to drop the insignificant coefficients. To offer a fair comparison, the tested wavelet-based coders used the following set of parameters: irreversible (9,7) filter-bank; 3 levels of dyadic wavelet decomposition; and a single layer not considering the progressive decoding mode. It is noteworthy that YUV format is the main raw video format used in video coding softwares. Particularly, the Main Still Picture profile in HEVC software only accepts YUV 4:2:0 format as input. So, before HEVC encoding, the digital holograms were converted to YUV 4:2:0 format to satisfy the unique input format supported by the reference software.

The Angular Spectrum Method (ASM) has been used to perform the reconstruction process for Girl and 3D Venus, whereas Fresnel diffraction Method (FDM) has been used to perform the reconstruction for Horse. It is worth noting that the ASM did not work on Horse because the object is substantially bigger than the hologram. Since angular spectrum preserves pixel pitch, it will not be able to fully show the object. On the other hand, by using numerical Fresnel diffraction (the Fourier form), the pixel pitch will change depending on the propagation distance, and it will grow so we can fully view the entire hologram. The softwares of both ASM and FDM methods have been kindly provided by Prof. Manuela Pereira from University of Beira Interior, Covilha in Portugal. In the proposed assessment methodology, the reconstructions associated to the original and decoded holographic data are compared by means of BD-PSNR using Bj∅ntegaard model [START_REF] Bj∅ntegaard | Calculation of average PSNR differences between RDcurves[END_REF]. In this case of study, the bit rate is the sum of the bit rates for the two compressed holographic data components. The BD-PSNR metric receives two sets of Rate-Distortion (RD) points defining two RD performance curves and measures the average difference in dB between the PSNR values for a same rate.

In order to evaluate the time complexity of the tested codecs, execution-time tests were conducted employing 32-bit C++ executables generated by Visual Studio 10.0 operating in release mode. The time-complexity evaluation for JPEG 2000, QT-L and SPIHT was carried out using the executables that are made publicly available by the authors.

B. Performance assessment in the reconstruction Domain

In this section, we evaluate the compression performance of JPEG-OPT, JPEG 2000, QT-L, SPIHT and JPEG in the reconstruction domain. The results of our experiments in terms of Bj∅ntegaard Delta-PSNR, when taking HEVC as a reference codec, are summarized in Table II.

Taking the Delta-PSNR to serve as a metric to compare the coding efficiency in the reconstruction domain, the gap between the average Delta-PSNR values using Real/Imaginary and shifted distance data representations for 3D Venus and Girl varies is 0.146 dB and 0.47 dB, respectively, which is relatively low. One may conclude that the performances of the Real/Imaginary components have results very close to those obtained for the shifted distances components. The results also reveal that HEVC significantly outperforms all the tested codecs, except for Real/Imaginary representation of Horse object, where HEVC yields lower performance than the other tested codecs. One may thus conclude that the speckle noise in Horse object has a greater impact on HEVC efficiency.

Examining the efficiency of JPEG-OPT, JPEG 2000, QT-L, SPIHT and JPEG, we observe that, for Girl object, JPEG 2000 performs best, followed by SPIHT. On the other hand, QT-L yields superior performance for 3D Venus whereas JPEG-OPT systematically yields superior reconstruction quality for Horse object. One may also notice that for the same object, both QT-L and SPIHT offer roughly similar performance despite the fact that QT-L is an intra-band coding method whereas SPIHT is an inter-band coding method. More generally, the reported average Delta-PSNR results show that JPEG-OPT produces satisfactory results which are competitive with those provided by JPEG 2000.

In Fig. 2, visual comparisons between JPEG-OPT and HEVC are provided using the reconstructed versions of 3D Venus and Girl object. Examining these renderings, we can clearly observe that for Real/Imaginary representation, JPEG-OPT architecture provides a visually appealing quality even for low bit rate. Besides, one may notice that both JPEG-OPT and HEVC methods are able to produce very similar results despite the difference in PSNR in favour of HEVC.

C. Computation times

Table III provides the coding and decoding execution time, in second (s), as a function of the tested representation formats, for 3D Venus, Girl and Horse holograms. Before discussing the results of the complexity assessments, we would like to emphasize the fact that the timings reported in Table III are very implementation dependant. But, they allow to know whether the order of complexity of the tested technique is acceptable or not. To prevent possible external processes from affecting the result, For each test condition, the runs were repeated 10 times and the average running time of repeated execution is retained. From the reported results, one may notice that, for the wavelet-based and HEVC codecs, the full processing time is highly influenced by the resolution of the test hologram. When comparing the coding delays of the four codecs, we can plainly discern that HEVC is very time-consuming compared to the other three tested codecs. This is not surprising since the quad-tree structured coding units and the Inter prediction with high number of modes, that characterise HEVC coding technique, cause high computational complexity. Note that in our simulation we used the HEVC Test Model 16.9 of HEVC. On the other hand, JPEG-OPT codec produces the fastest encoding and decoding times. On average, the JPEG-OPT encoder is 35 times faster than JPEG 2000, 83 times faster than SPIHT and 232 times faster than QT-L, which is more than fast enough for real-time applications.

D. Overview of the analysis of simulation experiments

To summarize, the results of the conducted experiments we can plainly discern that:

• In most cases, HEVC produces the best quality. But this is at the expense of extremely high coding and decoding time.

• JPEG 2000 requires the lowest encoding and decoding times among the wavelet-based codecs. • QT-L is the slowest among the wavelet-based codecs and often produces Delta-PSNR results comparable to those obtained with SPIHT. • JPEG-OPT yields very competitive compression performance and reconstruction quality that are sometimes even better than those of wavelet-based techniques. As a final note, it is obvious that fast encoding/decoding is a crucial task in modern low-complexity applications such as embedded systems. In addition, JPEG-OPT coding scheme inherently supports backward and forward compatibility with legacy JPEG formats, which is not the case for wavelet-based codecs.

IV. CONCLUSIONS

In this paper, we explore the efficiency of a DCT-based coding scheme for the compression of digital holograms. The proposed codec uses a Lagrangian-driven rate allocation mechanism while still compliant with the JPEG bitstream. Regarding compression efficiency and reconstruction quality, it yields competitive results with wavelet-based codecs that also include rate/distortion optimisation mechanism. Moreover, an empirical complexity evaluation showed that it is 37 times faster that JPEG 2000. Finally, we believe that contrarily to DH specific coding solutions in the literature, the optimized JPEG coding scheme offers a good balance of speed, file size and reconstruction quality. Since it guarantees backward compatibility with legacy JPEG decoders, it may serve as an attractive alternative to the current state-of-the-art and that the conclusions derived from the conducted studies are helpful for developing new JPEG-based hologram codec in line with the objectives of the upcoming JPEG Pleno standard.

The research done in this work could be extended by using perceptual metrics. Since the final destination of a compressed hologram is often to be displayed to a human observer, the perceived visual quality should be the criterion to take into account when evaluating a compression algorithm.

Fig. 1 .

 1 Fig. 1. Block diagram for the JPEG-OPT coding chain.

Fig. 2 .

 2 Fig. 2. Visual comparison of the reconstructed 3D Venus and Girl objects (real/imaginary format) after being compressed at 0.4 bpp, using (c-d) JPEG-OPT and (e-f) HEVC. Reconstructions of the original objects are shown in (a-b)

TABLE I PARAMETERS

 I OF CGH GENERATION SETUPS.

			Specifications of the different
	Database Hologram	generation setups
			Resolution (pixels)	Pixel pitch (µm)	Recons-truction distance (cm)	Wave length (nm)
	Télécom ParisTech	Girl	512 × 512	6.33	6.48	633
	Interfere-I	3D Venus 1920 × 1080	8	0.51	632.8
	EmergImg HoloGrail	Horse	968 × 968	4.4	14	632.8

TABLE II COMPRESSION

 II PERFORMANCE EVALUATION (IN BD-PSNR [DB] REGARDING HEVC INTRA) WHEN USING THE RECONSTRUCTED HOLOGRAMS.

TABLE III EXECUTION

 III TIMES (IN SECONDS) FOR CODING THE HOLOGRAMS IN THE TEST SET.

	Model		JPEG-OPT		QT-L	SPIHT	JPEG2000	HEVC
		Coding	Decoding	Coding Decoding Coding Decoding Coding Decoding	Coding	Decoding
	Girl	0.0018	0.0395×10 -3	0.1245	0.0907	0.0390	0.0225	0.0155	0.015	40.6500	0.3970
	3DVenus	0.0019	0.1390×10 -3	0.7578	0.7278	0.2783	0.2016	0.1205	0.0385	404.2500	3.5430
	Horse	0.0017	0.1287×10 -3	0.2867	0.2756	0.33	0.0892	0.0672	0.0542	125.664	1.2075
	Average	0.0018	0.1024 ×10 -3	0.3897	0.3647	0.2158	0.1044	0.0677	0.0359	190.188	1.7158

ACKNOWLEDGMENT

The authors are particularly grateful to David Blinder and Peter Schelkensa, (ETRO-VUB, Brussels, Belgium), Manuella Pereira (University of Beira Interior, Covilha, Portugal), José Peixeiro (Técnico Lisboa, Portugal) and Patrick Gioia (Orange Labs, Rennes, France) for their support, the very interesting discussions and for giving us access to their holographic display.