Andrea Federico

Sabattoli

Houssein El Sabattoli

Laurène Dirani

Francesco Youssef

Davide Garrisi

Luca Grassani

Camille Zatti

Erwine Petit-Etienne

J E Pargon

Sipe

Houssein El Dirani

Laurène Youssef

Francesco Garrisi

Davide Grassani

Luca Zatti

Camille Petit-Etienne

Erwine Pargon

J E Sipe

Marco Liscidini

Corrado Sciancalepore

Daniele Bajoni

Matteo Galli

Suppression of Parasitic Nonlinear Processes in Spontaneous Four-Wave Mixing with Linearly Uncoupled Resonators

come   L'archive ouverte pluridisciplinaire

Nonlinear optics has been widely exploited to produce nonclassical states of the electromagnetic field. There has been a recent focus on the use of CMOS compatible materials, such as silicon and silicon nitride, which are characterized by a strong third-order optical response that can be further enhanced by the light confinement attainable in integrated microresonators [1]. With the fabrication quality and technological maturity achieved in the last decade, integrated devices based on these materials offer the promise of efficient generation of entangled photons [2][3][4][5][6], heralded single-photons [7,8] and squeezed light [9,10].

A key third-order nonlinear process is dual-pump spontaneous four-wave mixing (d-SFWM), in which pump lasers at ω 1 and ω 2 are used to generate photons at ω T = (ω 1 + ω 2 )/2. This process can be used to generate single-mode squeezed light, a central resource in the development of continuous-variable quantum computing [11]. While conventional micrometric ring resonators have been used to enhance d-SFWM [START_REF] Guo | Optics InfoBase Conference Papers[END_REF], their multimode nature inevitably also leads to the amplification of unwanted nonlinear phenomena, such as singlepump spontaneous four-wave mixing (s-SFWM). This is a problem: as shown in Fig. 1(a), if the two pumps are resonant with the microring, the generation of photon pairs within a target resonance (ω T ) through d-SFWM occurs along with two s-SFWM processes. Both of these parasitic s-SFWM processes lead to the generation of two photons, with one, and only one, of them in the mode at ω T . This acts as a source of noise, leading to a degradation of the squeezed light [10]).

To quantitatively characterize these nonlinear parasitic processes, in Fig. 1(b) we plot the photon generation rates (measured as described in the Supplemental Material) due to the d-SFWM and s-SFWM processes in an integrated silicon microring resonator (length 372 µm and quality factor around 6 • 10 4 ). As expected, the generation rate for the d-SFWM is unchanged as long as the product of the two pump powers is held fixed, while the generation rate of the s-SFWM processes scale quadratically with the appropriate pump powers [START_REF] Helt | [END_REF]. The generation rates of the parasitic processes, given by the sum of the two s-SFWM rates, can even be one order of magnitude larger than the generation rate of the d-SFWM process, with a maximum signal-to-noise ratio (SNR) of only about 2 when the two pump powers are equal.

Strategies have been proposed to overcome this intrinsic limit. In one approach, the two pump lasers are slightly detuned from resonance, such that the s-SFWM rates are suppressed [10], but with a trade-off between SNR and generation rate. Alternatively, one can use systems composed of two or more linearly coupled resonators to engineer the spectral position of the resonant modes, with nonlinear phenomena selectively enhanced or suppressed [9,14,15]. While this can be very effective, the suppression of the parasitic processes is limited by the achievable resonant splitting of the coupled modes.

In 2019, Menotti et al. proposed a new strategy [16], in which two optical resonators are coupled solely through the nonlinear interaction. The system is characterized by two independent sets of resonances that can be tuned to selectively enhance or suppress different nonlinear phenomena. In Fig. 2(b) we represent a sketch of the realisation of d-SFWM in such a system. Two laser pumps are tuned to resonances of the first resonator, leading to the generation of photon-pairs within a target resonance of the second resonator. Yet s-SFWMs are expected to be suppressed, for the overall field enhancement All the processes are resonant and exploit a significant intensity (I.E.). (b) Measurement of the pair generation rate in a single silicon microring resonator when the product of the powers of the two pumps is fixed. The blue triangles and the red circles show the s-SFWM processes pumped by the lasers at ω1 and ω2, respectively. The black squares show the sum of the parasitic process rates, and the green stars show the d-SFWM generation rate. The measures are compared with the theoretical calculation (shown with the dotted blue line for the P1-pumped s-SFWM, the dashed red line for the P2-pumped s-SFWM, and the green line for d-SFWM.) See Supplemental Material for the experimental details, and Menotti et al. [16] for the theoretical calculations.

at ω S = 2ω 2 -ω T and ω I = 2ω 1 -ω T is strongly reduced by the absence of three resonances satisfying energy conservation. In this letter we experimentally demonstrate this scheme.

The device is presented in Fig. 2(a): two silicon-oninsulator (SOI) racetrack resonators (waveguide crosssection 600x220 nm 2 ) with different lengths L 1 = 372 µm and L 2 = 366 µm, are located side-by-side forming a directional coupler (DC) of length L DC = 92 µm, chosen to minimise cross-transmission and achieve isolation of the two resonators. Thus one set of modes is associated with Resonator 1 (R 1 ), and another with Resonator 2 (R 2 ). During the fabrication process, a H 2 thermal annealing was used to smooth the sidewalls of silicon waveguides [17]: this leads to the reduction of the propagation losses and, thus, to the realisation of resonators with high quality factors [START_REF] Sabattoli | 21st International Conference on Transparent Optical Networks (ICTON)[END_REF]. The optical characterization is shown in Fig. 2(c). The resulting two independent sets of resonances are clearly observed by injecting light from a tunable laser in the ports T 1 or T 2 and detecting the optical response at the ports T 6 and T 3 , respectively. The free spectral range (FSR) at 1550 nm is around 1.566 nm and 1.595 nm for R 1 and R 2 respectively. The loaded quality factors are around 6 • 10 4 and 3 • 10 5 , with the difference due to the different number of waveguides coupled to the racetracks (two to R 1 , and one to R 2 ). The optical isolation is at least 20 dB over a 40-nm bandwidth, as observed from the transmission spectrum T 1 → T 3 . Resonator 1

Resonator 2

