Wavefront Sensing II

1
Wavefront Sensing II

• Richard Lane
• Department of Electrical and Computer Engineering
• University of Canterbury

2
Contents

• Session 1 Principles
• Session 2 Performances
• Session 3 Wavefront Reconstruction

3
Session 2 Performances

• Geometrical wavefront sensing take 2
• The inverse problem
• The astronomical setting
• The basic methods

4
Geometric wavefront sensing(or curvature sensing
without curvature)
Plane 1
Image Plane
Plane 2
Improve sensitivity (signal stronger)
Improve the number of modes measurable (signal
weaker)
5
Geometric optics model

• Slopes in the wave-front causes the intensity
  distribution to be stretched like a rubber sheet

• Aim is to map the distorted
• distribution back to uniform

6
Geometric wavefront sensingTake 2
Intensity Plane 1
Plane 1
Image Plane
Intensity Plane 2
Plane 2
Intensity distribution gives the probability
distribution For the photon arrival
7
Recovering the phase
Intensity Plane 1
Integrate to Form CDF
Intensity Plane 2
Choose level
Probability density functions
Integrate slope to find the phase
Defocus!
8
Forward Problem
9
Inverse Problem
Performance is determined by amount of photons
entering the aperture and assumptions about the
object and turbulence
10
Imaging a star
 
 
11
Multiple layers
For wide angle imaging we need to know the
height of the turbulence
Layer 1
h1
h2
Layer 2
Aperture Plane
12
The fundamental problem

• How to optimally estimate the optical effects of
  turbulence from a minimal set of measurements

13
Limiting Factors

• Technological
• CCD read noise
• Design of wavefront sensor (Curvature,
  Shack-Hartmann, Phase Diversity)
• Fundamental
• Photon Noise
• Loss of information in measurements
• Quality of prior knowledge

14
In Its Raw Form the Inverse Problem Is Always
Insoluble

• There are always an infinite number of ways to
  explain data.
• The problem is to explain the data in the most
  reasonable way
• Example Shack-Hartmann sensing for estimating
  turbulence

15
Example fit a curve to known slopes

• Solution requires assumptions on the nature of
  the turbulence
• Use a limited set of basis functions
• Assume Kolmogorov turbulence or smoothness

16
Parameter estimation

• Essentially we need to find a set of unknown
  parameters which describe the object and/or
  turbulence
• The parameters can be in terms of pixels or
  coefficients of basis functions
• Solution should not be overly sensitive to our
  choice of parameters.
• Ideally it should be on physical grounds

17
Bayesian estimation 101An important problem

• Estimate

• And if you know that it models two people
• splitting the bill in a restaurant?

18
Possible phase functionsZernike basis
?
Zernike Polynomials Low orders are smooth
Pixel basis, highest frequency 1/(2?)
19
Estimation using Zernike polynomials
phase weighting Zernikie polynomial

• Measurement Interaction Zernike Polynomial
• vector matrix
  Coefficents

• ith column of T corresponds to the measurement
  that would occur if the phase was the ith Zernike
  polynomial

20
Extension to many modes

• Provided the set of basis functions is complete,
  the answer is independent of the choice
• The best functions are approximately given by the
  eigenfunctions of the covariance matrix C
• These approximate the low order Zernike
  polynomials, hence their use.
• Conventional approach is to use a least squares
  solution
• and estimate only the first M Zernikes when M
  N/2 (N is the number of measurements)

21
Ordinary least squares

• Minimise

Weighted least squares

• Not all measurements are equally noisy
• hence minimise

22
Conventional Results

• As the number M increases the wavefront error
  decreases then increases as M approaches N.
• Reason when MN there is no error and there
  should be as higher order modes exist and will be
  affecting the measurements

23
Phase estimation from the centroid

• Tilt and coma both produce displacement of the
  centroid
• According to Noll for Kolmogorov turbulence
• Variance of the tilt
• Variance of the coma

Ideally you should estimate a small amount of coma
24
Bayesian viewpoint

• The problem in the previous slide is that we are
  not modelling the problem correctly
• Assuming that the higher order modes are zero, is
  forcing errors on the lower order modes
