SINGLE-SHOT HYPERSPECTRAL WAVEFRONT SENSOR

Abstract

A method for determining wavefront shapes of a multi-spectral signal light beam from a single signal image acquisition of said multi-spectral signal beam, with a device including an optical assembly made at least of an optical mask and an imaging sensor, notably a matrix imaging sensor, for generating and recording intensity patterns of incident beams, by having these beams reflect on, or propagate through, the optical mask. The optical mask having the optical properties: i) to cause the intensity pattern to depend on the wavefront shape, so that a tilt applied to the wavefront shape results in a displacement amount of the intensity pattern, and ii) to produce uncorrelated intensity patterns over at least one surface area A of the imaging sensor, for a plurality of respective incident monochrome beams of different wavelengths having a same wavefront shape.

Claims

1. A method for determining each wavefront shape of N spectral channels of a multi-spectral signal light beam from a single signal image I(x,y) acquisition of said multi-spectral signal light beam, with a device comprising an optical assembly made at least of an optical mask and an imaging sensor for generating and recording intensity patterns of incident beams, by having these beams reflect on, or propagate through, the optical mask, the optical mask having the optical properties: i) to cause the intensity pattern to depend on the wavefront shape, so that a tilt applied to the wavefront shape results in a displacement amount of the intensity pattern, ii) to produce uncorrelated intensity patterns over at least one surface area A of the imaging sensor, for a plurality of respective incident monochromatic beams of different wavelengths having a same wavefront shape, two uncorrelated random intensity patterns being defined as statistically orthogonal relatively to a zero-mean cross-correlation product, the method comprising: a) recording several reference intensity patterns R.sub.L(x,y) using the device, each reference intensity pattern R.sub.L(x,y) being generated by a respective reference incident monochromatic beam L with wavelengths L, L varying from 1 to N, with N being the number of different reference incident monochromatic beams, x and y being coordinates; b) recording one single signal image I(x,y) of the intensity pattern generated by the said multi-spectral signal light beam which comprises at least the N wavelengths, using the device, the single signal image I(x,y) being representative of light impinging on the at least one surface area (A); c) computing intensity-weight data W.sub.L.sup.I(x,y) and deformation data T.sub.L.sup.I(x,y), for all L varying from 1 to N, the intensity-weight data W.sub.L.sup.I (x,y) and the deformation data T.sub.L.sup.I(x,y) being representative of an intensity modulation and a diffeomorphism, respectively, of each given reference intensity pattern R.sub.L(x,y), at wavelength λ.sub.L, for the single signal image I(x,y), all the N intensity-weight data W.sub.L.sup.I (x,y) and the N deformation data T.sub.L.sup.I(x,y) being computed, for L varying from 1 to N, so as to minimize, for all sampling points (x,y) of the surface area A, from the single signal image I(x,y): a difference D.sub.A between the single signal image I(x,y) on the one hand, and the sum of reference intensity patterns R.sub.L multiplied by intensity-weight data W.sub.L.sup.I (x,y) and deformed by deformation data T.sub.L.sup.I(x,y), on the other hand: $D_A={.Math.I\left(x,y\right)-\underset{L=1}{\overset{N}{.Math.}}{W}_L^I\left(x,y\right)R_L\left[\left(x,y\right)+T_L^I\left(x,y\right)\right].Math.}_A$ the symbol ∥.∥.sub.A designating a norm calculated for all (x,y) sampling points in the surface area A; for the surface A, each given reference intensity patterns R.sub.L(x,y) being orthogonal to each reference intensity pattern R.sub.K(x,y) relatively to the zero-mean cross-correlation product, when K entire natural different from L and chosen between [1; N]; d) generating data for each wavelength λ, representative of: the shape of the wavefront by integrating the deformation data T.sub.L.sup.I(x,y), the intensity map based on the weight W.sub.L.sup.I(x,y).

2. The method according to claim 1, wherein the imaging sensor is a matrix imaging sensor, and in the step c, the surface A=ΣA.sub.i, with A.sub.i being the surface of a macropixel, computing intensity-weight data W.sub.L.sup.I(A.sub.i) and deformation data T.sub.L.sup.I(A.sub.i) which are constant for the macropixel of surface A.sub.i, for all L varying from 1 to N, by minimizing the difference D.sub.Ai updating the intensity-weight data W.sub.L.sup.I(x,y) and the deformation data T.sub.L.sup.I(x,y) for all coordinates (x,y) on the surface A, with W.sub.L.sup.I(x,y)=W.sub.L.sup.I(A.sub.i) and T.sub.L.sup.I (x,y)=T.sub.L.sup.I(A.sub.i) for all (x,y) belongs to A.sub.i

3. The method according to claim 1, wherein the optical mask is an engineered pseudo-diffuser, the said engineered pseudo-diffuser having the property to produce intensity patterns that are exactly orthogonal relatively to the zero-mean cross-correlation product over at least one surface area A, for a plurality of respective incident monochromatic beams of different wavelengths having a same wavefront shape.

4. The method according to claim 1, wherein the optical mask is a diffuser or a diffractive optical element.

5. The method according to claim 1, wherein the optical mask is a dispersive optical fiber bundle, a metasurface, or a freeform optical element.

6. The method according to claim 2, wherein the sub-image is a macro-pixel whose size is larger or equal to N pixels of the matrix imaging sensor.

7. The method according to claim 3, where the surface area A is in the range from $A=\frac{1}{10}{N\left(\frac{\lambda }{\theta }\right)}^2\mathrm{to}A=10{N\left(\frac{\lambda }{\theta }\right)}^2,$ θ being a scattering or diffracting angle of the optical mask (14), N the number spectral channels and λ a central wavelength, an average wavelength or a median wavelength of the multispectral signal light beam.

