Genetic Algorithm Phase Retrieval for the Systematic Image‐Based Optical Alignment Test Bed

Since the early 1980s, requirements for telescopes and beam directors with greater light‐gathering capability and improved angular resolution have driven primary‐mirror designs to diameters in excess of 8 m. The prohibitive costs and risks associated with fabrication, metrology, and handling of such large mirrors have driven telescope designers to construct large primary mirrors from arrays of smaller mirror segments. For such a primary mirror, a closed‐loop control system is necessary to maintain relative alignment and phasing of the segments during operation. This control system must maintain the primary‐mirror figure while the telescope is subjected to disturbances such as thermoelastic deformations, gravity loads, wind loading, and vibrations in the telescope structure. Several different methods have been explored and demonstrated for controlling arrays of segmented primary mirrors.

The segmented mirror design was first realized in the 10 m Keck I telescope on Mauna Kea in Hawaii (Nelson et al. 1985). Keck I has a 10 m primary mirror comprising 36 segments in three rings of hexagons. Capacitive sensors (commonly known as edge sensors) are mounted on the edges of the mirror segments and allow the determination of the relative misalignments of the segments. Each mirror has three degrees of freedom for control: tip, tilt, and piston. The control computer converts the measurements into motion commands to the three actuators on each segment. This closed‐loop control system precisely maintains the figure of the primary mirror, enabling high angular resolution. In the Keck control system, the control commands come from a global optimization involving 168 sensors and 108 actuators; i.e., a 108 × 168 control gain matrix. McDonald Observatory's Hobby‐Eberly Telescope (HET) has an 11 m primary mirror comprising 91 hexagonal segments in five rings (Ramsey et al. 1998). HET also uses a global optimization control system that employs 273 actuators and 480 inductive edge sensors (Rakoczy et al. 2002, 2003). It is easy to see that for future large telescopes and beam directors in the 30–100 m class, the control gain matrix will increase in complexity very quickly, and computation speed and numerical round‐off will become potential problems with the global optimization systems.

A different method, demonstrated in the Phased Array Mirror Extendible Large Aperture (PAMELA), sought to decentralize the global control strategy (Rather et al. 1989; Ames et al. 1995; Rakoczy et al. 2000). In the PAMELA test bed (a 36 hexagonal segment, 0.5 m diameter primary mirror), a Shack‐Hartmann wave‐front sensor was employed to sense tip and tilt motions of each segment. Each subaperture of the Shack‐Hartmann wave‐front sensor mapped exclusively to an individual primary‐mirror segment. Each segment had its own local processing electronics to convert tip‐tilt errors to the three actuator commands. Piston control was achieved utilizing inductive edge sensors; three on each segment. The local segment electronics only processed the sensors mounted on that particular segment. This was dubbed a "nearest neighbor" approach to piston control, since the piston control command was based exclusively upon the motion of the segment relative to three of its adjacent neighbors. Processing at the segment level, rather than utilizing the global strategy, reduced the number of arithmetic operations and enabled high‐bandwidth control of the primary mirror, an impracticality for the Keck telescopes. However, the PAMELA control system involves the complexity of having many electronic components, considerable board‐level software development, and communication among several levels of electronics modules.

These more conventional approaches to segmented mirror control were reevaluated in the mid‐1990s with the emergence of the Next Generation Space Telescope, now named the James Webb Space Telescope (JWST). The JWST will be the first space telescope to exploit a segmented primary mirror. The design adopted by the JWST team sought to maximize the number of photons used for science. This meant minimizing the number of photons used by fiduciary instruments, such as Shack‐Hartmann wave‐front sensors. Since the JWST will be inaccessible for frequent calibration, maintenance, or repair, outfitting the mirrors with many edge sensors, cables, and electronics modules is not desirable. Rather, the JWST will perform wave‐front sensing at the image plane of the science camera. Such a system requires no dedicated wave‐front sensors, with their accompanying weight, power consumption, photon consumption, complexity, and risk. This new approach will exploit resources that are already there: the science camera and science image data. However, this design requires the maturation of a technique known as phase retrieval (Redding et al. 1998).

Phase‐retrieval techniques have been investigated since the 1970s (Gerchberg & Saxton 1972; Gonsalves 1976).¹ The objective is to extract the wave‐front error (the phase of a complex‐valued function at the aperture plane) from an exclusively real‐valued function (the intensity distribution at the image plane). Within the framework of Fraunhofer diffraction theory, the image plane intensity is approximately the squared magnitude of the Fourier transform of the aperture plane function. Thus, phase retrieval is a problem in nonlinear estimation, which involves ambiguity in converging to a global optimum. Many variations of the technique have been demonstrated over the years, including phase diversity, which attempts to resolve the ambiguity (Gonsalves 1982). Simulations and experiments at Marshall Space Flight Center (MSFC) with the PAMELA test bed revealed that the phase‐retrieval technique was computationally intensive and slow to converge. In fact, convergence was so slow that it could not provide adequate bandwidth to keep up with wave‐front perturbations in PAMELA. Furthermore, the technique was shown to be sensitive to noise in the detector of the image camera.