• Need to estimate the coefficients of all the
  modes as random variables

25
Example of Bayesian estimation for
underdetermined equations

• Measurement z is a linear function of two
  unknowns x,y

• The estimate (denoted by ) is a linear
• function of z

• We want to minimise the expected error

Statistical expectation
26
Minimisation of the error

• Key step, rewrite in terms of and

• Solution is a function of the covariance of the
  unknown
• parameters

27
Vector solution for the phase

• Express the phase as a sum of orthogonal basis
  functions
• Observed measurements are a linear function of
  the coefficients
• Reconstructor depends on the covariance of a

28
Simple example for tilt D/r04

• From Noll
• From Primot et al

29
Bayesian estimate of the wavefront
Minimizes
30
Summary Bayesian method

• When the data is noisy you need to put more
  emphasis on the prior.
• For example, if the data is very bad, dont try
  and estimate a large number of modes
• When done properly the result does not depend
  strongly on C being exact
• Error predicted to be
• where

31
Operation of a Bayesian estimator

• Minimizes
• When D becomes very large, the data is very noisy
  then more weight is placed on the prior
• data prior
• Ultimately as D?8, a?0 (for very noisy data no
  estimate is made)

32
Bayesian examination question

• You are on a game show.
• You can select one of three doors
• Behind one door is 10000, behind the others
  nothing
• After you select a door, the compere then opens
  one of the other doors revealing nothing.
• You are given the option to change your choice
• Should you?

33
Estimating the performance limits when it is
non-Gaussian

• The preceding analysis is fine when the
  measurement errors can be modelled as a Gaussian
  random variable
• On many equations you need to perform an analysis
  to work out the error in the analysis
• Cramer-Rao bounds

34
Cramer-Rao bound

• Linear unbiased estimators only
• Essentially the quality of the parameter estimate
  is given by the curvature of the pdf
• Doesnt tell you how to achieve the bound

35
Simple example

• Find the performance limit estimating the mean of
  a one-dimensional Gaussian from 1 sample

36
Points to note

• Limit is a lower bound. Clearly for 1 sample from
  the pdf it cannot be attained
• The variance decays as 1/N with more samples
• For a Gaussian asymptotically the centroid of the
  distribution can be shown to approach the
  Cramer-Rao bound

37
Estimation of a laser guidestar location,
Cramer-Rao bound
Small projection telescope
Large AO corrected projection telescope
Large uncorrected projection telescope
Key points In the presence of saturation a
focused spot may not be optimal Need to know the
pattern to reach the limit
38
Optimal estimation of a parameterwavefront tilt

• Important because the wavefront tilt is the
  dominant form of phase aberration
• A small error in estimating the tilt can be
  larger than the full variance of a higher order
  aberration.

39
Issues

• Displacement of the centroid of an image is
  proportional to the average tilt (not the least
  mean square) of the phase distortion
• Will discuss this issue later, but for the moment
  concentrate on estimating the mean square tilt.

40
How do you estimate the centre of a spot?

• The performance of the Shack-Hartmann sensor
  depends on how well the displacement of the spot
  is estimated.
• The displacement is usually estimated using the
  centroid (center-of-mass) estimator.
• This is the optimal estimator for the case where
  the spot is Gaussian distributed and the noise is
  Poisson.

41
Centroid estimation for a sinc2 function
42
Why Not Use the Centroid?

• In practice the spot intensity decays as
• This means that photons can still occur at points
  quite distant from the centre.
• Estimator is divergent unless restricted to a
  finite region in the image plane

43
Diffraction-limited spot

• For a square aperture, the distribution is

44
Photon arrival simulation
45
Solutions (1)

• Use a quad cell detector and discard the photons
  away from the centre
• The signal from the outer cells is discarded
  because it adds too much noise

46
Solutions 2

• Use an optimal estimator that weights the
  information appropriately
• Consider two measurements of an unknown parameter
  an estimate of a parameter with different
  variances
• A weighted sum is always a better estimator
• A non linear estimator is better still

47
Maximum-likelihood estimation

• If photons are detected at x1, x2, xN, the
  estimate is the value that maximizes the
  expression