8. The method according to claim 2, wherein an estimate of T.sub.L.sup.I(A.sub.i) and W.sub.L.sup.I(A.sub.i) and of the minimum of the difference D.sub.Ai is obtained by computing zero-mean cross-correlation product images for each macropixel of surface Ai, D.sub.Ai being defined relatively to the norm ∥.∥.sub.Ai, between the signal sub-image I(x,y) and each of the reference intensity patterns R.sub.L(x,y) with (x, y) coordinates of the portion A.sub.i, the zero-mean cross-correlation product image between the signal sub-image and each of the reference intensity patterns having a peak, the intensity-weight data W.sub.L.sup.I(A.sub.i) being the amplitude of the peak and the deformation data T.sub.L.sup.I(A.sub.i) being the displacement vector between the said peak or its centroid, from the center of the zero-mean cross-correlation product image.

9. The method according to claim 2, wherein an estimate of T.sub.L.sup.I(A.sub.i) and W.sub.L.sup.I(A.sub.i) and of the difference D.sub.A is obtained by computing a Wiener deconvolution for each macropixel of surface Ai, D.sub.Ai being defined relatively to the norm ∥.∥.sub.Ai, the Wiener deconvolution possibly relating to: a Wiener deconvolution of the signal sub-images by the reference intensity patterns R.sub.L(x,y), or a Wiener deconvolution of the reference intensity patterns R.sub.L(x,y) by the signal sub-images.

10. The method according to claim 2, wherein an estimate of T.sub.L.sup.I(A.sub.i) and W.sub.L.sup.I(A.sub.i) and of the difference D.sub.Ai is obtained for each macropixel of surface A.sub.i, by computing a matrix inversion algorithm, D.sub.Ai being defined relatively to the norm ∥.∥.sub.Ai, wherein the matrices to be inverted being related to at least a transform of a sub-image of an intensity pattern R.sub.L(x,y), such as a Fourier transform.

11. The method according to claim 2, wherein an estimate of the difference D.sub.A is obtained by computing intensity-weight data W.sub.L.sup.I(x,y) and deformation data displacement data T.sub.L.sup.I(x,y) thanks to an iterative optimization algorithm or thanks to a compressed sensing algorithm.

12. The method according to claim 1, wherein an estimate of T.sub.L.sup.I(A.sub.i) and W.sub.L.sup.I(A.sub.i) and of D.sub.Ai is computed thanks to deep learning approaches using neuronal networks which minimize D.sub.Ai with a collection of signal images, references images with respective N wavefronts known, D.sub.Ai being defined relatively to the norm ∥.∥.sub.Ai.

13. The method according to claim 1, wherein intensity-weight data W.sub.L.sup.I(x,y) and deformation data T.sub.L.sup.I(x,y) are computed thanks to the norm ∥.∥.sub.A which is a genetic algorithm with a collection of would-be solutions are tested, evaluated, selected for minimizing the value of D.sub.A.

14. The method according to claim 1, wherein the multi-spectral signal light beam is generated by a broadband laser system, notably a chirped pulse amplifier laser system.

15. The method according to claim 1, where an interferometer is added before the optical assembly in order to discretize the spectrum of the multispectral signal light beam, the said multispectral signal light beam potentially exhibiting a continuous spectrum.

16. The method according to claim 14, where a spectral phase of the light beam being known at least on a fraction of the surface area A of the imaging sensor, the method comprising the computing of the spectral phase for at least one surface fraction of the imaging surface other than the said fraction of the surface area A, in particular on all the surface area A, the method comprising the computation of a temporal profile of the multispectral signal light beam on the surface area A of the imaging sensor via a transform in the spectral/time domain, such as a Fourier transform, based on the spatial distribution of the spectral phase.

17. The method according to claim 1, wherein measured wavefront shapes are characteristic of some chromatic aberrations introduced by an optical system.

18. A wavefront sensor for determining each wavefront shape of N spectral channels of a multi-spectral signal light beam from a single signal image I(x,y) acquisition of said multi-spectral signal light beam, comprising an optical assembly made at least of: an optical mask and an imaging sensor for generating and recording intensity patterns of incident beams, by having these beams reflect on, or propagate through, the optical mask, the optical mask having the optical properties: i) to cause the intensity pattern to depend on the wavefront shape, so that a tilt applied to the wavefront shape results in a displacement amount of the said intensity pattern, ii) to produce uncorrelated intensity patterns over at least one surface area (A) of the imaging sensor, for a plurality of respective incident monochromatic beams of different wavelength having a same wavefront shape, two uncorrelated random intensity patterns being defined as statistically orthogonal relatively to a zero-mean cross-correlation product, the imaging sensor recording: (a) several reference intensity patterns R.sub.L(x,y), each reference intensity pattern R.sub.L(x,y) being generated by having a respective reference incident monochromatic beam L with wavelengths λ.sub.L, reflect on or propagate through the optical mask, L varying from 1 to N, with N being the number of different reference incident monochromatic beams, x and y being coordinates; b) one single signal image I(x,y) of the intensity pattern generated by having the multi-spectral signal light beam which comprises at least the N wavelengths λ.sub.L reflect on or propagate through the optical mask, the single signal image I(x,y) being representative of light impinging on the at least one surface area; computing means for: c) computing intensity-weight data W.sub.L.sup.I(x,y) and deformation data T.sub.L.sup.I(x,y), for all L varying from 1 to N, the intensity-weight data W.sub.L.sup.I(x,y) and the deformation data T.sub.L.sup.I(x,y) being representative of an intensity modulation and a diffeomorphism, respectively, of each given reference intensity pattern R.sub.L(x,y), at wavelength λ.sub.L, for the single signal image I(x,y), all the N intensity-weight data W.sub.L.sup.I(x,y) and the N deformation data T.sub.L.sup.I(x,y) being computed, for L varying from 1 to N, so as to minimize, for all sampling points (x,y) of the surface area A, from the single signal image I(x,y): a difference D.sub.A between the single signal image I(x,y) on the one hand, and the sum of reference intensity patterns R.sub.L multiplied by intensity-weight data W.sub.L.sup.I(x,y) and deformed by deformation data T.sub.L.sup.I(x,y), on the other hand: $D_A={.Math.I\left(x,y\right)-\underset{L=1}{\overset{N}{.Math.}}{W}_L^I\left(x,y\right)R_L\left[\left(x,y\right)+T_L^I\left(x,y\right)\right].Math.}_A$ the symbol ∥.∥.sub.A designating a norm calculated for all (x,y) sampling points in the surface area A; for the surface A, each given reference intensity patterns R.sub.L(x,y) being orthogonal to each reference intensity pattern R.sub.K(x,y) relatively to the zero-mean cross-correlation product, when K entire natural different from L and chosen between [1; N]; d) generating data for each wavelength λ.sub.L representative of: the shape of the wavefront by integrating the deformation data T.sub.L.sup.I(x,y), the intensity map based on the weight W.sub.L.sup.I (x,y).