ω 2 ω S ω T ω 1 ω I T 1 →T 6 T 2 →T 3 T 1 →T 3 T 5 T 4 T 3 T 1 T 2 FIG.
D 1 D 2 D 3 D 4 BS 1 PM 1 PM 2 + - R 1 R 2 C 2 C 3 C 4 T 1 T 3 FIG. 3.
Experimental set-up used to measure the coincidence rates. Pump light from two lasers, combined into the same fiber through a beam splitter BS1, is coupled to resonator R1, the resonances of which are tuned by changing the voltage applied to the heater. Photons generated through d-SFMW and s-SFWM are collected with a tapered fiber, and sorted by frequency using tunable fiber Bragg filters (TFB1-TFB3) and circulators (C1-C3). They are then detected using four superconductive single photon detectors (D1-D4). The correlations between the detection times are recorded and analysed through time tagging electronics. The spectral alignment of the lasers with the device is checked with two power meters (PM1 and PM2), with an additional fiber Bragg filter (TFB4) and circulator (C4) used to separate the two beams. Bandpass filters BPF1 and BPF2 are used to eliminate amplified spontaneous emission from the output of the pump lasers, and BPF3-BPF5 to filter out residual pump photons.

The two sets of resonances can be controlled independently by means of two electric heaters fabricated with a TiN layer on top of the silica cladding, and the electrical contact is realised with AlCu electric wires deposited on top of the chip.

Measurements of the nonlinear processes are realised with the experimental set-up presented in Fig. 3. The light from two tunable lasers (Santec TSL710 and TSL510) is coupled using a collimator and lens to match the mode field diameter of the edge coupling tip on the chip, after two band-pass filters BPF 1 and BPF 2 (FWHM around 8 nm) have been used to filter out the amplified spontaneous emission from the lasers. The spectral position of the resonances of R 1 is controlled by using an external voltage source (Keithley 2400 Source Meter) and a multi-contact wedge. The generated photon pairs are collected by a tapered fiber coupled to the port T 3 of the device, and the fiber is connected to the circulator C 1 . A tunable fiber Bragg filter (TFB 1 ) is used in reflection [65 pm 3dB-bandwidth (BW) and 40 dB of extinction ratio (ER)] to route photons at frequencies around ω T towards the detection system. Another band-pass filter (BPF 3 ) with more than 60 dB of ER filters out any residuals from the pump, and a fiber 50:50 beam splitter splits the photon flux towards the superconductive single-photon detectors D 1 and D 2 supplied by Photon Spot.

Photons at ω T can be generated through s-SFWM only if another photon is generated to guarantee the conservation of energy. In particular, s-SFWM pumped by the laser at ω 1 (P 1 -SFMW) can generate a photon pair with one photon at ω T and one at ω I , such that ω I = 2ω 1 -ω T . Similarly, s-SFWM associated with the laser at ω 2 (P 2 -SFMW) can generate one photon at ω T and one photon at ω S , where ω S = 2ω 2 -ω T . To detect these photons, we use two tunable fiber Bragg filters (TFB 2 and TFB 3 ) and two band-pass filters (BPF 4 and BPF 5 ) to route photons at ω I and ω S towards superconductive single-photon detectors D 3 and D 4 , respectively. Finally, the alignment of the two laser pumps is checked during the whole experiment using two power meters (PM 1 and PM 2 ), a circulator (C 4 ), and an additional band-pass filter (TFB 4 ).

The total insertion losses from the lasers to the waveguide coupled with R 1 at the port T 1 are around -7.8±0.1 dB for the light from both P 1 and P 2 . The losses from the waveguide coupled with R 1 at the port T 3 to D 1 and D 2 are around -10.7 ± 0.3 dB and -10.1 ± 0.3 dB, respectively, for the photons at ω T ; the losses from T 3 to D 3 are -11.9 ± 0.3 dB for photons at ω I , and -7.6 ± 0.3 dB from T 3 to D 4 for photons at ω S . The efficiency of the detectors is -1.00 ± 0.05 dB.

We begin by verifying the capability of our structure to suppress s-SFWM. To do that, we measure the vari-ation of the efficiency of an s-SFWM process by tuning the two sets of resonances, which leads to a modification of the nonlinear coupling of the modes. We start by adjusting the voltage applied to the heater of R 1 so that the frequency ω 2 of a resonance (vacuum wavelength λ 2 = 1560.916 nm) is in the middle of the frequencies of two resonances of R 2 (at vacuum wavelengths λ T = 1550.495 nm and λ U = 1571.479 nm). This configuration should guarantee a maximum efficiency of P 2 -SFWM, since ω U = ω S = 2ω 2 -ω T . Then, with laser P 1 turned off, we change the driving voltage of the heater of R 1 to shift its resonances, and we tune P 2 at each step to be on resonance with R 1 at each λ 2 by checking the transmission on the power meter PM 2 . The tunable fiber Bragg filter TFB 1 is aligned to the resonance at ω T at the beginning of the experiment, since thermal cross-talk leads to a tuning of the resonances of R 2 within the bandwidth of the filter. On the other hand, the filter TFB 2 is tuned to collect the photons generated at ω S = 2ω 2 -ω T , which requires the adjustment of its spectral alignment during the whole experiment. In Fig. 4 we show the attenuation of the single-pump SFWM process with respect to the detuning of the resonances involved, which is estimated from the driving voltage and the characterization of the heater (see Supplemental Material for details). In particular, that detuning is defined as ∆ω 2 -(∆ω T -∆ω U )/2. We find that the simulation and the data are in good agreement, and exhibit the expected behaviour. This result shows that by tuning the two sets of resonances one can manipulate the nonlinear coupling of the modes of the resonators, and, consequently, the generation of pairs due to the nonlinear processes occurring in the device. In particular, a detuning of around 65 pm (8.1 GHz, around 13 times the linewidth of the resonances of R 2 ) leads to a suppression of around three orders of magnitude with respect to perfect nonlinear coupling.

With the possibility of attenuating the parasitic processes demonstrated, we proceed to the dual-pump SFWM experiment. The laser pumps are aligned to two resonances of R 1 (vacuum wavelengths λ 1 = 1540.374 nm and λ 2 = 1560.742 nm), which are equally spaced in frequency about a resonance of R 2 at λ T = 1550.491 nm. This configuration is associated with a detuning of the s-SFWM resonances of around -160 pm (-20 GHz), and the estimation of the resulting attenuation of the s-SFWM processes is around -37.3 dB. The quality factors of the resonances of R 1 are around Q 1 = 6.7 • 10 4 and Q 2 = 5.7 • 10 4 , while the resonance of R 2 at λ T has a much higher quality factor, Q T = 3.2 • 10 5 . We align the tunable filter TFB 1 to λ T , while TFB 2 and TFB 3 are tuned to route photons at λ I = 1530.388 nm and λ S = 1571.129 nm to D 3 and D 4 , respectively. The analysis of the coincidences on the detectors allows for the recognition of the pairs generated through the different processes. In particular, the d-SFMW should give coincidence between the detection events at D 1 and D 2 . Similarly, the single-pump processes should result in coincident detection: P 1 -SFWM between the events at D 3 and D 1 (or D 3 and D 2 ), and P 2 -SFWM between the events at D 4 and D 1 (or D 4 and D 2 ).