We seek to investigate a genetic algorithm (GA) as a method for expediting and improving phase‐retrieval algorithms. The motivation for considering the GA as a solution technique is that their inherent parallelism could be exploited by specialized computing hardware to more rapidly calculate approximate solutions. A GA was applied to a simple test case targeted at the Systematic Image‐Based Optical Alignment (SIBOA) test bed at MSFC. The SIBOA test bed comprises a segmented primary mirror with seven circular segments (Fig. 1). The SIBOA configuration was selected in order to simplify the initial test case by minimizing the number of segments and consequently the number of parameters to be estimated. The test case was further simplified by neglecting the segment tip‐tilt degrees of freedom. With this assumption, only the relative piston errors among the seven segments needed to be sensed and corrected using a GA. In principle, tip‐tilt errors can be sensed and corrected as well. For instance, in an application such as the JWST, where photon loss is a high priority, tip‐tilt corrections could be made via an additional GA code in dedicated hardware. In a PAMELA‐type system, a Shack‐Hartmann wave‐front sensor could control the tip and tilt degrees of freedom in a decentralized manner, while the GA‐based phase‐retrieval scheme corrects the piston degrees of freedom. This architecture would be more favorable than the conventional PAMELA architecture, because no edge sensor hardware or software would be required, and the piston optimization would be global, not nearest neighbor.

Zoom In Zoom Out Reset image size

The problem of phase retrieval is that of estimating the real‐valued phase θ(u, v) of the pupil (or aperture) function

from measurements of intensity φ in the image plane, where

and the image function f(x,y) is given by

This must be done for each of the M mirrors in a segmented mirror system (M = 7 for the SIBOA test bed).

The image function and the pupil function make up a two‐dimensional Fourier pair:

In practice, only sampled values of the intensity distribution at the image plane are available, so the discrete Fourier transform (DFT) is used. Thus, using a DFT, we can estimate the image in the far field (i.e., at the detector) as

where N is the number of pixels on one side of an N × N pixel charge‐coupled device. Combining equations (1) and (6) and letting ω = e^2πi/N, the image function estimate can be written as

As a computational test case, uniform illumination was used. That is, the modulus of the pupil function was set to unity [F(u,v)| = 1] at each of the mirrors, and zero elsewhere. The independent correction of tip‐tilt provides a considerable simplification. Normally, θ(u, v) would be a continuously varying function of u and v for each mirror. Here, however, we are able to fix the phase differences θ₁,..., θ_M to a constant value across the surface of each mirror. Thus, the phase‐retrieval problem is one of estimating the M phases θ₁,..., θ_M so that the estimate of the image in the far field, given by equation (7), yields values approximating the measured image values f(x,y)|; i.e., the square roots of the intensity values. With these simplifications, the estimate of the image function takes the form

where

Here the index m determines the mirror of interest (where m = 1,...,M) and u and v are summed over the geometry of the mth mirror. The s^m_xy are completely determined by the characteristics of the optical system: the geometries of the aperture detector, the number of pixels, and wavelength of interest. The s^m_xy need be calculated only once. This technique results in approximately the same number of binary operations as a standard fast Fourier transform (FFT). It should be noted that while a GA was applied to the problem discussed here, this approximation technique is not specific to GAs; any optimization technique could use the same simplification.

As a demonstration of the technique, a real‐valued genetic algorithm as described in Coley (1999) was implemented using the MATLAB development environment and the Genetic and Evolutionary Algorithm Toolbox (GEATbx; Pohlheim 1999).² A real‐valued GA was chosen over a "canonical" binary‐valued algorithm (as described in Goldberg 1989 or Coley 1999), because it avoids the accuracy issues associated with the discrete nature of the binary representation and is in fact easier to implement. The genetic information for each "individual" consists of the M phases −π≤ θ_m≤π contained in a parameter vector (or "chromosome") of M elements ("genes"; recall that M = 7 for the SIBOA test bed). This GA is characterized as follows:

• 1.  
  Selection.—Rank‐based fitness selection is used in which the "most fit" parents from each generation are allowed to propagate their genetic information (the number of parents was varied to determine the optimal number of individuals allowed to propagate, as discussed below). Fitness is determined by generating an image using equation (8) and comparing it to a test image. Specifically, the absolute fitness is determined from ∥f(x,y)∣ - ∣f_est(x,y)∥, summed over the CCD pixels (smaller values are better). The individuals are then ranked according to this fitness criteria. In all cases, the highest ranked half of the population was allowed to propagate.
• 2.  
  Recombination.—Intermediate recombination is used for the exchange of genetic information as implemented in Pohlheim (1999). That is, two parents are randomly chosen from among the pool of fit individuals, and their chromosomes are allowed to exchange information. This is done by comparing the corresponding genes on each chromosome and generating a range for the random "offspring" that extends beyond the range between the two parents by 25% in each direction (as in Fig. 2).