19. The waveform sensor according to claim 18, wherein the imaging sensor is a matrix imaging sensor, the computing means are configured to: in the step c, compute intensity-weight data W.sub.L.sup.I(A.sub.i) and deformation data T.sub.L.sup.I(A.sub.i) which are constant for the macropixel of surface A.sub.i, for all L varying from 1 to N, by minimizing the difference D.sub.Ai, the surface A=ΣAi, with A.sub.i being the surface of a macropixel, update the intensity-weight data W.sub.L.sup.I(x,y) and the deformation data T.sub.L.sup.I(x,y) for all coordinates (x,y) on the surface A, with W.sub.L.sup.I(x,y)=W.sub.L.sup.I(A.sub.i) and T.sub.L.sup.I(x,y)=T.sub.L.sup.I(A.sub.i) for all (x,y) belongs to A.sub.i

20. The waveform a sensor according to a claim 1, wherein the said engineered pseudo-diffuser having the property to produce intensity patterns that are exactly orthogonal relatively to the zero-mean cross-correlation product over at least one surface area A, for a plurality of respective incident monochromatic beams of different wavelengths having a same wavefront shape.

21. The wavefront sensor according to claim 17, wherein the optical mask is a diffuser or a diffractive optical element.

22. The wavefront sensor according to claim 17, wherein the optical mask is a dispersive optical fiber bundle, a metasurface, or a freeform optical element.

23. The wavefront sensor according to claim 19, wherein the computing means are configured to estimate the minimum T.sub.L.sup.I(A.sub.i) and W.sub.L.sup.I(A.sub.i) and the minimum of the difference D.sub.Ai which is obtained by computing zero-mean cross-correlation product images for each macropixel of surface Ai, between the signal sub-image I(x,y) and each of the reference intensity patterns R.sub.L(x,y) with (x,y) coordinates of the portion A.sub.i, the zero-mean cross-correlation product image between the signal sub-image and each of the reference intensity patterns having a peak, the intensity-weight data W.sub.L.sup.I(A.sub.i) being the amplitude of the peak and the deformation data T.sub.L.sup.I(A.sub.i) being the displacement vector between the said peak or its centroid, from the center of the zero-mean cross-correlation product image.

24. The wavefront sensor according to claim 17, wherein the computing means are configured to estimate of the difference D.sub.A is obtained by computing intensity-weight data W.sub.L.sup.I(x,y) and deformation data T.sub.L.sup.I(x,y) with: a Wiener deconvolution of the signal sub-images by the reference intensity patterns R.sub.L(x,y), or a Wiener deconvolution of the reference intensity patterns R.sub.L(x,y) by the signal sub-images.

25. The wavefront sensor according to previous claim 17, wherein the computing means are configured to estimate of the difference D.sub.A is obtained by computing intensity-weight data W.sub.L.sup.I(x,y) and deformation data T.sub.L.sup.I(x,y) with matrix inversion algorithm, where the matrices to be inverted being related to at least a transform of a sub-image of an intensity pattern R.sub.L(x,y), such as a Fourier transform.

26. The wavefront sensor according to claim 17, wherein the computing means (17) are configured to estimate the difference D.sub.A is obtained by computing intensity-weight data W.sub.L.sup.I(x,y) and deformation data displacement data T.sub.L.sup.I(x,y) with an iterative optimization algorithm or thanks to a compressed sensing algorithm.

27. An optical device with: a wavefront sensor according to claim 17; a light emitting source for generating a multi-spectral signal light beam with specific wavelengths which is a broadband laser system, notably a chirped pulse amplifier laser system.

28. Use of a multi-spectral wavefront sensor for a multi-spectral signal light beam according to any one of claim 17, in: Optical metrology; Laser metrology; Quantitative phase microscopy (Biology, chemistry . . . ); Adaptive optics; Ophthalmology; Optical thermometry

Description

BRIEF DESCRIPTION OF THE DRAWINGS

[0134] The accompanying drawings, which are incorporated in and constitute a part of the specification, illustrate presently preferred embodiments of the invention, and together with the description, serve to explain the principles of the invention. In the drawings:

[0135] FIG. 1a, 1b, 1c and 1d are schematic illustrations of different embodiments of an optical system according to the invention,

[0136] FIGS. 2 to 6 are schematic illustrations of optical mask properties according to the invention,

[0137] FIG. 7 is a schematic of one embodiment of a method according to the invention, and

[0138] FIGS. 8 and 9 show examples of reconstruction of the shapes of wavefront using a method according to the invention.