As discussed earlier [16], the efficiency of the d-SFWM process in our device is reduced from what could be achieved with a simple microring, since here the nonlinear interaction occurs only in the directional coupler region shared by the two rings. This decreases both the overlap integral of the fields (1/4) and the interaction length (L DC ≈ L/4). This would be a limitation in the framework of squeezed-light sources, where high generation efficiency is required, but it could be overcome by working with ultra low-loss waveguides with high modal confinement in silicon nitride [START_REF] Dirani | [END_REF]20]. In such a device one could also increase the pump power, since the twophoton absorption that affects silicon waveguides would not be present [21].

In Fig. 5 we show the result of a coincidence experiment with the current device, where 0.95 ± 0.05 mW of optical power is coupled into the waveguide from each pump laser. The data are acquired for 20 minutes, and the histograms have a bin-width equal to 35 ps, which is comparable with the time jitter of the detectors. The red peak shows that photons arriving on D 1 and D 2 are emitted at the same time, since their arrival time is correlated: this clearly demonstrates they are emitted through d-SFWM, since any other process that can generate photons at λ T cannot be characterized by this temporal correlation. We measure a coincidence rate equal to 164.2±0.4 Hz, which corresponds to an internal generation rate in the second resonator equal to 62 ± 6 kHz, and the CAR is 1190 ± 10. The coincidence rate is estimated by integrating the peak within the whole window (and subtracting the noise), while the CAR is measured by taking into account the FWHM of the peak (16 bins of the histogram).

From the histograms in Fig. 5 it is also possible to evaluate the photon pairs generated through s-SFWM, which lead to coincidence events between D 1 and D 2 with D 3 or D 4 . As can be seen, no peaks are clearly visible. In fact, the dark counts, the noise from the environment, and other parasitic processes occurring in the set-up (such as Raman emission from optical fibers) hide the signal given by time-correlated photons. This is due to the severe attenuation of the processes. Nevertheless, we can estimate a lower bound for the SNR by analysing the coincidences between D 1 and D 4 . We integrate the peak due to d-SFWM over its FWHM (16 bins), and we divide it by the noise, which is calculated as the average on the black histogram in Fig. 5 multiplied by 16 bins. By taking into account the losses, we obtain a SNR equal to (11.3 ± 3.2) • 10 3 . From the model, we should expect that the SNR is around 21 • 10 3 , an improvement over that of a single ring system by over four orders of magnitude. Our experimental estimation is lower, since we are overestimating the noise due to s-SFWM.

These results demonstrate that the use of linearly un- coupled resonators allows for the effective control of the nonlinear interaction between optical modes, leading to the ability to enhance or suppress the generation of photon pairs through spontaneous four wave mixing. This kind of device will be a central component in strategies for the implementation of continuous-variable quantum computing. It will enable the possibility of realising sources of highly squeezed light, since the main sources of noise, due to single-pump spontaneous four wave mixing, can be effectively suppressed. In our experiment, we could confirm the improvement of the signal-to-noise ratio by over four orders of magnitude. This work has been supported by Ministero dell'Istruzione, dell'Università e della Ricerca (Dipartimenti di Eccellenza Program (2018-2022)). J.E.S. acknowledges support from the Natural Sciences and Engineering Research Council of Canada. The device has been designed using the open source Nazca design TM framework.

FIG. 1 .

 1 FIG. 1. (a) Sketch of the s-SFWM and d-SFWM processes occurring in a single microring resonator when two laser pumps P1 and P2 are aligned to two resonances (ω1 and ω2). All the processes are resonant and exploit a significant intensity (I.E.). (b) Measurement of the pair generation rate in a single silicon microring resonator when the product of the powers of the two pumps is fixed. The blue triangles and the red circles show the s-SFWM processes pumped by the lasers at ω1 and ω2, respectively. The black squares show the sum of the parasitic process rates, and the green stars show the d-SFWM generation rate. The measures are compared with the theoretical calculation (shown with the dotted blue line for the P1-pumped s-SFWM, the dashed red line for the P2-pumped s-SFWM, and the green line for d-SFWM.) See Supplemental Material for the experimental details, and Menotti et al. [16] for the theoretical calculations.

  FIG. 2. (a) Optical image of the device with two linearly uncoupled racetracks. A second image of the device without the electrical circuits has been superimposed to show the position of waveguides, which are covered by the metal. (b) Sketch of the s-SFWM and d-SFWM processes occurring in the device when two laser pumps are aligned to two resonances (ω1 and ω2). Single-pump processes which generate photons at ωT are suppressed due to the weak field enhancement experienced at the spectral position where energy conservation is preserved. (c) Linear characterization of the device.

FIG. 4 .

 4 FIG. 4. Plot showing the attenuation of s-SFWM as a function of the detuning of the resonances of the two resonators.

3 FIG. 5 .

 35 FIG.5. Histograms representing the coincidences events detected. The red, black, and blue histograms represent, respectively, the coincidences due to d-SFWM (events on D1 and D2), P2-SFWM (events on D1 and D4), and P1-SFWM (events on D1 and D3).The width of each bin is 35 ps. The delay is calculated by taking into account the different optical lengths of the fibers from the device to the detectors.


Bruno Lecouat

Jean Ponce

Julien Mairal

Lucas-Kanade Reloaded: End-to-End Super-Resolution from Raw Image Bursts

Output of Dcraw

Camera's ISP (hi-quality jpeg) Joint demosaick+single-image SR Proposed method

: ×4 super-resolution results obtained from a burst of 30 raw images acquired with a handheld Panasonic Lumix GX9 camera at 12800 ISO for the top image and 25600 for the bottom image. Dcraw performs basic demosaicking.

Introduction