• 3.  
  Mutation.—Mutation is the practice of introducing new, random genetic material. The mutation is such that each parameter has a probability of 1/M of being randomly changed to a value that is its current value ±10%. The probability distribution of the mutated parameter is uniform within this range.
• 4.  
  Reinsertion.—Each mating results in two offspring, and mating continues until the number of offspring in the new population is equal to the number of individuals in the original population.
• 5.  
  Elitism.—Elitism is the practice of allowing a few of the best individuals in each generation to move into the next generation as potential parents. In all cases, the fittest two individuals were kept in the population.

Zoom In Zoom Out Reset image size

This problem has a few unique characteristics in the application of the genetic algorithm technique. It must be noted that the overall solution for the M phases is unique only up to an overall phase shift (that is, an identical phase shift for all M phases). Thus, the genetic algorithm must really only find the relative phases of the mirrors. It is tempting in a case like this to fix the value of the phase of one of the mirrors to be zero; thus, the GA must only find M - 1, rather than M variables. However, this dramatically slows down the convergence of the algorithm, since it must find a unique solution, as opposed to simply finding the "shape" of the solution.

Another interesting characteristic is the periodic nature of the phase solutions. This can lead to difficulties in the convergence of the GA to the correct solutions if the GA attempts to converge to a phase outside the defined range (−π≤ θ_m≤π). One way to allow for this possibility is to not absolutely restrict convergence to these solutions. This can be implemented in the population by generating random phase values with a Gaussian distribution whose mean is zero and whose standard deviation is π. Thus, 68% of the phase values will fall in the desired range, while values outside this range are not forbidden; 95% of the values will fall in the range (- 2π,2π).

In order for a parallel GA phase‐retrieval technique to be a viable alternative to a more traditional iterative serial approach using Fourier transform methods such as the one used by Gonsalves (1976), it must meet two basic criteria. First, it must converge as fast or faster than a serial method to the correct set of relative phase shifts, and it must do this consistently for a variety of test cases. Second, it must be implementable using a minimum of additional hardware above that required for a serial technique (i.e., one processor), or else the GA becomes cost prohibitive.

The amount of additional hardware necessary for implementation of the GA in parallel is essentially determined by the number of individuals in the population. Evaluation of the fitness criterion is the most computationally demanding stage in the algorithm, since here it requires the evaluation of an FFT; thus, to speed calculation, each processor should handle just one individual. However, it should be noted that an increase in the number of processors in a parallel cluster represents an increase in the amount of system and communications overhead for the cluster, and thus a reduction in its computational efficiency.

To compare the computational efficiency of the GA to that of an iterative serial technique, an estimate of the GA's computational speed must be made. As a target value, the serial iterative technique described in Gonsalves (1976) typically requires on the order of 30 iterations for convergence; however, it requires one standard and one inverse two‐dimensional FFT per iteration, while the GA technique described here requires only one inverse two‐dimensional FFT per individual per generation (iteration). Because of the increased computational overhead of a parallel cluster, the GA technique will have to converge in a similar number of iterations to be competitive. If this could be accomplished with a population size in the neighborhood of 16 individuals (requiring 16 processors), then the GA might be a viable alternative to conventional techniques.

In optimizing the GA, the parameters that were considered were the mutation rate, the presence (or lack) of elitism, and the population size/number of mating parents. It was found that the standard mutation rate of 1/M worked well for the largest percentage of trial cases, as did a small amount of elitism (two individuals per generation); however, small changes in these parameters did not have a significant affect on the rate of convergence or the correctness of the solution. (Of course, large changes in these parameters invariably produced a detrimental effect).

The GA was found to be considerably more sensitive to the number of individuals in the population and the number of parents allowed to propagate. This observation is not surprising, since these values determine the amount of genetic information available in each generation: the number of individuals in the population determines the likelihood that the needed genetic information is present; too few parents, and the solution will converge quickly but with poor accuracy; too many parents, and the solution converges too slowly. It was determined that allowing 50% of the population to propagate was most effective. Again, small deviations from this value had little effect, but large ones were always detrimental.

Population sizes of 16, 32, 64, and 128 individuals were tested for several different simulated phase distributions. As might be expected, larger populations generally produced more accurate results. Unfortunately, a population size of 16 individuals almost never converged to the correct solution in our test cases. Most test cases required at least 32 individuals for consistent, accurate convergence within approximately 30 generations; a few required 64 or even 128 individuals.

While the genetic algorithm outlined here does successfully converge to the correct relative phases for each of the seven mirrors in the SIBOA system for several test cases, it does not do so more quickly or with greater computational efficiency than a traditional serial technique. Implementation of the GA would require a greater investment in hardware even to be competitive with a serial technique, without a significant increase in speed of the phase‐retrieval process.