sensors Article

Dynamic Aberration Correction for Conformal Window of High-Speed Aircraft Using Optimized Model-Based Wavefront Sensorless Adaptive Optics Bing Dong 1, *, Yan Li 2 , Xin-li Han 1 and Bin Hu 2 1 2

*

School of Optoelectronics, Beijing Institute of Technology, Beijing 100081, China; [email protected] Beijing Institute of Space Mechanics & Electricity, Beijing 100090, China; [email protected] (Y.L.); [email protected] (B.H.) Correspondence: [email protected]; Tel.: +86-10-6891-3794

Academic Editor: Vittorio M. N. Passaro Received: 21 June 2016; Accepted: 22 August 2016; Published: 2 September 2016

Abstract: For high-speed aircraft, a conformal window is used to optimize the aerodynamic performance. However, the local shape of the conformal window leads to large amounts of dynamic aberrations varying with look angle. In this paper, deformable mirror (DM) and model-based wavefront sensorless adaptive optics (WSLAO) are used for dynamic aberration correction of an infrared remote sensor equipped with a conformal window and scanning mirror. In model-based WSLAO, aberration is captured using Lukosz mode, and we use the low spatial frequency content of the image spectral density as the metric function. Simulations show that aberrations induced by the conformal window are dominated by some low-order Lukosz modes. To optimize the dynamic correction, we can only correct dominant Lukosz modes and the image size can be minimized to reduce the time required to compute the metric function. In our experiment, a 37-channel DM is used to mimic the dynamic aberration of conformal window with scanning rate of 10 degrees per second. A 52-channel DM is used for correction. For a 128 × 128 image, the mean value of image sharpness during dynamic correction is 1.436 × 10−5 in optimized correction and is 1.427 × 10−5 in un-optimized correction. We also demonstrated that model-based WSLAO can achieve convergence two times faster than traditional stochastic parallel gradient descent (SPGD) method. Keywords: conformal optics; dynamic aberration; adaptive optics; wavefront sensorless; deformable mirror

1. Introduction Currently, targeting and sensing on airborne platforms is mostly realized by imaging through a spherical or flat window. Although easy to design and fabricate, these windows show poor aerodynamic performance and can substantially increase air drag and radar cross section (RCS). For aircraft with high Mach number, an aspheric conformal window that follows the general contour of surrounding surface is preferred. However, the local shape of conformal windows changes across the field-of-regard (FOR), leading to large amounts of FOR-dependent wavefront aberrations and degraded images [1]. A variety of aberration correction methods using fixed or dynamic correctors have been proposed for conformal optics in recent years [2–4]. We are particularly interested in dynamic correction approach using a deformable mirror (DM) which can be taken as a freeform element in optical design [5]. Li et al. designed a missile seeker with an ellipsoidal conformal dome and a DM for correction [6,7]. The DM can be driven in open-loop or closed-loop mode. In open-loop mode, DM is shaped in response to parameters like look angle, window temperature and other flight conditions. The commands of DM for a given set of parameters are predetermined and stored in

Sensors 2016, 16, 1414; doi:10.3390/s16091414

www.mdpi.com/journal/sensors

Sensors 2016, 16, 1414

2 of 13

a lookup table, so a lot of simulations and experiments must be done to calibrate the control system before use. In closed-loop mode, DM is integrated into an adaptive optics (AO) control system [8]. The control signal of DM is derived from a dedicated wavefront sensor (i.e., classical AO) or from the final image detector (i.e., wavefront sensorless adaptive optics—WSLAO). Classical AO uses a wavefront sensor to measure aberration and a DM to make conjugated correction. The correction bandwidth of classical AO is usually higher than WSLAO. However, typical wavefront sensors like Hartmann-Shack may be inapplicable because of insufficient photons, non-common path error or structural constraints. A dedicated wavefront sensor is omitted in WSLAO and the control signal of DM is derived from the intensity information on the image detector. Employing WSLAO in an airborne imaging system can make most of the imaging photons and benefit from its simple structure. The control algorithms of WSLAO can be divided into two categories: model-free algorithms and model-based algorithms. “Model-free” means the exact relationship between the metric function and aberration is unknown. Several model-free algorithms like stochastic parallel gradient descent (SPGD) have been proposed and demonstrated [9,10]. The common disadvantage of model-free algorithms is that a large number of iterations are usually needed and a global convergence is not always guaranteed. Several hundreds of iterations are typically needed to achieve convergence with SPGD, which severely limits its temporal bandwidth and dynamic correction performance. In model-based WSLAO, wavefront aberration is decomposed into specific modes like Zernike [11], Lukosz [12] or DM modes [13]. The choice of mode is based on the model used to describe the mathematical relationship between modal coefficients and metric function. The modal coefficients can be calculated directly from that relationship. The metric function can be derived from Strehl ratio (SR), image spectral density or mean-square spot radius [11–13]. Model-based WSLAO can achieve convergence much faster than model-free algorithms and avoid dropping into local optima, which are both crucial for dynamic aberration correction. The model-based WSLAO must be calibrated before correction in order for the DM to generate a specific mode, which is not required in model-free algorithm like SPGD. The influence functions of DM actuators should be measured by a wavefront sensor or interferometer in a calibration process. In this paper, model-based WSLAO is applied for correction of dynamic aberration induced by conformal window of an infrared remote sensor on high-speed aircraft. Two approaches are proposed to accelerate the correction. Numerical simulations and experiments are made to demonstrate the feasibility and performance of our method. 2. Optical System The infrared remote sensing system coupled with a conformal window and a scanning mirror is depicted in Figure 1. The optical system is modeled and optimized in ZEMAX. The secondary mirror of the two-mirror reflective telescope is replaced by a DM that can be shaped to compensate the FOR-dependent aberrations. The system’s FOR is ±30◦ and the field of view (FOV) at any instant is ±0.25◦ . The conformal window is toroidal type that is formed by rotating a curve in the Y-Z plane about an axis parallel to Y axis. The working wavelength is in long-wave infrared from 7.7 to 10.3 µm. The system parameters are summarized in Table 1. Table 1. System Parameters. Parameter

Value

Parameter

Value