[0139] Whenever possible, the same reference numbers will be used throughout the drawings to refer to the same or like parts.

[0140] In accordance with the invention, and as broadly embodied, in FIGS. 1a, 1b and 1c, an optical device 15 is provided.

[0141] The optical device 15 comprises a light source 11 for generating a multi-spectral signal light beam 12. Said multi-spectral signal light beam may potentially exhibits a continuous spectrum.

[0142] The optical device 15 further comprises a wavefront sensor 10. The latter includes device 16 comprising an optical mask 14 and an imaging sensor 18 for generating and recording intensity patterns of incident beams. For calculation purpose, computing means are also provided. The computing means may comprise a memory storing, a set of instructions and processor for executing said set of instructions, such as a PC computer or dedicated microcontroller bound. For purposes of illustration, the computing means are shown as block 17, which may be part of the wavefront sensor 10 or separate from the latter.

[0143] As shown in FIGS. 1a, 1b, 1c and 1d, the optical mask 14 is illuminated by the light beam 12 presenting a wavefront 13 and originated from a light source 11. Depending on the application of the optical device 10, the source of light 12 may be a broadband laser system, notably a chirped pulse amplifier laser system. In some embodiments, an interferometer may be added on the path of the beam 12, between the light source 11 and the wavefront sensor 10, in order to discretize the spectrum of the multispectral light beam 12.

[0144] FIG. 1a displays an embodiment of the wavefront sensor 10 where the optical mask is working in transmission preferably at the close vicinity of the imaging sensor 18, for example at a distance d, ranging from 0 to 10 L/θ, where L is the size of the imaging sensor and θ the scattering angle of the optical mask. In the illustrated example, the axis of incidence of the light beam is perpendicular to the optical mask. In a variant, the axis of incidence of the light beam forms an angle smaller than 90° with the optical mask.

[0145] When the incident light beam 12 with input wavefront 13 passes through the optical mask, the latter is scattered. Optical mask crossing causes the intensity pattern to depend on the wavefront shape, so that a tilt applied to the wavefront shape results in a displacement amount of the said intensity pattern, as can be seen in FIG. 2. Such property is also commonly referred to as the “memory effect” of the optical mask. This assumption may be valid until a so-called memory effect angle. The latter may range from 0° to 45°.

[0146] FIGS. 1b and 1d show an alternative embodiment where the optical mask 14 is working in reflection, preferably at a close vicinity of the imaging sensor. In this case the light beam 12 is incident along the parallel direction to an axis of incidence Z′ of the reflective optical mask 14, and is deviated towards the imaging sensor along an axis Z.

[0147] The light beam 22 emerging from the optical mask 14 is then captured using the imaging sensor 18, yielding an intensity pattern image I(x,y). The intensity pattern of emerged light beam 22 is a superposition of individual intensity patterns weighted by their respective contributions to the signal beams constituting light beam 12. An example of a signal intensity pattern generated with the optical mask 14 is shown in FIG. 3. In FIG. 3, the signal intensity pattern corresponds to a speckle pattern generated by a diffuser.

[0148] The imaging sensor 18 may be a monochromatic sensor.

[0149] In a preferred embodiment, the imaging sensor 18 is a matrix imaging sensor.

[0150] The imaging sensor 18 may be, for example, a CCD, a CID or a CMOS imaging sensor.

[0151] The wavefront sensor 10 may include more than one imaging sensor. In a such case, the intensity pattern image may be produced by a single imaging sensor or several imaging sensors.

[0152] Referring to FIGS. 1c and 1d, other embodiments of a wavefront sensor 10 are shown. In these embodiments, an intermediate optical system 19 is placed between the optical mask 14 and the imaging sensor 18. The optical mask 12 is placed at a distance d from the conjugate plan 21 defined by the intermediate optical system.

[0153] The intermediate optical system 19 may comprise at least one lens. For example, as illustrated in FIGS. 1c and 1d, the intermediate optical system consists of a a set of relay lenses 19a; 19b.

[0154] The device 16 comprises the optical mask 14, the imaging sensor 18 and the intermediate optical system 19 in FIGS. 1c and 1 d.

[0155] Preferably, the optical mask 14 is a diffuser or an engineered diffuser. An example of a suitable diffuser is described in the article (Berto, P., Rigneault, H., & Guillon, M. (2017). Wavefront sensing with a thin diffuser. Optics letters, 42(24), 5117-5120).

[0156] Otherwise, the optical mask 14 may be any optical mask: [0157] i) causing the intensity pattern to depend on the wavefront shape as mentioned above, and

[0158] producing uncorrelated intensity patterns over at least one surface area A of the imaging sensor, for a plurality of respective incident monochromatic beams of different wavelength having a same wavefront shape. The property ii) is evaluated relatively to a given measure of similarity.

[0159] The optical mask 14 can be: [0160] a ground glass or holographic diffuser [0161] a engineered pseudo-diffuser, the said engineered pseudo-diffuser having the property to produce intensity patterns that are exactly orthogonal relatively to the zero-mean cross-correlation product over at least one surface area A, for a plurality of respective incident monochromatic beams of different wavelengths having a same wavefront shape; [0162] a diffractive optical element; [0163] a dispersive optical fiber bundle, a metasurface, or a freeform optical element.

[0164] A measure of similarity between two intensity patterns can be mathematically characterized using correlation tools. For example, a measure may consist in localizing critical points of the intensity such as intensity maxima and to look for correlations between intensity maxima localization of the two intensity patterns. This can be quantified using a Pearson correlation coefficient. It can also be quantified using the mean average distance to the nearest neighbor weighted by the critical point density.