• The Cramer-Rao lower bound for the variance is
• For a large number of photons, N, the variance
  approaches the Cramer-Rao lower bound.

48
Centroid location by model fitting

• Technique relies on finding a model of the object
• Not sensitive to the size of window (unlike the
  centroid)
• Centroid is a closed form
  solution for fitting a Gaussian of variable width

49
Tilt estimation in curvature sensing

• The image is displaced by the atmospheric tilt,
  how well you can estimate it is determined by the
  shape of the image formed.

50
Tilt estimation in the curvature actual
propagated wavefronts
51
Performance versus detector positionfor a
curvature sensor
52
Actual Wavefront sensor data

• Observation at Observatoire de Lyon
• SPID instrument on 1-m telescope
• 20x20 Shack-Hartmann lenslet array
• Exposure time 2ms
• Objects Pollux, point object 2500 frames
• Castor3 arc second binary 2500
  frames

53
Centroiding issues

• Accuracy required to a fraction of a pixel
• Sampling rate 60 of Nyquist

54
Finding the model

• Need a good model of the object
• In each lenslet of the Shack-Hartmann acts like a
  small telescope the dominant effect is one of
  tilt.
• gt We have a large number of images of the same
  object shifted before they are sampled.

55
Solution approach

• Use blind deconvolution to find model
• MAP framework (Hardie et al, FLIR)
• Data-capturing process
• Choose initially so that
• Prior information
• Laplacian smoothness for the optics
• Maximum entropy for CCD

56
Typical SPID data frames
Single Wavefront Sensor Frame
Long term WSF
Blow up of a spot
Movie of a spot
57
Simulations

• Inputs
• Object f point source
• Optics diffraction-limited pattern of square
  aperture
• CCD structure Gaussian-like
• Random displacements
• White Gaussian noise dB, 30dB, 15dB

58
Simulation result 15dB noise
Optics reconstruction
CCD reconstruction
59
Traditional centroiding

• Centre of gravity of spot image
• Problems
• Finite pixel size
• Finite window size
• Readout noise (more pixels more noise)
• Bias
• Problems become worse with extended objects

60
Model-fitting

• Full blind deconvolution computationally
  unreasonable
• Fit a model estimated by blind deconvolution
• Use model to determine centroids

61
Error in centroid calculation
62
Blind deconvolution results
Optics reconstruction
CCD reconstruction
63
Results from speckle image deconvolution
(narrowband)
Binary estimated with model fitted centroids
Binary estimated with traditional centroids
64
Phase reconstructions of Binary
Model based centroiding
Traditional centroiding
65
Conclusions

• Bayesian approaches provide a logical framework
  for filling in missing data
• Make sure of what you are assuming
• Cramer-Rao bound can provide a performance limit
• You need to look at the whole process when
  deriving an algorithm

66
And the answer is(ref Stark and Woods)

• Yes change the door

?
?
?
67
Actual Wavefront sensor data

• Observation at Observatoire de Lyon
• SPID instrument on 1-m telescope
• 20x20 Shack-Hartmann lenslet array
• Exposure time 2ms
• Objects Pollux, point object 2500 frames
• Castor3 arc second binary 2500
  frames

68
Subpixel displacement estimation
Wavefront sensing is based on estimating the
tilts produced by atmospheric distortion, the
accuracy of displacement estimation is critical.
Estimated optics psf
Data from SPID 2500 frames undersampled by 40
Spot displacements
Estimated CCD pixel sensitivity
69
Explanation of the terms

• Results from

70
Possible phase functionsZernike basis
71
The inverse problem
72
Alternatively
73
Prior information

• Infinite number of unknowns, but a finite number
  of centroid measurements from the sensor
• Conventional approach is to choose the basis
  functions and estimate M coefficients, where M lt
  N the number of measurements

74
Using real data binary star

• 14x14 Shack-Hartmann lenslet array
• Exposure time 3.2ms
• Object Castor, a binary star
• Intensity ratio 2.1
• Separation 3.1 arcseconds

75
Blind deconvolution results

• Intensity ratio 2.4
• Separation 3 arc seconds