Window surface type Window material Field of view (◦ ) Primary mirror diameter (mm) Working f-number

Toroid Germanium ±0.25 200 2.8

Window thickness (mm) Window diameter (mm) Field of regard (◦ ) DM diameter (mm) Working wavelength (µm)

20 500 ±30 48 7.7–10.3

Sensors 2016, 16, 1414 Sensors Sensors 2016, 2016, 16, 16, 1414 1414

3 of 13 33 of of 13 13

Figure 1. Optical system at 00°◦ and 30° 30◦ look angle. Figure 1. Optical system at 0° and 30° look angle.

RMS Wavefront Error(rad) RMS Wavefront Error(rad)

The conformal aberrations that vary with looklook angle. The wavefront rootThe conformal window windowgives givesrise riseto aberrations that vary with angle. The wavefront The conformal window gives rise toto aberrations that vary with look angle. The wavefront rootmean-square (RMS) value changing with look angle is shown in Figure 2. The window is symmetrical root-mean-square (RMS) value changing with look angle isinshown in The Figure 2. The window is mean-square (RMS) value changing with look angle is shown Figure 2. window is symmetrical to the center to ofthe FOR, so we only so show positive degrees. we use an ideal to make ideal symmetrical center of FOR, we only show positiveIf degrees. If we useDM an ideal DM an to make to the center of FOR, so we only show positive degrees. If we use an ideal DM to make an ideal wavefront correction at a wavelength of 9 μm, errorserrors at other wavelengths are are shown in an ideal wavefront correction at a wavelength of 9residual µm, residual at other wavelengths shown wavefront correction at a wavelength of 9 μm, residual errors at other wavelengths are shown in Figure 3. The dotted line denoting diffraction limit corresponds to RMS error of 0.45 rad or Strehl in Figure 3. The dotted line denotingdiffraction diffractionlimit limitcorresponds correspondstotoRMS RMSerror errorof of 0.45 0.45 rad rad or or Strehl Figure 3. The dotted line denoting Strehl Ratio of 0.8. 0.8. Wavefront than this criterion has has little little impact impact on on image image quality. quality. So the Ratio of Wavefront error error lower lower than this criterion So the Ratio of 0.8. Wavefront error lower than this criterion has little impact on image quality. So the system’s chromatic aberration aberration which which cannot cannot be be compensated by DM is negligible. system’s chromatic compensated by DM is negligible. In In our our simulations, simulations, system’s chromatic aberration which cannot be compensated by DM is negligible. In our simulations, the wavefront aberrations aberrations are are estimated estimated and and corrected corrected at at 99 µm. μm. the the wavefront wavefront aberrations are estimated and corrected at 9 μm. be insufficient in correction of large A critical parameter parameter of DM is actuator stroke which might A critical of DM is actuator stroke which might be in correction of A critical parameter of DM is actuator stroke which might be insufficient insufficient in correction of12large large aberrations. As our system is working at infrared, the stroke requirement can be as large as μm aberrations. As our system is working at infrared, the stroke requirement can be as large as 12 µm aberrations. As our system is working at infrared, the stroke requirement can be as large as 12 μm which is estimated from the peak-to-valley value of aberrations. Thanks to the development of DM which is estimated from the peak-to-valley value of aberrations. Thanks to the development of DM which is estimated from DMs the peak-to-valley value ofand aberrations. Thanks toavailable the development of DM technology, large-stroke have been developed are commercially now [14–16]. technology, technology, large-stroke large-strokeDMs DMshave havebeen beendeveloped developedand andare arecommercially commerciallyavailable availablenow now[14–16]. [14–16].

Figure 2. 2. RMS RMS of of wavefront wavefront aberration aberration as as aa function function of of look Figure look angle angle (at (at 99 µm). μm). Figure 2. RMS of wavefront aberration as a function of look angle (at 9 μm).

Sensors 2016, 16, 1414 Sensors 2016, 16, 1414

4 of 13

4 of 13

Figure 10.3 μm µm after aftercorrection correctionatat9 9μm. µm. Figure3.3.Residual Residualerror errorat at 7.7 7.7 µm μm and 10.3

PrincipleofofModel-Based Model-Based WSLAO WSLAO 3. 3. Principle extended objects,low lowspatial spatialfrequency frequency content content of (m) can be ForFor extended objects, of the the image imagespectral spectraldensity densitySJS J (m) can be used as a metric function in model-based WSLAO: used as a metric function in model-based WSLAO: 1 Z2 MZ2 M  −1 G   2π SJ (2m)mdmd   G =  0 M1 S J (m)mdmdξ 

0

(1)

M1

(1)

where M1 and M2 are normalized spatial frequency, ξ is the angle of spatial frequency, and where M2 are is normalized frequency, ξ is relationship the angle of between spatial frequency, m = (mcosξ, 1 and msinξ) m =M (mcosξ, the spatial spatial frequency vector. The wavefrontand aberration Φ msinξ) is the spatial frequency vector.by: The relationship between wavefront aberration Φ and the metric and the metric function G is given function G is given by: 1x 2 2 G q q+ q 1 Φ ds (2) (2) G≈ 1 1 q22   pp |∇ Φ | ds π s where q1 qand q2qare constants structure. ∬p ds denotes overthe thepupil pupil 1 and 2 are constantsdepending depending on on object object structure. denotes an an integral integral over where area p. Wavefront aberration Φ can be represented by a series of Lukosz modes L whose derivatives i derivatives area p. Wavefront aberration Φ can be represented by a series of Lukosz modes Li whose areare orthogonal with each other: orthogonal with each other: N

x p