[0165] A measure of similarity between two intensity patterns may also be estimated by computing a zero-mean cross-correlation product. A statistical average of the zero-mean cross-correlation product between two uncorrelated random signals is zero. In contrast, a statistical average of the zero-mean cross-correlation product of an intensity pattern with itself is a function admitting an extremum, in particular a maximum. An illustration of the zero-mean cross-correlation product of the speckle of FIG. 3 with itself (also called auto-correlation product) is shown in FIG. 4. The result displayed in FIG. 4 exhibits an intense central bright spot demonstrating the high degree of correlation of an intensity pattern with itself.

[0166] Alternatively, a Wiener deconvolution can be used rather than a zero-mean cross-correlation product. The result of a Wiener deconvolution applied to an intensity pattern with itself is similar to what is displayed in FIG. 4.

[0167] An example illustrating property ii) is shown in FIG. 5.

[0168] FIG. 5 shows an illustration of correlation images computed between reference intensity patterns R.sub.L(x,y) obtained with the optical mask 14 at different wavelengths, ranging from 485 nm to 680 nm. Here, by way of illustration, the correlation is evaluated for a signal sub-image extracted from each intensity pattern R.sub.L(x,y). Each sub-image is a macro-pixel whose size is 32×32 pixels. Here, Wiener deconvolution is used as the correlation tool.

[0169] FIG. 5 shows a 6×6 matrix of images. Each image at position (i,j) represents the Wiener deconvolution between the sub-images at wavelength λ.sub.1 and λ.sub.1. As shown, the diagonal of the matrix exhibits bright spots demonstrating the high degree of correlation of an intensity pattern with itself, whereas out of diagonal terms are almost perfectly zero, illustrating the absence of correlation between patterns obtained at different wavelengths.

[0170] The degree of correlation between spectral channels may also be represented in a “correlation matrix” as illustrated in FIG. 6. The amplitude of the coefficient at position (j, k) in the matrix represents the amplitude of the maximum signal in the image located at position (j,k) in FIG. 5. Characterizing the correlation matrix enables to reveal the number of independent spectral channels, 6 in this example. As shown, this matrix is almost diagonal, which demonstrates the feasibility of efficiently uncoupling the contribution of every individual spectral channel once a multi-spectral signal image is recorded.

[0171] The wavefront sensor 10 mentioned above may be used to determine the wavefront shapes of the multi-spectral signal light beam 22 from a single signal image acquisition of the said multi-spectral signal light beam.

[0172] FIG. 7 illustrates a method according to the present invention.

[0173] First, at step 101, several reference intensity patterns R.sub.L(x,y) are recorded. Each reference intensity pattern R.sub.L(x,y) is generated by sending a respective reference incident monochromatic beam L with wavelengths λ.sub.L onto the wavefront sensor 10, L varying from 1 to N with N the number of different reference incident monochromatic beams, x and y are coordinates. Reference incident monochromatic beams may have a same known wavefront shape, for example a planar or a spherical wavefront.

[0174] Then, at step 103, a single signal image I(x,y) of an intensity pattern is recorded using the imaging sensor 18. The latter is generated by sending the said multi-spectral signal light beam 12 onto the optical mask 14. The light beam 12 comprises at least the N wavelengths I.sub.L.

[0175] In a variant, reference intensity patterns R.sub.L(x,y) are recorded after recording the single signal image I(x,y) of the intensity pattern generated by the said multi-spectral signal light beam 12.

[0176] In order to determine the wavefront shape at the N wavelengths λ.sub.L, deformation data T.sub.L.sup.I(x,y) are computed at step 109 using computing means. Such data is representative of a diffeomorphism of each reference intensity pattern R.sub.L(x,y).

[0177] At step 109, intensity-weight data W.sub.L.sup.I(x,y) is also computed. This data is representative of an intensity modulation of each reference intensity pattern at wavelength λ.sub.L.

[0178] All the N intensity-weight data W.sub.L.sup.I(x,y) and the N deformation data T.sub.L.sup.I(x,y) being computed, for L varying from 1 to N, so as to minimize, for all sampling points (x,y) of the surface area A, from the single signal image I(x,y):

[0179] a difference D.sub.A between the single signal image I(x,y) on the one hand, and the sum of reference intensity patterns R.sub.L multiplied by intensity-weight data W.sub.L.sup.I(x,y) and deformed by deformation data T.sub.L.sup.I(x,y), on the other hand:

[00010]$D_A={.Math.I\left(x,y\right)-\underset{L}{.Math.}{W}_L^I\left(x,y\right)R_L\left[\left(x,y\right)+T_L^I\left(x,y\right)\right].Math.}_A$

[0180] the symbol ∥.∥.sub.A designating a norm calculated for all (x,y) sampling points in the surface area A;

[0181] The quantity may comprise a regularization term.

[0182] Finally, at step 111 data are generated for each wavelength λ.sub.L representative of: the shape of the wavefront by integrating the deformation data T.sub.L.sup.I(x,y) over at least one direction of the intensity pattern image, preferably over the two direction of the intensity pattern image; and the intensity map based on the intensity-weight data W.sub.L.sup.I (x,y).

[0183] In some applications, for instance for complex distorted wavefronts, it is preferable to work at a local scale rather than to perform the estimation on the global intensity pattern. One possibility consists in splitting the intensity pattern into signal sub-images, each representative of light impinging on a portion A′ and to estimate local deformation data T.sub.L.sup.I(x,y) for each sub-image and for every wavelength λ.sub.L.

[0184] In this case, the method may comprise two additional steps 105 and 107 preceding step 109. At step 105, the intensity pattern is split into several sub-images, each representative of light impinging on a portion A′ of the at least one surface area A. At step 107, the intensity-weight data W.sub.L.sup.I(x,y) and the deformation data T.sub.L.sup.I(x,y) are calculated between the signal sub-image and each of the reference so as to minimize the differences D.sub.A′ for all the surface areas; and at step 109, deformation data T.sub.L.sup.I(x,y) are obtained by updating the value of T.sub.L.sup.I(x,y) at at least one point (x,y) inside the said surface area A′.