The problem of reconstructing high-resolution (HR) images from lower-resolution (LR) ones comes in multiple flavors, that may significantly differ from each other in both technical detail and overall objectives. When a single LR image is available, the corresponding inverse problem is severely ill-posed, requiring very strong priors about the type of picture under consideration [18,[START_REF] Yang | Image super-resolution: historical overview and future challenges[END_REF]. For natural images, data-driven methods based on convolutional neural networks (CNNs) have proven to be very effective [START_REF] Ledig | Photorealistic single image super-resolution using a generative adversarial network[END_REF][START_REF] Wang | ESRGAN: Enhanced super-resolution generative adversarial networks[END_REF]. Generative adversarial networks (GANs) have also been used to synthesize impressive HR images that may, how-ever, contain "hallucinated" high-frequency details [9,28].

In the true super-resolution setting [START_REF] Milanfar | Super-resolution imaging[END_REF][START_REF] Dong | Image super-resolution using deep convolutional networks[END_REF][START_REF] Yang | Image super-resolution: historical overview and future challenges[END_REF], 1 where multiple LR frames are available, HR details are present in the data, but they are spread among multiple misaligned images, with technical challenges such as recovering subpixel registration, but also the promise of recovering veridical information in applications ranging from amateur photography to astronomy, biological and medical imaging, microscopy imaging, and remote sensing.

Videos are of course a rich source of multiple, closelyrelated pictures of the same scene, with several recent approaches to super-resolution in this domain, often combining data-driven priors from CNNs with self-similarities between frames [START_REF] Jo | Deep video super-resolution network using dynamic upsampling filters without explicit motion compensation[END_REF][START_REF] Li | Mucan: Multi-correspondence aggregation network for video super-resolution[END_REF][START_REF] Wang | Edvr: Video restoration with enhanced deformable convolutional networks[END_REF]. However, most digital videos are produced by a complex pipeline mapping raw sensor data to possibly compressed, lower-resolution frames, resulting in a loss of high-frequency details and spatially-correlated noise that may be very difficult to invert [12]. With the ability of modern smartphone and prosumer cameras to record raw image bursts, on the other hand, there is a new opportunity to restore the corresponding frames before the image signal processor (ISP) of the camera produces irremediable damage [4,[START_REF] Wronski | Handheld multi-frame superresolution[END_REF]. This is the problem addressed in this presentation, and it is challenging for several reasons: (i) images typically contain unknown motions due to hand tremor, 2 making subpixel alignment difficult; (ii) converting noisy raw sensor data to full-color images is in itself a difficult problem known as demosaicking [START_REF] Kimmel | Demosaicing: image reconstruction from color ccd samples[END_REF][START_REF] Lecouat | Fully trainable and interpretable non-local sparse models for image restoration[END_REF]; and (iii) effective image priors are often data driven, thus requiring a differentiable estimation procedure for end-to-end learning.

In this paper, we jointly address these issues and propose a new approach that retains the interpretability of classical inverse problem formulations while allowing end-to-end learning of models parameters. This may be seen as a bridge between the "old world" of signal processing and the "brave new one" of data-driven black boxes, without sacrificing interpretability: On the one hand, we address an inverse problem with a model-based optimization procedure alternating motion and HR image estimation steps, directly building on classical work from the 1980s [START_REF] Baker | Lucas-kanade 20 years on: A unifying framework[END_REF][START_REF] Bruce | An iterative image registration technique with an application to stereo vision[END_REF] and 1990s [16]. On the other hand, we also fully exploit modern technology in the form of a plug-and-play prior [START_REF] Stanley H Chan | Plugand-play admm for image restoration: Fixed-point convergence and applications[END_REF][START_REF] Singanallur V Venkatakrishnan | Plug-and-play priors for model based reconstruction[END_REF] that gracefully mixes deep neural networks with variational approaches. In turn, unrolling the optimization procedure [START_REF] Chen | Theoretical linear convergence of unfolded ISTA and its practical weights and thresholds[END_REF][START_REF] Lecouat | Fully trainable and interpretable non-local sparse models for image restoration[END_REF][START_REF] Zhang | Deep unfolding network for image super-resolution[END_REF] allows us to learn the model parameters end to end by using training data with synthetic motions [4].

Since aliasing produces low-frequency artefacts associated with undersampled high-frequency components of the original signal, it is typically considered a nuisance, mo-1 "Single-image super-resolution" has become a popular nickname for single-image upsampling under strong priors; here, we use the classical definition of super-resolution from multiple LR snapshots [START_REF] Dong | Image super-resolution using deep convolutional networks[END_REF][START_REF] Yang | Image super-resolution: historical overview and future challenges[END_REF]. 2 Image bursts acquired on a tripod may also present subpixel misalignments in practice due to floor vibrations, as observed in our experiments.

tivating camera manufacturers to add anti-aliasing (optical) filters in front of the sensor. 3 Yet, aliased images carry highfrequency information, which may be recovered from multiple shifted measurements. Perhaps surprisingly, aliasing is thus an ally in the context of super-resolution, a fact already noted in earlier references, see [41]. As shown in the rest of this presentation, our approach to raw burst superresolution also exploits this insight, and it achieves a new state of the art on several standard benchmarks that use synthetic motion for ground truth. It also gives excellent qualitative results on real data obtained with smartphone and prosumer cameras. Interestingly, as illustrated by Figure 1, our method has turned out to be surprisingly robust to noise given the particularly challenging setting of raw image super-resolution, which involves simultaneous blind denoising, demosaicking, registration, and upsampling.

Summary of contributions.

• To the best of our knowledge, we propose the first modelbased architecture learnable end to end for joint image alignment and super-resolution from raw image bursts. • We introduce a new differentiable image registration module that can be applied to images of different resolutions, is readily integrable in neural architectures, and may find other uses beyond super-resolution. • We show that our approach gives excellent results on both real image bursts (with up to ×4 upsampling for raw images) and synthetic ones (up to ×16 for RGB images).

Related Work

Classical multiframe super-resolution. Tsai and Huang wrote the seminal paper in this setting [START_REF] Dong | Image super-resolution using deep convolutional networks[END_REF], with a restoration model in the frequency domain assuming known translations between frames. Most latter approaches have focused on the spatial domain, and they generally fall into two main categories [START_REF] Köhler | Multi-frame super-resolution reconstruction with applications to medical imaging[END_REF]: In interpolation-based methods, LR snapshots aligned with sub-pixel precision are jointly interpolated into an HR image [15,38]. Impressive results have recently been obtained for hand-held cameras using the variant of this method proposed by Wronski et al. [START_REF] Wronski | Handheld multi-frame superresolution[END_REF], whose insight of exploiting aliasing effects has been one of the inspirations of our work. However, due to the sequential nature of their algorithm, errors may propagate from one stage to the next, leading to sub-optimal reconstructions [START_REF] Park | Superresolution image reconstruction: a technical overview[END_REF]. In contrast, iterative spatial domain techniques iteratively refine an estimate for the super-resolved image so as to best explain the observed LR frames under some image formation model. Variants of this approach include the early iterated backprojection algorithm of Irani et al. [START_REF] Irani | Improving resolution by image registration[END_REF], the maximum likelihood technique of Elad and Feuer [START_REF] Elad | Restoration of a single superresolution image from several blurred, noisy, and undersampled measured images[END_REF],

and the model regularized by bilateral total variation of Farsiu et al.

[13]. The image formation parameters are either be assumed to be known a priori through calibration, or estimated jointly with the HR image. In general, inter-frame motion can either be estimated separately, or be treated as an integral part of the super-resolution problem [START_REF] Bercea | Confidence-aware levenberg-marquardt optimization for joint motion estimation and super-resolution[END_REF]16], thus avoiding motion estimation between LR frames, whose accuracy may be affected by undersampling [START_REF] Vandewalle | Super-resolution from unregistered aliased images[END_REF]. The method proposed in the rest of this paper combines the best of both worlds since it performs joint estimation while aligning the LR frames with the reconstructed HR image.

Learning-based approaches. In this context, the multiframe case has received less attention than its single-image counterpart, for which several loss functions and architectures have been proposed [9, 28,[START_REF] Zhang | Deep unfolding network for image super-resolution[END_REF]. Most multi-frame algorithms focus on video super-resolution. Model-based techniques learn non-uniform interpolation or motion compensation using convolutional neural networks [START_REF] Shi | Real-time single image and video super-resolution using an efficient sub-pixel convolutional neural network[END_REF] but the most successful approaches so far are model free, leveraging instead diversity with 3D convolutions or attention mechanisms [START_REF] Jo | Deep video super-resolution network using dynamic upsampling filters without explicit motion compensation[END_REF][START_REF] Wang | Edvr: Video restoration with enhanced deformable convolutional networks[END_REF]. Learning-based methods have also been used in remote sensing applications, using 3D convolutions [32] or joint registration/fusion architectures [8] for example. Finally, and closer to our work, Bhat et al.

[4] have recently proposed a network architecture for raw burst super-resolution, together with a very interesting dataset featuring both synthetic and real images for training and testing. It is important to note that learning-based approaches to super-resolution are typically trained on synthetically generated LR images [START_REF] Lugmayr | Ntire 2020 challenge on real-world image super-resolution: Methods and results[END_REF], a strategy that may not generalize well to real photographs unless great care is taken in modeling the image corruption process [START_REF] Brooks | Unprocessing images for learned raw denoising[END_REF]. Learning super-resolution models from real LR/HR image pairs is quite challenging since it requires in general using separate cameras with different lenses and spatial resolution, with inevitable spatial and spectral misalignments. As shown by our experiments, our method, although trained from synthetic LR images, gives excellent results with real bursts taken from different smartphones and cameras. Leveraging real images at training time is, for now, left for future work.

Proposed Approach

This section presents the three main components of our approach: its image formation model, an optimization procedure for solving the corresponding inverse problem, and its unrolled implementation in a feedforward architecture whose parameters can be learned end to end.

Image formation model

Image acquisition in a digital camera starts from an instantaneous irradiance function f γ,t : [0, 1] 2 → R + defined on a continuous retinal domain with nonnegative values, such that f γ,t (u) is the spectral irradiance value at point u, time t, and wavelength γ, accounting for blur due to optics, atmospheric effects, etc. The camera sensor integrates f γ,t in the spatial, time, and spectral domains to construct a raw digital image y : [1, . . . , n] 2 → R + , where each pixel's spectral response is typically dictated by the 2 × 2 RGGB Bayer pattern, with twice as many measurements for the green channel than for the red and blue ones [START_REF] Kimmel | Demosaicing: image reconstruction from color ccd samples[END_REF]. Modern cameras turn the raw image y into a full blown, three-channel RGB image x with the same spatial resolution through an interpolation process called demosaicking.

In practice, we do not have access to f γ,t to use as ground truth for learning an image restoration process, even when an accurate model of the f γ,t → x map is available. Thus, we model instead the process x → y k , where x is a latent high-resolution (HR) image we wish to recover, and the low-resolution (LR) images y k (k = 1, . . . , K) have been observed in a burst of length K. We assume that x is sharp, without any blur, and noiseless. The burst images are obtained through the following forward model (Figure 3):

y k = DBW p k x + ε k for k = 1, . . . , K, (1)
where ε k is some additive noise. Here, both the HR image x and the frames y k of the burst are flattened into vector form. The operator W p k parameterized by p k warps x to compensate for misalignments between x and y k caused by camera or scene motion between frames, assumed here to be a 6parameter affine transformation of the image plane, then resamples the warped image to align its pixel grid with that of y k . Finally, the corresponding HR image is blurred to account for integration over both space (the LR pixel area, using either simple averaging or, as in the figure, a Gaussian filter) and time (accounting for camera and/or scene motion during exposure), and it is finally downsampled in both the spatial and spectral domains by the operator D, with an (a priori) arbitrary choice of where to pick the sample from (pixel corner or center for example), the spectral part correponding to selecting one of the three RGB values to assemble the raw image. It will prove convenient in the sequel to rewrite (1) as y = U p x + ε, where

U p =    DBW p1 . . . DBW p K   , y =    y 1 . . . y K   , p =    p 1 . . . p K   , ε =    ε 1 . . . ε K   . (2)
Before closing this section, let us note that simpler motion models with two (translation) or three (rigid motion) parameters, or (much) more complex piecewise-affine or elastic models could be considered depending on the application. We focus here on the scenario where a user wishes to zoom in on a relatively small crop (say, between 100 × 100 to 800 × 800 pixels) of a multi-megapixel image, and the affine model has proven effective with real handheld cameras in this setting. This implicitly corresponds to a globally piecewise-affine motion model.

Inverse problem and optimization

Given the image formation model of Eqs. ( 1)-( 2), recovering the HR image x from the K LR frames y k in the burst can be formulated as finding the values of x and p that minimize

1 2 y -U p x 2 + λφ θ (x), (3)
where φ θ is a parameterized regularizer, to be detailed later, and λ is a parameter balancing the data-fidelity and regularization terms. Many methods are of course available for minimizing this function. Like others (e.g., [24]), and mainly for simplicity, we choose here a quadratic penalty method [33, Sec. 17.1] often called half-quadratic splitting (or HQS) [START_REF] Geman | Nonlinear image recovery with half-quadratic regularization[END_REF]: the original objective is replaced by

E µ (x, z, p) = 1 2 y -U p z 2 + µ 2 z -x 2 + λφ θ (x), ( 4
)
where z is an auxiliary variable, and µ is a parameter increasing at each iteration, such that, as µ → +∞, the minimization of (4) with respect to x, z and p becomes equivalent to that of (3) with respect to x and p alone. Each iteration of HQS can be viewed as one step of a blockcoordinate descent procedure for minimizing E, changing one of variables z, x and p at a time while keeping the others fixed, with the value of µ increasing after each iteration.

Convergence guarantees for quadratic penalty methods require an approximate minimization of Eq. ( 4) with increasing precision over time [START_REF] Nocedal | Numerical Optimization[END_REF]. Following common practice in computer vision (e.g. [24]), we use HQS without formally checking that its precision indeed increases with iterations. This very simple procedure turns out to work well in practice. Its steps are detailed in the next three paragraphs, the exponent t being used to designate the value of the variables at iteration t. The sequence of weights (µ t ) t≥0 is learned end-to-end as explained in Section 3.3.

Updating z. Several strategies are possible for minimizing Eq. ( 4) with respect to z. Given the dimension of the problem, one may choose for instance a fast iterative minimization procedure such as conjugate gradient descent.

Since an approximate minimization is sufficient for our needs, we have chosen to use instead a single step of plain gradient descent, which converges more slowly in theory, but is also simpler and more easily amenable to the unrolled optimization strategy for end-to-end learning that will be presented next. The update at iteration t is given by

z t ← z t-1 -η t U p t-1 (U p t-1 z t-1 -y)+µ(z t-1 -x t-1 ) , (5)
where η t > 0 is some step size, also learned end to end.

Updating the motion parameters p. Let p k denote the part of the parameter vector p responsible for the alignment of z t and y k in (4). The corresponding optimization problem can be rewritten as

min p k 1 2 y k -DBW p k z t 2 . (6)
This is a non linear least-squares problem, which can once again be solved using many different techniques. Here, we pick a Gauss-Newton approach, which corresponds to a variant of the Lucas-Kanade algorithm [START_REF] Baker | Lucas-kanade 20 years on: A unifying framework[END_REF][START_REF] Bruce | An iterative image registration technique with an application to stereo vision[END_REF], showing again that a 40-year old technique can still be relevant today. Specifically, we perform one Gauss-Newton step at each iteration t for each p k in parallel:

p t k ← p t-1 k -J t k J t k -1 J t k r t k , (7)
where r t k = U p t-1 k z ty k is the residual of the non-linear least-squares problem (6), and J t k = (∂U p t-1 k /∂p k )z t is the Jacobian of the DBW p k operator. The only difference with a Lucas-Kanade iteration is the presence of a highresolution frame z t and the downsampling operator DB. This is similar to [16], or more recently [START_REF] Bercea | Confidence-aware levenberg-marquardt optimization for joint motion estimation and super-resolution[END_REF][START_REF] He | A nonlinear least square technique for simultaneous image registration and super-resolution[END_REF], which align high-resolution images with low-resolution ones.

Estimating the HR image x. The x update is obtained as

x t ← arg min x µ t-1 2 z t -x 2 + λφ θ (x),
which amounts to computing the proximal operator of the prior φ θ . In practice, we follow a "plug-and-play" approach [START_REF] Stanley H Chan | Plugand-play admm for image restoration: Fixed-point convergence and applications[END_REF]35,[START_REF] Singanallur V Venkatakrishnan | Plug-and-play priors for model based reconstruction[END_REF], and replace the proximal operator by a parametric function f θ (z t ) (here, a CNN, see implementation details). Using such an implicit prior has proven very effective in our setting. More traditional image priors such as total variation could of course have been used as well.

Unrolled optimization and backpropagation

The optimization procedure described so far requires choosing hyper-parameters such as the sequence (µ t ) t≥0 , and its implicit prior also involves model parameters θ. By using a training set of n LR burst/HR image pairs, we propose to learn all these parameters in a supervised fashion. We denote the training set by

(Y i , x i ) n i=1
, where Y i = {y i j } K j=1 is the i-th burst of LR images associated to the HR image x i . We then unroll the optimization procedure for T steps and, denoting by xT (Y i ) the HR image estimated from burst Y i , we consider the objective function

1 n n i=1 L(x T (Y i ), x i ), (8)
where L is the 2 or 1 loss (in practive we have observed that the 1 loss performs slightly better). Because every step of our estimation procedure is differentiable, we minimize (8) by stochastic gradient descent.

Learned data prior. Good image priors are essential for solving ill-posed inverse problems. As noted earlier, instead of using a classical one, such as total variation (TV) or bilateral total variation (BTV) [13], we learn an implicit prior parameterized by a convolutional neural network f θ in a data-driven manner. We use the ResUNet architecture introduced in [START_REF] Zhang | Deep unfolding network for image super-resolution[END_REF] in practice. It involves four scales, each of which has an identity skip connection between downscaling and upscaling operations.

Implementation details and variants

Downsampling and blurring operators D, B. We have tried different variants of downsampling/blurring strategies such as Gaussian smoothing. In practice, we have observed that simple averaging, which is differentiable and parameter-free, gives good results in all our experiments. As a consequence, we do not assume any knowledge about the blur used to generate data, corresponding to an operator B that only captures blur due to photon integration on the sensor without addressing optical blur. We argue that this limited model is relevant because modern cameras and smartphone are aliased [START_REF] Wronski | Handheld multi-frame superresolution[END_REF], which may explain the generalization to real images, as soon as the scene is static.

Initialization by coarse alignement. To initialize the motion parameters p, we cannot minimize [START_REF] Stanley H Chan | Plugand-play admm for image restoration: Fixed-point convergence and applications[END_REF] as in the previous section, because no good estimate of the HR image is available. Therefore, we align each LR frame to an arbitrary one from the burst (e.g., the first one) by using the Lucas-Kanade forward additive algorithm [START_REF] Baker | Lucas-kanade 20 years on: A unifying framework[END_REF]36] which is known to be robust to noise. Note that another difficulty lies in the raw format of images. To overcome this issue, we simply convert raw images into grayscale images by using bilinear interpolation. This is of course sub-optimal, but sufficient for obtaining coarse motion parameters.

Initialization via coarse-to-fine strategy. For extreme upsampling factors (×16), we found a coarse-to-fine initialization strategy to be useful: We initialize the motion parameters p 0 j and high-resolution image z 0 by using the output of the algorithm trained at a lower upsampling factor. For instance, ×16 can be obtained by applying twice a ×4 algorithm, or four times ×2 algorithm.

Experiments

Experiments were conducted on synthetic and real raw image bursts. We also provide experiments on RGB bursts in the appendix, allowing easier comparison with earlier approaches that cannot handle raw data.

Training procedure and data. For synthetizing realistic raw bursts from groundtruth RGB images, we follow the approach described in [4], using the author's publicly available code4 on the training split of the Zurich raw to RGB dataset [19]. The approach consists of applying the inverse RGB to raw pipeline introduced in [START_REF] Brooks | Unprocessing images for learned raw denoising[END_REF]. Displacements are randomly generated with Euclidean motions and frames are downscaled with bilinear interpolation in order to simulate LR frames containing aliasing. Synthetic, yet realistic, noise is added to the frames, and color values are discarded according to the Bayer pattern. Then, we train our models for minimizing the loss (8). We perform 100 000 iterations of the ADAM optimizer with a batch size of 10, a burst size of 14 and with a learning rate of 3×10 -5 decaying by a factor 2 after 50 000 iterations. Our approach is implemented in Pytorch and takes approximately 1.5 days to train on an Nvidia Titan RTX GPU. We evaluate our models in all our experiments with a burst size of 14 unless specified.

Extreme ×16 upsampling on RGB images. As a proofof-concept, we also perform experiments for an unusual ×16 super-resolution task, using the coarse-to-fine strategy of Sec. 3.4. A result is presented in Fig. 2, showing impressive reconstruction and additional ones can be found in the appendix. Even though not realistic, we believe the exper-iment to be of interest, as it demonstrates the effectiveness of our approach in an idealistic, yet extreme, setting. Evaluation on synthetic RAW images. The evaluation protocol of [4] allows us to perform quantitative comparison with their state-of-the-art method for processing raw image bursts. An additional comparison with [START_REF] Wronski | Handheld multi-frame superresolution[END_REF] would have been interesting but this method is part of a commercial product that could not be shared with us. with synthetic affine motions. Reconstruction error in average PSNR and geometrical registration error in pixels for our models. "known p" is the oracle performance our model could achieve, if motion estimation was perfect.

Demosaic+SISR ETH [4] Ours

We provide a quantitative comparison in Table 1 with the model introduced in [4], as well as a single-image upsampling baseline based on the ResUNet architecture [START_REF] Zhang | Deep unfolding network for image super-resolution[END_REF], which we use as a plug-and-play prior in our model.

To that effect, we first use the validation set of [4] available online (with no overlap with the training set), for which motions are unknown, allowing us to compare with their method, which we outperform by more than 2dBs. In order to perform further comparison and conduct the ablation study, we also build an additional validation set by randomly extracting 266 images from the Zurich raw to RGB dataset, allowing us to generate validation data with known motion. We evaluate variations of our model in the same table, notably comparing the registration accuracy achieved by these variants by using the geometrical error presented in [36]. More precisely, we perform a small ablation study by introducing a simpler baseline that does not perform joint alignement and only exploits the coarse registration module (no refine baseline). Performing motion refinement significantly improves the registration accuracy and subsequently the image reconstruction quality. Last, we also report the oracle performance of our model with known motions.

We provide a visual comparison in Figure 4 with singleimage SR baselines and the state-of-the-art method [4] for processing raw image bursts. Only the two approaches processing bursts are able to recover high-frequency details, demonstrating their ability to leverage and remove aliasing artefacts, which are very present in the top image. Significantly better quality results are obtained with our approach.

Impact of burst length and cropping size. The dataset Zurich rgb-to-raw [19] was very useful for training our models, but it unfortunately features relatively small image crops of size 96 × 96 without giving access to the original megapixel images. By experimenting with real raw data, it became apparent to us to our method was performing better with larger crops (e.g., more than 200 × 200 pixels), achieving better registration and visually better results. To study the impact of the crop size and burst length, we have thus synthesized additional raw bursts from the DIV2K dataset, and report our experimental results in Figure 5, confirming our findings. Note that this does not appear to be a strong limitation of our approach, since in real-life scenarios, we can always assume that the original megapixel image is available. As expected, the performance of our approach is also increasing with the burst size, even though our models were trained with bursts of size 14.

Results on real raw image bursts, dataset of [4]. In Figure 6, we show a comparison with [4] using their dataset featuring small crops of size 96 × 96. As discussed previously, this setup is suboptimal for our approach, but still produces visually pleasant results. Choosing which method performs best here is however very subjective and we found conclusions hard to draw on this dataset. Whereas the im- ages produced by [4] may sometimes look slightly sharper, one may argue that our approach seems to recover more reliable details, e.g., the text is perhaps easier to read. Note that our models were trained on synthetic data only and we leave fine-tuning with real data on this dataset for future work. There is an attempt in [4] to address the open problem of quantitative evaluation with real data using a custom metric, but, like any other attempt so far, it is flawed since (i) it is based on the alignment method of [4], with an unavoidable slight bias in its favor, and (ii) it assumes ground truth from a particular Canon camera. Interestingly, this score improvement does not always correlate with visual quality, as shown by Figure 6. This is by no means a criticism of [START_REF] Bhat | Ntire 2021 challenge on burst super-resolution: Methods and results[END_REF]: we believe instead that quantitative evaluation on real images is an extremely challenging problem, far from being solved. Since the submission of our paper, the results of the NTIRE 2021 burst super-resolution challenge have been published [START_REF] Bhat | Ntire 2021 challenge on burst super-resolution: Methods and results[END_REF]. Our method ranked third quantitatively in the "synthetic data" part of the challenge that we entered. Results on real raw image bursts from various devices. Finally, we demonstrate the effectiveness of our approach on real raw bursts acquired by different devices. We consider a Panasonic Lumix GX9 camera, which is interesting for SR as it does not feature an optical anti-aliasing filter, a Canon Powershot G7X camera, a Samsung S7 and a Pixel 4a smartphones. Results obtained in high noise regimes have already been presented in Figure 1, showing that our approach is surprisingly robust to noise. We believe that the result is of interest since it may allow photographers to use high ISO settings in low-light conditions, without sacrificing image quality. Other results are presented in Figure 6 on low-noise outdoor conditions with bursts of 20 to 30 raw images. In all cases, the method succeeds at recovering high-frequency details. Many more examples and comparisons with other multiframe methods are provided in the supplementary material. We also present failure cases, corresponding in large parts to scene motion. Last, we remark that our method is relatively fast at inference time.

Processing a burst of 20 raw 300 × 300 images takes for instance about 1s on an Nvidia Titan RTX GPU, producing an upsampled image of size 1200 × 1200.

Conclusion

We have presented a simple but effective method for superresolution that combines the interpretability of modelbased approaches to inverse problems with the flexibility of data-driven architectures and can be learned from pairs of synthetic LR and real HR images. We plan several extensions, including using multiple cameras to add real LRburst/HR-image pairs to the training mix, and at test time to take advantage of the multiplicity of imaging devices now available on high-end smartphones. This will open the door to wide-baseline super-resolution applications, such as the construction of high quality panoramas and finely detailed texture maps in multi-view stereo reconstructions. Finally, we plan to explore several other extensions of our approach, including tackling blurry bursts, extending super-resolution to reconstruct HDR images, and pursuing applications in the astronomy and microscopy domains

Figure 2 :

 2 Figure 2: Proof of concept for extreme ×16 upsampling. The right image is obtained by processing a burst of 20 LR images presented on the left obtained with synthetic random affine movements and bilinear downsampling.

Figure 3 :

 3 Figure 3: Image formation:The HR image x is warped then resampled to align it with the LR image y using the operateur Wp k . It is then blurred by the operator B to account for integration over LR pixels and finally downsampled in the spatial and spectral domains by the operator D (the spectral downsampling from RGB to R, G or B is not illustrated here for simplicity).

Figure 4 :

 4 Figure 4: Visual comparison on synthetic raw image bursts used in [4]. Demosaic+SISR is our single-image baseline based on the ResUNet architecture [48] (see main text). The two right columns are produced by methods dedicated to raw burst processing, respectively [4] and ours.

Figure 5 :

 5 Figure 5: Left: Impact of the crop size on the registration and reconstruction performance. Right: Impact of the burst length, see main text for details.

Figure 6 :

 6 Figure 6: Results from real raw bursts from dataset of [4] including Aligned PSNR score (see main text).

Figure A6 :Figure A7 :

 A6A7 Figure A6: Visual differences caused by merging a different number of frames in the case of low SNR scenes. With a larger number of frames we can observe a quality increase and better denoising.

Table 1 :

 1 Results with synthetic raw image bursts of 14 images generated from the Zurich raw to RGB dataset[19]

	Method	PSNR (db) Geom (pix) SSIM
	Scores on public validation set
	ETH [4]	39.09	-	-
	Ours (refine)	41.45	-	0.95
	Scores on our own validation set to conduct the ablation study
	Bicubic Single Image	33.45	-	-
	Multiframe L2 only	34.21	-	-
	Multiframe L2 + TV prior 34.48	-	-
	Single Image	36.80	-	-
	Ours (no refinements)	40.38	0.55	0.958
	Ours (refinements)	41.30	0.32	0.963
	Ours (known motion)	42.41	0.00	0.971

There is, however, a trend today toward removing these filters, as in the prosumer camera used in some of our experiments with real images.

https://github.com/goutamgmb/NTIRE21_BURSTSR.

Proceedings of Imaging Understanding Workshop, 1981. 2,

Acknowledgments

We thank Frédéric Guichard for useful discussions and comments. This work was funded in part by the French government under management of Agence Nationale de la Recherche as part of the "Investissements d'avenir" program, reference ANR-19-P3IA-0001 (PRAIRIE 3IA Institute). JM and BL were supported by the ERC grant number 714381 (SOLARIS project) and by ANR 3IA MIAI@Grenoble Alpes (ANR-19-P3IA-0003). JP was supported in part by the Louis Vuitton/ENS chair in artificial intelligence and the Inria/NYU collaboration. This work was granted access to the HPC resources of IDRIS under the allocation 2020-AD011011252 made by GENCI. Panasonic Panasonic Canon G7X Samsung s7 Pixel3a Pixel4a Source image ISP Camera (jpeg output)

Joint demosaic+single-image SR Our burst processing method Figure 7: Results from real raw image bursts obtained with various cameras. We provide comparisons with single image and multiframe baselines. Finest restored details can be seen by zooming on a computer screen. The last three digits of the phone number, only legible in our reconstruction, are masked in the figure for privacy concerns.

Supplementary material

This supplementary material presents additional qualitative and quantiative results. In Figure A1 we present additional visual comparison with two burst denoising methods on real images. In Table A1 we present additional experiments on RGB images. Figures A2 andA3 are devoted to super-resolution experiments from real raw data from different smartphones (Google Pixel 3a and 4a, Samsung S7 and S10) and cameras (Panasonic Lumix GX9 and Canon Powershot G7X) and comparison with additional baselines. In Figures A4, we present extreme upsampling results by using synthetic RGB image bursts. In Figure A5, we present restoration results obtained from real images with very low SNR to illustrate the efficiency of our method to perform blind denoising. In Figures A6 andA7, we study the effect of the number of frames in the burst on the reconstruction, both in the low SNR and high SNR settings. Finally, we present failure cases in Figure A8, where fast moving objects are present in the scene.

Comparison with burst denoising methods. We perform additional qualitative comparison on a real image with two burst denoising methods. We compare our method with [START_REF] Ehret | Joint demosaicking and denoising by fine-tuning of bursts of raw images[END_REF] which performs joint denoising and demosaicking on a burst of raw images. We use the code and the pretrained model made available online. We also use the code and pretrained model of [START_REF] Xia | Basis prediction networks for effective burst denoising with large kernels[END_REF]. However the model is only designed to perform grayscale burst denoising, so we perform denoising independently on each RGB channel and then perform demosaicking to get an RGB image. Despite our best efforts for tuning the parameters of these methods to maximize visual quality, the results obtained are not as good as our method (see Figure A1 below). We believe this is not surprising since each one of these methods only addresses a subset of our problem. Adapting them successfully to our general setting is not trivial. Evaluation on RGB images. We compare our approach on the BSD68 dataset against state-of-the-art single-image and video superresolution algorithms (considering a burst as a video sequence) and report the HR image reconstruction accuracy in terms of average PSNR in Table A1. For the training with RGB data, we perform 80 000 iterations of the ADAM optimizer with a batch size of 10, a burst size of 14 and with a learning rate of 3 × 10 -5 decaying by a factor 2 after 40 000 iterations. For evaluating the model VSR-DUF [START_REF] Jo | Deep video super-resolution network using dynamic upsampling filters without explicit motion compensation[END_REF], we use the code and the pretrained models made available online by the authors. Other single-image reconstruction results are from [START_REF] Zhang | Deep unfolding network for image super-resolution[END_REF]. In the present setting, we consistently outperform other baselines, notably demonstrating that burst SR cannot simply be addressed effectively by current video SR approaches. We also note that our models perfom better with less blurring (and more aliasing). Finally, we evaluate variations of our model in the same table, notably comparing the registration accuracy achieved by these variants by using the geometrical error presented in [36]. More precisely, we perform a small ablation study by introducing a simpler baseline that does not perform joint alignement and only exploits the coarse registration module (no refine baseline). Performing joint alignment and image estimation systematically improves motion estimation. Last, we also report the oracle performance of our model with known motions. + EDSR [28] groupsc [START_REF] Lecouat | Fully trainable and interpretable non-local sparse models for image restoration[END_REF] + VSR-DUF [START_REF] Jo | Deep video super-resolution network using dynamic upsampling filters without explicit motion compensation[END_REF] Joint mosa & SR Ours Full image Camera ISP groupsc [START_REF] Lecouat | Fully trainable and interpretable non-local sparse models for image restoration[END_REF] + EDSR [28] groupsc [START_REF] Lecouat | Fully trainable and interpretable non-local sparse models for image restoration[END_REF] + VSR-DUF [START_REF] Jo | Deep video super-resolution network using dynamic upsampling filters without explicit motion compensation[END_REF] Joint mosa & SR Ours

Method