N a L Φ=∑ i i Φ  i =4 ai Li i 4 (

1 i=j ∇ Li · ∇ L j ds =1 i  j  L  L ds   p i j 0 0 i i 6=j j

(3)

(4)

(3)

(4)

The expression and aberration definition of Lukosz modes can be found in [12]. Using Equations The expression and aberration definition of Lukosz modes can be found in [12]. Using (3) and (4), the metric function can be rewritten in terms of Lukosz coefficients: Equations (3) and (4), the metric function can be rewritten in terms of Lukosz coefficients: NN

 q22 G G=qq11+ ∑aai 2i 2

(5) (5)

ii =44

After introducing iLi by DM, The initial metric functionisisassumed assumedas asGG0.. After The initial metric function introducingaapositive positivemodal modalbias bias+b+b 0 i Li by DM, the metric function becomes G+. Applying a negative modal bias –biLi corresponds to metric function the metric function becomes G+ . Applying a negative modal bias –bi Li corresponds to metric function of G–. Then we can get an equation set as: of G− . Then we can get an equation set as:

Sensors 2016, 16, 1414 Sensors 2016, 16, 1414

5 of 13 5 of 13

   qq22aai 2i 2 G0 q1 q+1 q2q2∑aakk22 + G0=    k 6=k i i    2 bi )b2i )2 G+G= q1q+ aakk 2 +qq2 2( a( ai i+  1 qq 22 ∑ k ii  k6=   2  G   q q bi )b2 )2   1 2  G = q + q aakk 2 +qq2 2( a( ai i−  − 2∑ 1 i k i 

(6)(6)

k 6 =i

required correction amount each modecan canbebeestimated estimatedby bysolving solving Equation Equation (6): (6): TheThe required correction amount forfor each mode b (G  G ) ai,corr   ai   bii( G+ − G− ) 2 G ai,corr = − ai = −   4G0  2G

(7) (7) 2G+ − 4G0 + 2G− From Equation (7), to estimate one modal coefficient, at least three metric function From Equation (7), to estimate one modal coefficient, at least three metric function measurements measurements are required. To correct N modes, two strategies might be used. One can apply arecorrection required. for Toall correct modes, two strategies be used. One can apply for+ all modesNsimultaneously after theirmight coefficients are estimated, whichcorrection requires 2N 1 modes simultaneously after their coefficients are estimated, which requires 2N + 1 images in total. images in total. Alternatively, the correction can be applied in turn for each mode immediately after Alternatively, the can be applied in turn each in mode after its coefficient its coefficient is correction estimated. This strategy requires 3Nfor images totalimmediately as the metric function with zero is estimated. This strategy requires 3N images in total as the metric function with zero ismodal modal bias (i.e., G0) must be recalculated for each mode. The flow chart of 3N algorithm shownbias is (i.e., G0 ) must be recalculated each mode. correction The flow chart of 3N algorithm shown is Figurein4. Figure 4. We will compareforthe dynamic performance of these istwo algorithms Wenext will sections. compare the dynamic correction performance of these two algorithms in next sections.

Figure Algorithmflow flowchart chartof of3N 3Nalgorithm. algorithm. Figure 4.4.Algorithm

Sensors 2016, 16, 1414

6 of 13

Sensors 2016, 16, 1414

6 of 13

4. Algorithm Optimization and Simulation

4.1.Algorithm SimulationOptimization Method 4. and Simulation The model-based algorithm is programmed in the MATLAB environment. The wavefront 4.1. Simulation Method aberration induced by conformal window shape at any look angle is obtained from the optical design The model-based algorithm is programmed in the MATLAB environment. wavefront software (ZEMAX). These aberrations are extracted into MATLAB through its The dynamic data aberration induced by conformal window shape at any look is obtained from the exchange (DDE) interface with ZEMAX. The compensated DMangle deformation calculated by optical modeldesign software (ZEMAX). aberrations are extracted into MATLAB its dynamic data based algorithm is sent backThese to ZEMAX to evaluate the image quality afterthrough correction. exchange (DDE) that interface with ZEMAX. Assuming sampling width ofThe thecompensated object is L DM and deformation diameter ofcalculated the pupilby ismodel-based D, then the algorithm is sent back to ZEMAX to evaluate the after sampling interval of spatial frequency is equal toimage 1/L. Toquality improve thecorrection. sampling rate and prevent the sampling width ofadopted the object L and diameter of thethe pupil is D,Here then γ the aliasAssuming effect, L = that γD (γ > 1) is usually byispadding zero outside pupil. is sampling the zerointerval spatial frequency is equaltotothe 1/L. To improve the sampling rate and prevent the alias paddingoffactor which corresponds number of samples along radial direction within λ/D effect, at the Limage = γDplane. (γ > 1)The is usually adopted by padding zero outside the pupil. Here γ is the zero-padding factor wavefront aberration for any field of view obtained from ZEMAX is sampled by a which corresponds to the number of along radial direction within λ/D atobject the image plane. 128 × 128 grid and then embedded to samples a 256 × 256 matrix (i.e., γ = 2). The extended is chosen as The wavefront aberration any field view obtained sampled a 128 × 128 grid a typical remote sensing for image with of256 × 256 pixels.from TheZEMAX detectedis image is by generated by twoand then embedded to convolution a 256 × 256 matrix (i.e.,and γ = 2). extended is chosen as a typical remote dimensional discrete of object theThe point spreadobject function. As proven in [12], the sensing image with 256 256 pixels. The detected image is generated by two-dimensional discrete wavefront estimation of × our algorithm is independent of object structure. convolution of object and the 52 point spread (red function. Asisproven in [12], the wavefront estimation of our The layout of DM with actuators circles) depicted in Figure 5. All influence functions algorithm is independent of object not structure. are assumed identical. However, all actuators are fully illuminated. The dashed inner-circle in The layout ofthe DM with 52 actuators (red is depicted 5. All influence functions Figure 5 denotes clear aperture of DM. Thecircles) influence functioninisFigure assumed to be Gaussian-type as are assumed identical. However, not all actuators are fully illuminated. The dashed inner-circle in given in Equation (8). The coupling coefficient of DM with continuous faceplate is usually between Figure 5 denotes thethe clear aperture of DM. The is assumed be Gaussian-type 10% and 20%. Here coupling coefficient c isinfluence set as 0.2.function The actuator stroke to is larger than 12 μm as given in Equation (8).2:The coupling coefficient of DM with continuous faceplate is usually between estimated in Section 10% and 20%. Here the coupling coefficient c is set as 0.2. The actuator stroke is larger than 12 µm as 2 2   estimated in Section 2: " ln c ( x  xi )  ( y  yi ) # Ii (x, y)  exp (8)  ( x − x )2 +2 (y − y )2 i d i  Ii ( x, y) = exp lnc (8) d2 Here  xi , yi  denotes actuator position and d is actuator pitch. Here( xi , yi ) denotes actuator position and d is actuator pitch.

Figure 5. Layout of DM with 52 actuators. Figure 5. Layout of DM with 52 actuators.

4.2. Optimized Dynamic Correction 4.2. Optimized Dynamic Correction For dynamic aberration correction, the temporal bandwidth of the AO system must be higher For dynamic aberration correction, the temporal bandwidth of the AO system must be higher than that of aberration. Besides using more powerful hardware (i.e., image detector and DM with than that of aberration. Besides using more powerful hardware (i.e., image detector and DM with high frame rate) [17], the wavefront control algorithm must also be optimized to accelerate the high frame rate) [17], the wavefront control algorithm must also be optimized to accelerate the

Sensors Sensors2016, 2016,16, 16,1414 1414 Sensors 2016, 16, 1414

13 77ofof13 7 of 13

Mode coefficient Mode coefficient

convergence. Two approaches are investigated here to improve the efficiency of the model-based convergence. Two Two approaches are investigated investigated here improve the of model-based convergence. here to toused improve the efficiency efficiency of the the model-based algorithm. The firstapproaches one is to reduce the mode number in correction. As shown in Equation (7), algorithm.The Thefirst first one one is is to to reduce reduce the mode number used in correction. As shown ininEquation (7), algorithm. used in correction. As shown Equation (7), the estimation of modal coefficient is independent of each other. If we have a priori knowledge of the the estimation of modal coefficient is independent of each other. If we have a priori knowledge of the the estimation of modal coefficient ison independent of each other. If we have a priori knowledge of the aberration, making corrections only those dominant modes with large amplitude can significantly aberration,making makingcorrections corrections only only on on those those dominant modes with large amplitude can significantly aberration, dominant modes with large amplitude can significantly improve the correction speed. improvethe thecorrection correctionspeed. speed. improve The wavefront aberrations induced by conformal window and subsequent imaging system can The wavefront aberrations induced by by conformal conformal window subsequent imaging system can The wavefront aberrations window and and subsequent imaging can be obtained from ZEMAX. Theninduced the aberrations are decomposed into Lukosz modes. Thesystem dominant be obtained from ZEMAX. Then the aberrations are decomposed into Lukosz modes. The dominant be obtained from ZEMAX. Then the aberrations are decomposed into Lukosz modes. The dominant Lukosz modes of our system are 4th (defocus), 6th (astigmatism), 7th (coma), 9th (trefoil) and their Lukoszmodes modesof ofour oursystem system are are 4th 4th (defocus), (defocus), 6th 6th (astigmatism), (astigmatism), 7th 9th (trefoil) and Lukosz 7th (coma), (coma), (trefoil) andtheir their coefficients varied across the FOR are shown in Figure 6. The coefficients of other9th modes are less than coefficients varied across the FOR are shown in Figure 6. The coefficients of other modes are less coefficients varied across the FORhere are is shown in Figure coefficients of other modes lessthan than 0.1. The modal coefficient used normalized by 6. itsThe RMS value, which means a unitare amplitude 0.1.The Themodal modalcoefficient coefficientused usedhere hereisisnormalized normalizedby byitsitsRMS RMSvalue, value,which whichmeans meansa aunit unitamplitude amplitudeof 0.1. of mode corresponding to a phase aberration with RMS of 1 rad. The correction results of four of mode corresponding to a phase aberration withofRMS of The 1 rad. The correctionofresults of four mode corresponding to a phase aberration with RMS rad. correction four dominant dominant Lukosz modes and first 20 Lukosz modes are1both shown in Figureresults 7 for comparison. The dominant Lukosz modes and first 20 Lukosz modes are both shown in Figure 7 for comparison. The Lukosz modes andcorrection first 20 Lukosz modes is are both shown in Figure for comparison. The Strehl ratio Strehl ratio after of 20 modes just slightly better than 7correction of 4 dominant modes. Strehl ratio after correction of 20 modes is just slightly better than correction of 4 dominant modes. after correction of 20 modes is just slightly better than correction of 4 dominant modes. However, However, 32 more images have to be captured and then evaluated when using “2N + 1” algorithm. However, 32 more images have to be captured and then evaluated when using “2N + 1” algorithm. 2) [18]. Here Strehl 32 more images have to be captured andRMS then(σ) evaluated whenaberration using “2Nby + 1” algorithm. Here Strehl ratio is calculated from the of wavefront exp(–σ Here Strehl ratio is calculated from the RMS (σ) of wavefront aberration by exp(–σ2) [18]. 2 ratio is calculated from the RMS (σ) of wavefront aberration by exp(–σ ) [18].

Strehl Ratio Strehl Ratio

Figure 6. Coefficients Coefficients of dominant dominant Lukosz Lukosz modes varying across FOR. Figure modes varying varying across acrossFOR. FOR. Figure 6. 6. Coefficients of of dominant Lukosz modes

Figure or four four dominant dominantmodes. modes. Figure7.7.Strehl Strehlratio ratio variation variation curve curve in correction of first 20 modes or Figure 7. Strehl ratio variation curve in correction of first 20 modes or four dominant modes.

Sensors 2016, 16, 1414 Sensors 2016, 16, 1414

8 of 13 8 of 13

Strehl Ratio

In only arise arise from from the the asymmetry asymmetry of of the the conformal conformal Inpractice, practice, wavefront wavefront aberrations aberrations may not only window but also be induced by aero-optical effects. The advantage of using DM on-board is the that window but also induced by aero-optical effects. The advantage of using DM on-board is that the mirror shape can adjustedaccordingly accordinglytotoaccommodate accommodatevarious variousdisturbances. disturbances.The Theaberrations aberrations mirror shape can bebe adjusted caused complicated than than conformal conformal window window caused by by aero-optical aero-optical effects effects are are random random and much more complicated shape. aberrations depends dependson onmany manyfactors factorssuch suchasas shape.The Themagnitude magnitude and and distribution distribution of aero-optical aberrations flightspeed speedand andattitude, attitude, thermal thermal environment, environment, window shape flight shape and and material, material,etc. etc.[19]. [19].However, However,ititisis stillpossible possibleto topredict predictthe the aero-optical aero-optical aberrations aberrations by simulation simulation and still and experiment experimentand andthen thenextract extract dominant aberration modes, which is our future work [20,21]. Alternatively, low order DM modes dominant aberration modes, which is our future work [20,21]. Alternatively, low order DM modes can can be selected for correction we have no prior knowledge of aberration the aberration be selected for correction if weifhave no prior knowledge of the [13].[13]. Thesecond secondapproach approach to to accelerate accelerate the correction is to minimize The minimize the the image image size size while whilekeeping keeping satisfactory corrected accuracy. The most time-consuming part of our algorithm is the image satisfactory corrected accuracy. The most time-consuming part of our algorithm is the image acquisition acquisition and subsequent of metric that function that involves two-dimensional Fast process and process subsequent calculationcalculation of metric function involves two-dimensional Fast Fourier Fourier Transform (FFT) of captured image. The computing complexity (Mwhere log M)M where Transform (FFT) of captured image. The computing complexity of FFT isofOFFT (M is logOM) is the M is the image size,aso using aimage smaller a lot of computational thusthe improve the image size, so using smaller canimage save acan lot save of computational time thus time improve correction correction speed. Fortargets, point-like low order only affect the intensity nearofthe center of speed. For point-like lowtargets, order modes onlymodes affect the intensity near the center image spots. image spots. In this paper, we use extended targets and the metric function is low spatial frequency In this paper, we use extended targets and the metric function is low spatial frequency content of the content of the image spectral density. Adopting a smaller leads to a larger sampling interval image spectral density. Adopting a smaller image leads toimage a larger sampling interval in the spatial in the spatial frequency domain and reduces the correction accuracy. Simulations can be made to frequency domain and reduces the correction accuracy. Simulations can be made to determine the determinevalue the minimum minimum of image value size. of image size. For 100 samples of random aberrations composed For 100 samples of random aberrations composed of of first first 15 15Lukosz Lukoszmodes modeswith withRMS RMSofof11rad, rad, the mean Strehl ratio after correction is calculated for different image size in Figure 8. The correction the mean Strehl ratio after correction is calculated for different image size in Figure 8. The correction results of of using using 2N 2N ++ 11 images images and and 3N 3N images images (here results (here mode mode number number N N == 15) 15) are are both both given givenfor for comparison. To obtain Strehl ratio better than 0.8, images containing at least 50 × 50 pixels are comparison. To obtain Strehl ratio better than 0.8, images containing at least 50 × 50 pixels are required required for “3N algorithm” and at least 150 × 150 pixels for “2N + 1 algorithm”. From Figure 8, “3N for “3N algorithm” and at least 150 × 150 pixels for “2N + 1 algorithm”. From Figure 8, “3N algorithm” algorithm” is superior to “2N + 1 algorithms” in correction accuracy for same image size. In “3N is superior to “2N + 1 algorithms” in correction accuracy for same image size. In “3N algorithm”, algorithm”, each modal aberration is corrected by DM immediately after the estimation of each modal each modal aberration is corrected by DM immediately after the estimation of each modal coefficient, coefficient, which is more appropriate for fast-changing dynamic aberration and can alleviate the which is more appropriate for fast-changing dynamic aberration and can alleviate the crosstalk problem crosstalk problem in coefficient estimation. in coefficient estimation.

Figure with image image size. size. Figure 8. 8. Strehl ratio variation with

4.3.Dynamic DynamicCorrection CorrectionSimulation Simulation 4.3. Inour ourapplication, application, speed at which aberrations change depends the mirror scan rate, In thethe speed at which aberrations change depends on theon mirror scan rate, without without consideration of aero-optical effects. In practice, the scan rate depends on many factors like consideration of aero-optical effects. In practice, the scan rate depends on many factors like the aircraft’s the aircraft’s speed and height, rate of radiation the detector, radiation intensity of observed targets andof speed and height, frame rate of frame the detector, intensity of observed targets and the speed the speed of the data processing unit. The scan rate may vary widely during the flight mission, from

Sensors 2016, 16, 1414 Sensors 2016, 16, 1414

9 of 13 9 of 13

Strehl Ratio

the datatoprocessing unit.per Thesecond. scan rate may during the flight mission, from several to several tens of degrees What wevary care widely about most is the AO system’s performance under tens of degrees per second. weconstantly care aboutchanging most is the AO system’s under different different scan rates. BecauseWhat of the aberration and performance inevitable time delay of AO scan rates. Because of create the constantly changingwavefront aberrationfor and inevitable timeangle delaybased of AO systems, DM can only a complementary the current look onsystems, images DM can only create a look complementary current look angle on images acquired at previous angles. Thewavefront time delayfor inthe each correction cycle,based typically on the acquired order of at previous look angles. The time delay in each correction cycle, typically on hardware the order of milliseconds milliseconds level, depends on the specific correction algorithm and AO configuration. level, depends on the specific correction algorithm and aberration AO hardware configuration. However, partial However, partial correction is still possible if the changing is continuous with no correction is still possible if the aberration changing is continuous with no abrupt disturbance. abrupt disturbance. wewe dodo notnot need to know thethe exact timetime delay value and To predict predictthe thedynamic dynamicperformance performanceofofAO, AO, need to know exact delay value the exact scanning rate. rate. We only needneed to know how much the aberration changed after the time and the exact scanning We only to know how much the aberration changed after thedelay. time The time estimation and correction of a single (i.e., one (i.e., correction cycle) is assumed delay. Thedelay timefor delay for estimation and correction of amode single mode one correction cycle) is as t seconds includes of image metric function DM response. assumed as twhich seconds which time includes time acquisition, of image acquisition, metric calculation, function calculation, DM After t seconds, the change look angle assumed as ∆θ. as If Δθ. ∆θ is small, the correction is still response. After t seconds, thein change in lookisangle is assumed If Δθ is small, the correction is effective because the aberration changing is usually continuous without abrupt disturbance. Based on still effective because the aberration changing is usually continuous without abrupt disturbance. the DDE between MATLAB and ZEMAX, we can get correction result for any Based on interface the DDE interface between MATLAB and ZEMAX, we dynamic can get dynamic correction result given Given as 1◦Δθ , 2◦ as and ratio after each correction cycle is shown Figurein 9. for any∆θ. given Δθ.∆θ Given 1°,4◦2°, the andStrehl 4°, the Strehl ratio after each correction cycle isinshown ◦ From Figure 9,Figure the variation of look angle time t should be lessbe than to keep ratio Figure 9. From 9, the variation of lookduring angle during time t should less2than 2° to Strehl keep Strehl betterbetter than than 0.8. The fluctuation region of Strehl ratioratio curve is caused by the change of the ratio 0.8. rapid The rapid fluctuation region of Strehl curve is caused by fast the fast change of aberration. From Figure 2, this region roughly corresponds to look angle fromfrom 25◦ to the aberration. From Figure 2, this region roughly corresponds to look angle 25°30to◦ .30°.

Figure Figure 9. 9. Strehl Strehl ratio ratio variation variation in in dynamic dynamic correction for different scanning rate.

5. 5. Experimental ExperimentalDemonstration Demonstration 5.1. System System Setup Although the optical system in Section 2 is working in the infraredregion, the algorithm can be demonstrated in the visible range. The experimental system is depicted in Figure 10. The The observation observation target illuminated by a white light source is a remote sensing image printed on a transparent illuminated by a white light source is a remote sensing image printed on a transparent film. film. DM1, a 37-channel micromachined micromachined membrane membrane DM, DM, is used to mimic dynamic aberrations induced by conformal mirror’s scanning raterate can can be reflected by the rate ofrate DM1. conformalwindow. window.The The mirror’s scanning be reflected byrefresh the refresh of DM2, DM1. aDM2, large-stroke DM with 52 electromagnetic actuators, is used forfor wavefront a large-stroke DM with 52 electromagnetic actuators, is used wavefrontcorrection correctionas aswell well as providing appropriate modal bias in aberration estimation. estimation. The The target target is is finally finally imaged imaged by a 12-bit charge-coupled charge-coupled device device (CCD) (CCD) camera. camera. The The camera camera and and DM2 DM2 are are both both connected connected to the computer to form a closed-loop control system.

Sensors 2016, 16, 1414 Sensors 2016, 16, 1414 1414

10 of 13 10 of 13

Figure 10. Experimental system layout. Figure 10. Experimental system layout. Figure 10. Experimental system layout.

5.2. System Error Clearance 5.2. System Error Clearance 5.2. System Errorsystem Clearance The initial errors induced by figure error of DMs and misalignment should be removed The initial system errors induced figure error of DMs and misalignment should be removed by DM itself before adding any specificby aberrations. Both and model-based WSLAO were by The initial system errors induced by figure error ofmodel-free DMs and misalignment should be removed DM itselfcalibrate before adding any specific Both model-free andThe model-based WSLAO were used used theadding system andspecific theaberrations. results are shown in Figure 11. model-free algorithm used by DMtoitself before any aberrations. Both model-free and model-based WSLAO were to calibrate the system and the results are shown in Figure 11. The model-free algorithm used here here is SPGD and image sharpness defined by Equation 9 is used as its metric function. In modelused to calibrate the system and the results are shown in Figure 11. The model-free algorithm used isbased SPGDalgorithm, and image sharpness defined byare Equation 9 isand usedcorrected as its metric In model-based first 15 Lukosz modes estimated usingfunction. “3N algorithm”. The here is SPGD and image sharpness defined by Equation 9 is used as its metric function. In modelalgorithm, 15 Lukosz modes are estimated and corrected using algorithm”. variation offirst image sharpness in SPGD and model-based algorithm are“3N shown in Figure The 12: variation of based algorithm, first 15 Lukosz modes are estimated and corrected using “3N algorithm”. The image sharpness in SPGD and model-based algorithm are shown in Figure 12: variation of image sharpness in SPGD and model-based algorithm are shown in Figure 12: 2 s I 2( x, y)dxdy  J (9) I ( x, y)dxdy 2 J = (I I2 ((xx,, yy)) (9) ) dxdy 2  J  (∑ I ( x, y))2 (9) ( I ( x, y)) whereII(x, (x,y) y) is is the the intensity intensity distribution where distributionat atimage imageplane. plane.

where I (x, y) is the intensity distribution at image plane.

5 ; (b) Image after SPGD correction, −5;−(b) Figure (a) Image Image with with initial initial aberration, aberration,J J==1.385 1.385× ×1010 Figure 11. 11. (a) Image after SPGD correction, − 5 5 . Image size: 128 × 128. −5− JJ== 1.459 10−5; (c) ; (c)Image Imageafter aftermodel-based model-basedWSLAO WSLAOcorrection, correction,J J==1.457 1.457××1010 1.459 × × 10 . Image size: 128 × 128.

Figure 11. (a) Image with initial aberration, J = 1.385 × 10−5; (b) Image after SPGD correction,

FromFigures Figures 11 and12, 12,the thecorrection correctionresult resultofofSPGD SPGD slightly better than that of modelFrom is is slightly better than that of thethe model-based −5. Image J = 1.459 × 10−5;11 (c)and Image after model-based WSLAO correction, J = 1.457 × 10 size: 128 × 128. based WSLAO, however, at the cost of more DM deformation number and more computation time. WSLAO, however, at the cost of more DM deformation number and more computation time. −5 From Figure 12, the image sharpness in the 3N algorithm increases from 1.385 × 10 to 1.457 × 10−5 − 5 FromFrom Figure 12, the image sharpness in the 3N algorithm increases from 1.385 × 10 to 1.457 × 10−5 Figures 11 and 12, the correction result of SPGD is slightly better than that of the modelafter 90 DM deformations. However, 200 DM deformations are required in SPGD to achieve the same after DM deformations. 200more DM DM deformations are required in SPGD achieve the same based90WSLAO, however, atHowever, the cost of deformation number and moretocomputation time. improvement, so the convergence speed of the -based WSLAO (3N algorithm) is roughly two times −5 improvement, sothe theimage convergence speed of the WSLAO (3N algorithm) two ×times From Figure 12, sharpness in the 3N-based algorithm increases from 1.385 is × roughly 10−5 to 1.457 10 faster than SPGD. In Figure 12, we use the DM deformation number to evaluate the convergence faster than SPGD. In Figure 12, we use the DM deformation number to evaluate the convergence speed after 90 DM deformations. However, 200 DM deformations are required in SPGD to achieve the same speed instead of using iteration number which can be ambiguous for different algorithms. In modelinstead of using iteration number which canthe be ambiguous for different algorithms. In model-based improvement, the is convergence speed WSLAO (3N algorithm) based WSLAO,soDM deformed 2N + 1 of or 3N -based times in one iteration. However,isinroughly SPGD, two DM times is WSLAO, DM is deformed 2N + 1 or 3N times in one iteration. However, in SPGD, DM is deformed faster than SPGD. In Figure 12, we use the DM deformation number to evaluate the convergence only in one iteration. Thenumber correction accuracy “3N” algorithm is alsoalgorithms. higher thanIn“2N + 1” speedtwice instead of using iteration which can beofambiguous for different modelalgorithm from Figure based WSLAO, DM is12. deformed 2N + 1 or 3N times in one iteration. However, in SPGD, DM is

Sensors 2016, 16, 1414 Sensors 2016, 16, 1414

11 of 13 11 of 13

ImageSharpness Sharpness(J) (J) Image

deformed only twice deformed twice in in one one iteration. iteration. The The correction correction accuracy accuracy of of “3N” “3N” algorithm algorithm is is also also higher higher11than than Sensors 2016, only 16, 1414 of 13 “2N + 1” algorithm from Figure 12. “2N + 1” algorithm from Figure 12.

Figure Figure 12. Image sharpness sharpness vs. vs. DM deformation deformation number number in in SPGD SPGD and model-based WSLAO. Figure 12. 12. Image Image sharpness vs. DM DM deformation number in SPGD and and model-based model-based WSLAO. WSLAO.

5.3. 5.3. Dynamic Results Dynamic Correction Correction Results After After removing removing the the system’s system’s initial initial error, error, dynamic dynamic aberrations aberrations of of the the conformal conformal window window across across full FOR as given in Section 2 are sampled with 1 degree interval and then reproduced by DM1 full FOR as given in Section 2 are sampled with 1 degree interval and then reproduced by DM1 with with update update rate rate of of 10 10 Hz, Hz, corresponding corresponding to to aa scanning scanning rate rate of of 10 10 degrees degrees per second. The time time delay delay in in degrees per per second. second. The correction of one mode is about 90 milliseconds, which is limited by our image acquisition hardware. correction of one mode is about 90 milliseconds, which is limited limited by by our our image image acquisition acquisition hardware. hardware. DM2 for the aberration and the results are in 13. DM2 is isused usedto compensate aberration the correction results are shown in Figure 13. is used totocompensate compensate forfor thethe aberration andand the correction correction results are shown shown in Figure Figure 13. Here Here the “3N” algorithm is used and the image sharpness as defined in Equation (9) is used to evaluate Here the “3N” algorithm is used image sharpness defined Equation(9) (9)isisused used to to evaluate the “3N” algorithm is used andand the the image sharpness as as defined ininEquation the the image image quality quality after after correction correction of of each each mode. mode. Dominant are Dominant Lukosz Lukosz modes modes as as shown shown in in Figure Figure 66 are corrected in optimized correction. First 15 Lukosz modes are corrected in un-optimized corrected in optimized correction. First First 15 15 Lukosz Lukosz modes modes are corrected corrected in un-optimized un-optimizedcorrection. correction. correction. Images of 128 ×× 128 Images are are captured captured with with the the same same size 128 as as aa tradeoff tradeoff between between correction correction speed speed and and same size size of of 128 128 × accuracy. From Figure 13, using the optimized method can achieve convergence much faster than unaccuracy. From FromFigure Figure13, 13,using usingthe theoptimized optimized method can achieve convergence much faster method can achieve convergence much faster thanthan un−5 optimized method. The mean value of image sharpness during dynamic correction is 1.436 × 10 −5 − un-optimized method. mean value image sharpness during dynamiccorrection correctionisis1.436 1.436×× optimized method. TheThe mean value of of image sharpness during dynamic 1010 in in5 −5 optimized correction and is ×× 10 in correction. 5unoptimized in optimized correction and is 1.427 ×−510 in unoptimized correction. optimized correction and is 1.427 1.427 10 in−unoptimized correction. 1.44 1.44

10 -5-5 10

1.43 1.43

ImageSharpness Sharpness Image

1.42 1.42 1.41 1.41 1.4 1.4 1.39 1.39 1.38 1.38 Optimized Correction Optimized Correction Un-optimized Correction Un-optimized Correction Before Correction Before Correction

1.37 1.37 1.36 1.360 0

5 5

10 10

15 15

20 20

25 25

Correction Number Correction Number

30 30

35 35

40 40

45 45

Figure Figure 13. 13. Dynamic Dynamiccorrection correction results results using using the optimized method and the unoptimized method. Figure 13. Dynamic correction results using the the optimized optimized method method and and the the unoptimized unoptimized method. method.

Sensors 2016, 16, 1414

12 of 13

6. Conclusions An infrared remote sensing system with conformal window and secondary deformable mirror was fully established. Model-based WSLAO was proposed and optimized for dynamic correction of aberrations induced by conformal window. To optimize the correction, only some dominant Lukosz modes are corrected and the image size is minimized to save the computing time of metric function. In experiment, dynamic aberration of conformal window with scanning rate of 10 degrees per second is reproduced by a 37-channel DM. A 52-channel DM is used for correction. For a 128 × 128 image, the mean value of image sharpness during dynamic correction is 1.436 × 10−5 in optimized correction and is 1.427 × 10−5 in un-optimized correction. We also demonstrated that the model-based WSLAO can achieve convergence two times faster than traditional stochastic parallel gradient descent (SPGD) method. Dynamic aberrations caused by aero-optical effects is not taken into account in this paper. Testing the performance of model-based WSLAO in correction of aero-optical aberration will be the subject of our future work. The temporal bandwidth of aero-optical aberrations surrounding the conformal window can be much higher than that of dynamic aberrations induced by window shape and scanning mirror. The temporal bandwidth of aero-optical aberrations is on the order of KHz for boundary layer turbulence. It is much more challenging to compensate aero-optical aberrations by WSLAO system. The wavefront control algorithm must be optimized and more powerful hardware should be used. To accelerate the correction of model-based WSLAO, dominant aberration modes induced by aero-optical effects should be determined in advance by simulation or experiment. Micro-electro-mechanical-system (MEMS) DM with high response speed (up to 50 KHz) [22] can be used to further decrease the time delay. Acknowledgments: This work was supported by the National Natural Science Foundation of China (Nos. 61505008, 11304012). We thank the anonymous reviewers for their valuable comments which significantly improved the paper. Author Contributions: B.D. and Y.L. conceived and designed the experiments; X.H. performed the experiments; B.D. and X.H. analyzed the data; B.H. contributed analysis tools; B.D. wrote the paper. Conflicts of Interest: The authors declare no conflict of interest.

References 1. 2. 3. 4. 5. 6. 7. 8. 9. 10.

Mills, J.P. Conformal Optics: Theory and Practice. In Proceedings of the Novel Optical Systems Design and Optimization IV, San Diego, CA, USA, 29 July 2001; pp. 101–107. Zhao, C.; Cui, Q.; Mao, S. Aberration and boresight error correction for conformal aircraft windows using the inner window surface and tilted fixed correctors. Appl. Opt. 2016, 55, 2626–2633. [CrossRef] [PubMed] Zhang, W.; Chen, S.; Hao, C.; Wang, H.; Zuo, B.; Fan, Z. Conformal dome aberration correction with gradient index optical elements. Opt. Express. 2014, 22, 3514–3525. [CrossRef] [PubMed] Zhang, W.; Zuo, B.; Chen, S.; Xiao, H.; Fan, Z. Design of fixed correctors used in conformal optical system based on diffractive optical elements. Appl. Opt. 2013, 52, 461–466. [CrossRef] [PubMed] Hao, Q.; Cheng, X.; Kang, J.; Jiang, Y. An image stabilization optical system using deformable freeform mirrors. Sensors 2015, 15, 1736–1749. [CrossRef] [PubMed] Li, Y.; Li, L.; Huang, Y.; Liu, J. Conformal optical design with combination of static and dynamic aberration corrections. Chin. Phys. B 2009, 18, 565–570. Li, Y.; Li, L.; Huang, Y.; Du, B. Correcting dynamic residual aberrations of conformal optical systems using AO technology. Chin. Phys. B 2009, 18, 2769–2773. Lombardo, M.; Serrao, S.; Devaney, N.; Parravano, M.; Lombardo, G. Adaptive optics technology for high-resolution retinal imaging. Sensors 2013, 13, 334–366. [CrossRef] [PubMed] Vorontsov, M.A.; Carhart, G.W. Adaptive optics based on analog parallel stochastic optimization: Analysis and experimental demonstration. J. Opt. Soc. Am. A 2000, 17, 1440–1453. [CrossRef] Gao, Q.; Jiang, Z.; Yi, S.; Xie, W.; Liao, T. Correcting the aero-optical aberration of the supersonic mixing layer with adaptive optics: concept validation. Appl. Opt. 2012, 51, 3922–3929. [CrossRef] [PubMed]

Sensors 2016, 16, 1414

11. 12. 13. 14. 15. 16. 17. 18. 19.

20. 21. 22.

13 of 13

Booth, M.J. Wave front sensor-less adaptive optics: A model-based approach using sphere packings. Opt. Express. 2006, 14, 1339–1352. [CrossRef] [PubMed] Débarre, D.; Booth, M.J.; Wilson, T. Image based adaptive optics through optimisation of low spatial frequencies. Opt. Express. 2007, 15, 8176–8190. [CrossRef] [PubMed] Dong, B.; Wang, R. Laser beam cleanup using improved model-based wavefront sensorless adaptive optics. Chin. Opt. Lett. 2016, 14, 031406. [CrossRef] Wei, H.; Su, G.J. Design and fabrication of a large-stroke deformable mirror using a gear-shape ionic-conductive polymer metal composite. Sensors 2012, 12, 11100–11112. [CrossRef] [PubMed] Hossain, M.M.; Bin, W.; Kong, S.H. Large-stroke convex micromirror actuated by electromagnetic force for optical power control. Opt. Express. 2015, 23, 28358–28368. [CrossRef] [PubMed] ALPAO Deformable Mirrors. Available online: http://www.alpao.com/adaptive-optics/deformablemirrors.html (accessed on 20 June 2016). Tyson, R.K. Introduction to Adaptive Optics; SPIE Press: Bellingham, WA, USA, 2000. Conrad, R.A.; Wilcox, W.E.; Williams, T.H.; Michael, S.; Roth, J.M. Emulation of dynamic wavefront disturbances using a deformable mirror. Opt. Express. 2009, 17, 3447–3459. [CrossRef] [PubMed] Hao, C.; Chen, S.; Zhang, W.; Ren, J.H.; Li, C.; Pang, H.J.; Wang, H.H.; Liu, Q.; Wang, C.; Zou, H.Y.; et al. Comprehensive analysis of imaging quality degradation of an airborne optical system for aerodynamic flow field around the optical window. Appl. Opt. 2013, 52, 7889–7898. [CrossRef] [PubMed] Vogel, C.R.; Tyler, G.A.; Wittich, D.J. Spatial–temporal-covariance-based modeling, analysis, and simulation of aero-optics wavefront aberrations. J. Opt. Soc. Am. A 2014, 31, 1666–1679. [CrossRef] [PubMed] Goorskey, D.J.; Drye, R.; Whiteley, M.R. Dynamic modal analysis of transonic Airborne Aero-Optics Laboratory conformal window flight-test aero-optics. Opt. Eng. 2013, 52, 071414. [CrossRef] BMC Deformable Mirrors. Available online: http://www.bostonmicromachines.com/deformable-mirrors. html (accessed on 20 June 2016). © 2016 by the authors; licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC-BY) license (http://creativecommons.org/licenses/by/4.0/).