[0185] The intensity-weight data W.sub.L.sup.I(x,y) and the deformation data T.sub.L.sup.I(x,y) may be estimated by computing zero-mean cross-correlation product images, for each signal sub-image, between the signal sub-image and each of the reference intensity patterns, the plurality of reference incident monochromatic beams L being uncorrelated relatively to the zero-mean cross-correlation product, the zero-mean cross-correlation product image between the signal sub-image and each of the reference intensity patterns R.sub.L(x,y) having a peak, the intensity-weight data being the amplitude of the peak and the deformation data T.sub.L.sup.I(x,y) being the displacement vector between the said peak or its centroid from the center of the zero-mean cross-correlation product image.

[0186] In other embodiments, the intensity-weight data W.sub.L.sup.I(x,y) and displacement data T.sub.L.sup.I(x,y) are computed thanks to a Wiener deconvolution of the signal sub-images by the reference intensity patterns R.sub.L(x,y), or a Wiener deconvolution of the reference intensity patterns R.sub.L(x,y) by the signal sub-images.

[0187] The intensity-weight data W.sub.L.sup.I(x,y) and displacement data T.sub.L.sup.I(x,y) may be computed thanks to a matrix inversion algorithm, where the matrices to be inverted being related to at least a transform of a sub-image of a reference intensity pattern R.sub.L(x,y), such as a Fourier transform.

[0188] The intensity-weight data W.sub.L.sup.I(x,y) and displacement data T.sub.L.sup.I(x,y) may be computed thanks to an iterative optimization procedure, for example a mean-squared-difference minimization algorithm as a steepest descent optimization algorithm.

[0189] Intensity-weight data W.sub.L.sup.I(x,y) and displacement data T.sub.L.sup.I(x,y) may be computed thanks to a compressed sensing algorithm.

[0190] Other reconstruction algorithms may be used among which all techniques relying on matrix inversion, for example, Moore-Penrose pseudo-inversion, singular value decomposition, Tikhonov regularization.

[0191] The intensity-weight data W.sub.L.sup.I(x,y) and displacement data T.sub.L.sup.I(x,y) may be computed thanks to a stochastic optimization method, for example, a genetic algorithm or a Markov chains Monte-Carlo (MCMC) method.

[0192] FIG. 8 displays an example of a reconstruction of a single line beam at 485 nm illuminated the imaging sensor 18. In the example illustrated herein, the light beam is modulated by a cylindrical lens. The imaging sensor 18 is a matrix imaging sensor whose size is 1024×768 pixels.

[0193] The full signal intensity pattern is split in macro-pixels of size 64×64 pixels, so resulting in images of size 15×12 macro-pixels. For each macro-pixel the signal intensity pattern is Wiener deconvolved by each reference intensity pattern R.sub.L(x,y) at wavelengths 485 nm, 530 nm, 570 nm, 610 nm, 645 nm, and 680 nm.

[0194] The first line represents a map of intensity weight data W.sub.L.sup.I(x,y). Lines 2 and 3 show displacement data T.sub.L.sup.I(x,y) over horizontal and vertical direction, respectively. Line 4 display data representative of the wavefront shape. In this example, data displayed are representative of beam phases.

[0195] As can be seen, only in the first column, corresponding to the wavelength 485 nm, a significant signal can be detected, so demonstrating that no signal is present in the other spectral channels. Since the wavefront of the signal beam was modulated by a cylindrical lens, a cylindrical phase delay was measured in this spectral channel. Phase maps computed in other spectral channels look like noise since associated with no intensity signal.

[0196] FIG. 9 illustrates an example of performing simultaneous two-color spectro-wavefront measurements. In this case, a synthetic signal intensity pattern is obtained by summing the experimental signal intensity patterns obtained from two different wavelengths, namely at 530 and 610, each of signal intensity patterns being obtained with two different wavefront modulations. Wiener deconvolution macro-pixel per macro-pixels allows retrieving the wavefronts at both wavelengths from this single signal intensity pattern.

[0197] Like in FIG. 4, only spectral channels with significant enough intensity maps are to be considered.

Additional Example of Invention

[0198] In an embodiment, a method according to the present invention relies on the specific properties of random intensity patterns generated by diffusers or more generally by an optical mask having the optical properties:

[0199] i) to cause the intensity pattern to depend on the wavefront shape, so that a tilt applied to the wavefront shape results in a displacement amount of the intensity pattern,

[0200] ii) to produce uncorrelated intensity patterns over at least one surface area A of the imaging sensor, for a plurality of respective incident monochromatic beams of different wavelengths having a same wavefront shape,

[0201] The principle is described below.

[0202] If two beams exhibit the same wavefront and the same intensity distribution, and only differ in their wavelength, the intensity patterns on the camera will be different. The similarity degree between two intensity patterns{s.sub.j}.sub.j=1,2, obtained at two different wavelengths depends on the spectral difference δA=λ.sub.1−λ.sub.2.

[0203] The degree of similarity between two patterns can be mathematically characterized using several correlation tools. One solution to characterize two intensity patterns consists in localizing critical points of the intensity such as intensity maxima and to look for correlations between intensity maxima localization of two patterns. The localization correlation can be quantified using a Pearson correlation coefficient. It can also be quantified using the mean average distance to the nearest neighbor weighted by the critical point density. Such correlation characterization between two different populations of critical points in different field patterns was achieved in Gateau et al. ArXiv 1901:11497.

[0204] Alternatively, correlation between two intensity patterns is characterized by computing the zero-mean cross-correlation product. The statistical average of the cross-correlation product between two uncorrelated random signals (of zero mean) is zero.

[0205] For intensity patterns {R.sub.L}.sub.L=1,2:

[00011]$.Math.\left(R_1-.Math.R_1.Math.\right)*\left(R_2-.Math.R_2.Math.\right).Math.=0$

[0206] Where . stands for a statistical average. This statistical average may be achieved by integrating the intensity pattern over a given surface area of the imaging sensor under ergodic hypothesis.

[0207] Conversely, the zero-mean cross-correlation product of a signal intensity pattern with itself, also called auto-correlation function, is a peaked function:

[00012]$C\left(R_1,R_2\right)=\left(R_1-.Math.R_1.Math.\right)*\left(R_1-.Math.R_1.Math.\right)=\delta \left(r\right)$

[0208] Where δ designates the Dirac distribution, for example if s.sub.1 is assumed to be of variance l, and r is a spatial coordinate on the imaging sensor.

[0209] An optical mask may be designed so that to satisfy orthogonality properties between patterns obtained at different wavelengths.

[0210] Alternatively, a diffuser generating wavelength-dependent random intensity patterns, speckles.

[0211] Mathematically, two uncorrelated random intensity patterns (or speckles) are thus statistically orthogonal relatively to the zero-mean cross-correlation product.

[0212] In cases where the intensity pattern is a speckle, the integration surface area should at least contain one “speckle grain” and the larger the integration surface area, the better the validity of the orthogonality approximation. The zero-mean cross-correlation product thus provides a useful inner product.

[0213] Noteworthy, the orthogonality property satisfied by speckles generated by a diffuser is in particular true in a statistical sense, or in other words, when averaged over a large amount of speckle grains. For a given number of speckle grains, a small spectral difference 8R will not result in any significant change in the intensity pattern. The larger the number of speckle grains, the higher the spectral sensitivity.

[0214] In some embodiments, optical masks with tailored properties may be designed to optimize the orthogonality of patterns obtained at different wavelengths. In this case, regular diffuser would be replaced by specifically engineered diffractive optical elements or pseudo-diffusers. Optimizing the orthogonality of patterns not only reduces the cross-talk between spectral channels but also provides a simple mean to retrieve the wavefronts at each spectral channel thanks to the simple zero-mean cross-correlation product or a Wiener deconvolution.

[0215] Mathematical Grounding, Theoretical Principle

[0216] Given a polychromatic signal beam composed of a discrete spectrum {λ.sub.1, . . . , λ.sub.n}.

[0217] If illuminating the diffuser with such a polychromatic light beam, every individual wavelength will produce its own proper signal intensity pattern. The total signal intensity pattern of the polychromatic beam is then the sum of all the individual signal intensity patterns.

[0218] In a calibration step, the reference intensity patterns s.sub.j({right arrow over (r)}) are obtained at each individual wavelength λ.sub.j, where {right arrow over (r)} is the spatial coordinate vector on the imaging sensor.

[0219] In an embodiment, the wavefront changes are, in particular, beam-tilts, hence resulting in global shifts of the monochromatic signal intensity patterns at the imaging sensor. The total single signal intensity pattern at the imaging sensor 18 will then be:

[00013]$\begin{array}{cc}I\left(\stackrel{.fwdarw.}{r}\right)=\underset{j=1}{\overset{n}{.Math.}}{\alpha }_j{R}_j\left(\stackrel{.fwdarw.}{r}-{\stackrel{.fwdarw.}{r}}_j\right)& \left(1\right)\end{array}$ [0220] where α.sub.j is the weight of the contribution to the signal beam at wavelength λ.sub.j and {right arrow over (r)}.sub.j a translation vector resulting from the beam-tilt.

[0221] Mathematical treatments allow to retrieve every individual α.sub.j and {right arrow over (r)}.sub.j.

[0222] To achieve so, the zero-mean cross-correlation product on R.sub.n may be performed:

[00014]$.Math.\left(I-.Math.I.Math.\right)*\left(R_n-.Math.R_n.Math.\right).Math.=\underset{j=1}{\overset{n}{.Math.}}{\alpha }_j.Math.\left(R_j\left(\stackrel{.fwdarw.}{r}-{\stackrel{.fwdarw.}{r}}_j\right)-.Math.R_j.Math.\right)*\left(R_n-.Math.R_n.Math.\right).Math.$$.Math.\left(I-.Math.I.Math.\right)*\left(R_n-.Math.R_n.Math.\right).Math.={\alpha }_n\delta \left(\stackrel{.fwdarw.}{r}-\stackrel{.fwdarw.}{r_n}\right)$

[0223] Where . stands for a statistical average. The result of this mathematical treatments may be a translation image exhibiting a peaked function of amplitude α.sub.j,whose peak is centered at {right arrow over (r)}.sub.j. The parameter {right arrow over (r)}.sup.j gives access to the uniform deformation data T.sub.L (x,y)=−{right arrow over (r)}.sub.j while α.sub.j gives access to the uniform intensity weight data W.sub.L(x,y)=α.sub.j of the beam at wavelength λ.sub.j.

[0224] In the following, several possible algorithmic implementations are cited for retrieving these parameters. However, whatever the algorithm used, the orthogonality of the wavelength-dependent reference intensity patterns relatively to the zero-mean cross-correlation product is the common theoretical grounding requirement.

[0225] Practical Algorithmic Implementation

[0226] In some embodiments, typical wavefront distortions are more complex than simple tilts. Integration of the displacement vector maps will then provide the wavefronts, as described above.

[0227] Algorithms such as Demon Algorithm may be used in order to compute direcTLIy a distortion map (Berto, P., Rigneault, H., & Guillon, M. (2017), Wavefront sensing with a thin diffuser. Optics letters, 42(24), 5117-5120). Deep-learning computational approaches based on neuronal networks can also be implemented in order to find the distortion maps.

[0228] A more simple approach consists in splitting the original intensity map into macro-pixels. At a “local scale”, the signal intensity pattern is just translated and pattern deformation may be neglected, making algorithms more simple. The “local” translation of the signal intensity pattern is proportional to the “local tilt” of the wavefront. One algorithmic difficulty consists defining what the size of the “local scale” is, what should the size of the macro-pixel be. In practice, it depends on geometrical and physical parameters of the WFS. This point is thus discussed now.

[0229] As a first simple algorithmic implementation, the zero-mean cross-correlation product, extensively discussed above, may be used. At the defined “local scale”, the signal intensity pattern is just a translation of the reference intensity pattern. A crop in the pixilated full field camera image of the signal intensity pattern must thus be performed. The crop size now defines a macro-pixel (composed of several camera pixels), which will correspond, in the end, to the phase pixel; i.e. the pixel of the final phase image. Although the simplest solution consists in considering a square macro-pixel of constant size, the macro-pixel size may vary depending on the camera region. The size of the macro-pixel is a balance between the spatial resolution of the final phase image (the larger the macro-pixel, the fewer they are) and the spectral resolution of the WFS (the larger the macro-pixel, the larger the number of speckle grains per macro-pixel and the larger the number of orthogonal spectral modes over a given spectral width).

[0230] Once the signal intensity pattern is split into macro-pixels, every macro-pixel can be projected on the corresponding macro-pixels of the reference intensity pattern, thanks to the zero-mean cross-correlation product, as described above.

[0231] Alternatively, a Wiener deconvolution can be used rather than a zero-mean cross-correlation product. Wiener deconvolution is very similar to a zero-mean cross-correlation product. This similarity appears when performing these operations in the Fourier domain.

[0232] Considering a signal image I and a reference image R, their two-dimensional Fourier transform are then written Î and {circumflex over (R)}. The Fourier transform of their cross-correlation is then F=*{circumflex over (R)}*Î where {right arrow over (R)}* is the complex conjugate of R. By comparison, the Fourier transform of the Wiener deconvolution of S by R is:

[00015]$W=\frac{R^{*}I}{{\sigma }^2+R^2}.$

Here the term σ.sup.2 is an average signal to noise ratio. Note that both the cross-correlation product and Wiener deconvolution are slightly dissymmetric in R and I. For our data treatment exchanging the role of R and I does not change significantly the result. Our preliminary data presented below are obtained using Wiener deconvolution.

[0233] In practice, intensity patterns obtained at different wavelengths are only orthogonal in a statistical sense. Consequently, there may be some cross-talk between spectral channels. To reduce these cross-talk Iterative Wiener deconvolution can be implemented to suppress the cross talk effect. In this implementation, an estimate of α.sub.j and {right arrow over (r)}.sub.j are first deduced from the translation image obtained thanks to a Wiener deconvolution. Second, assuming the expected signal intensity pattern is numerically rebuilt according to equation (1). As a third step, the expected signal intensity pattern is compared to the actual experimental signal intensity pattern and differences are sent as an input for step 1. Steps 1 to 3 can then be iterated.

[0234] A more elaborated compressed sensing algorithms can also be used, taking into account that a given wavelength can only be responsible for a single peak. Such algorithms are optimized to simultaneously minimize two quantities, one of which being the root mean squared error between the experimental data (the signal intensity beam) and the rebuilt data. The other quantity to be minimized here is the number of non-zero coefficients in each translation image.

[0235] More elaborated reconstruction algorithms may be used among which all techniques relying on matrix inversion: Moore-Penrose pseudo-inversion, singular value decomposition (also called principal component analysis), Tikhonov regularization etc. For instance, principal component analysis was used in N. K. Metzger et al., Nat Commun 8:15610 (2017) to make a very sensitive spectrometer. Such a matrix pseudo-inversion can be achieved the following way. In the Fourier domain, equation (1) can be written:

[00016]$F\left(I\right)\left(\stackrel{.fwdarw.}{k}\right)=\underset{j=1}{\overset{n}{.Math.}}{\alpha }_j{R}_j\left(\stackrel{.fwdarw.}{k}\right)e^{-i\stackrel{.fwdarw.}{k}.Math.{\stackrel{.fwdarw.}{r}}_j}$

[0236] Where R.sub.n({right arrow over (k)}) are the 2D Fourier transforms of reference intensity patterns, we then have:

[00017]$R_n^{*}\left(\stackrel{.fwdarw.}{k}\right)F\left(I\right)\left(\stackrel{.fwdarw.}{k}\right)=\underset{j=1}{\overset{n}{.Math.}}{\alpha }_j{R}_n^{*}\left(\stackrel{.fwdarw.}{k}\right)R_j\left(\stackrel{.fwdarw.}{k}\right)e^{-i\stackrel{.fwdarw.}{k}.Math.{\stackrel{.fwdarw.}{r}}_j}$

[0237] Which can be simply re-written as a matrix equation to invert:

[00018]$V_I\left(\stackrel{.fwdarw.}{k}\right)=M\left(\stackrel{.fwdarw.}{k}\right)V_O\left(\stackrel{.fwdarw.}{k}\right)$$\mathrm{Where}{V}_I=R_n^{*}F\left(I\right),M_{n,j}=R_n^{*}{R}_j\mathrm{and}{V}_O\left(\stackrel{.fwdarw.}{k}\right)={\alpha }_j{e}^{-i\stackrel{.fwdarw.}{k}.Math.{\stackrel{.fwdarw.}{r}}_j}.$

[0238] The present invention is not limited to these embodiments, and various variations and modifications may be made without departing from its scope.