Aspheric Surface Generator


Revision: 0.1

Last Update: 20131205

This notebook generates even asphereic surfaces using the definition in optical shop testing. The higher order terms use a normalized coordinate system such that the departure at the edge for the term is equal to the coefficient.

The resulting surface is witten in a 4sight compatible csv file.

Governing Equation

$$z = \dfrac{cS^2}{1+[1-(\kappa+1)c^2S^2]^{1/2}}+A_1S'^4+A_2S'^6+A_3S'^8+A_4S'^{10}$$

$~~~~~~~~~~\text{z} = \text{sag}$

$~~~~~~~~~~\text{c} = \text{curvature} = \frac{1}{\text{Radius of Curvature}}$

$~~~~~~~~~~\kappa = \text{Conic Constant}$

$~~~~~~~~~~S = \text{Radial distance in real units}$

$~~~~~~~~~~S' = \text{Normalized radial distance such that }S' = 1\text{ when }S = S_{max}$

$~~~~~~~~~~A_n = \text{Higher order term coefficient. Value is equal to contribution at the edge of the aperture.}$

Define Constants


In [1]:
# array size and center offset in pixels
arraySizeX = 1000.0
arraySizeY = arraySizeX
arrayOffsetX = -0.5
arrayOffsetY = -0.5

# Array length in mm
arrayLength = 1.1 * 20

# Wavelength (nm)
wavelength = 633.0

# Physical pupil diameter in mm
diameter = 20

# Radius of curvature in mm
radiusOfCurvature = 20

# Conic Constant
k = -1

# Higher Order Coeficients in Waves
a1 = 0.0
a2 = 0.0
a3 = 0.0
a4 = 0.0

# Radius of currvature in mm for test wavefront
radiusOfCurvatureTest = 1.0015 * radiusOfCurvature

Assign Base Plot Style


In [2]:
import json
s = json.load( open("../styles/bmh_matplotlibrc.json") )
matplotlib.rcParams.update(s)

Ensure entered values are in correct cast


In [3]:
arraySizeX = float(arraySizeX)
arraySizeY = float(arraySizeY)
arrayOffsetX = float(arrayOffsetX)
arrayOffsetY = float(arrayOffsetY)
arrayLength = float(arrayLength)
wavelength = float(wavelength)
diameter = float(diameter)
radiusOfCurvature = float(radiusOfCurvature)
radiusOfCurvatureTest = float(radiusOfCurvatureTest)
k = float(k)
a1 = float(a1)
a2 = float(a2)
a3 = float(a3)
a4 = float(a4)

Generate Output File Name


In [4]:
if k < 0:
    fileName = 'R' + str(int(radiusOfCurvature)) + '_D' + str(int(diameter)) + \
                   '_A' + str(int(arrayLength)) + '_Km' + str(int(abs(k))) + '_a1' + \
                       str(int(a1)) + '_a2' + str(int(a2)) + '_a3' + str(int(a3)) + \
                           '_a4' + str(int(a4)) + '_N' + str(int(arraySizeX)) + '_waves.csv'

    fileNameBFS = 'R' + str(int(radiusOfCurvature)) + '_D' + str(int(diameter)) + \
                   '_A' + str(int(arrayLength)) + '_Km' + str(int(abs(k))) + '_a1' + \
                       str(int(a1)) + '_a2' + str(int(a2)) + '_a3' + str(int(a3)) + \
                           '_a4' + str(int(a4)) + '_N' + str(int(arraySizeX)) + '_waves_BFS.csv'
                
else:
    fileName = 'R' + str(int(radiusOfCurvature)) + '_D' + str(int(diameter)) + \
                   '_A' + str(int(arrayLength)) + '_K' + str(int(k)) + '_a1' + \
                       str(int(a1)) + '_a2' + str(int(a2)) + '_a3' + str(int(a3)) + \
                           '_a4' + str(int(a4)) + '_N' + str(int(arraySizeX)) + '_waves.csv'
                
    fileNameBFS = 'R' + str(int(radiusOfCurvature)) + '_D' + str(int(diameter)) + \
                   '_A' + str(int(arrayLength)) + '_K' + str(int(k)) + '_a1' + \
                       str(int(a1)) + '_a2' + str(int(a2)) + '_a3' + str(int(a3)) + \
                           '_a4' + str(int(a4)) + '_N' + str(int(arraySizeX)) + '_waves_BFS.csv'

Define functions for calculating the asphereic surface and generating a normalized array


In [5]:
def GetAsphericValue(c, k, s, sNorm, a1, a2, a3, a4):
    return c * s**2 / (1 + sqrt(1 - (k + 1) * c**2 * s**2)) \
                 + a1 * sNorm**4 + a2 * sNorm**6 + a3 * sNorm**8 + a4 * sNorm**10

In [6]:
def GenerateNormalizedArray(arraySizeX, arraySizeY, arrayOffsetX, arrayOffsetY):
    radius = (min([arraySizeX, arraySizeY]) / 2.0)
    
    x = linspace((arrayOffsetX - arraySizeX / 2.0) / radius,
                      (arrayOffsetX + arraySizeX / 2.0) / radius, arraySizeX)
    y = linspace((arrayOffsetY - arraySizeY / 2.0) / radius,
                      (arrayOffsetY + arraySizeY / 2.0) / radius, arraySizeY)
    yy, xx = meshgrid(y,x)
    return sqrt(xx**2 + yy**2)

In [7]:
def EstimatedPower(zSmax, Smax):
    return 2 * zSmax / (Smax**2 + zSmax**2)

Generate the aspheric surface and the base power surface without the aspheric terms


In [8]:
# Generate array normalized to the edge of the array
s = GenerateNormalizedArray(arraySizeX, arraySizeY, arrayOffsetX, arrayOffsetY)

# Generate array normalized to the pupil
sNorm = s * arrayLength / diameter
sNorm = ma.masked_greater(sNorm, 1.0)

# Generate array of real lateral coordinates
s = s * arrayLength / 2.0
#s = ma.masked_greater(s,diameter / 2.0)
s = ma.masked_array(s, dtype='float',mask=sNorm.mask)

# Calculate the curvature in 1/mm
c = 1.0 / radiusOfCurvature

# Calcuate curvature in 1/mm for test wavefront
cTest = 1.0 / radiusOfCurvatureTest

# Convert coefficients from waves to mm
a1 = a1 * wavelength / (10**6)
a2 = a2 * wavelength / (10**6)
a3 = a3 * wavelength / (10**6)
a4 = a4 * wavelength / (10**6)

# Get aspheric surface in units of mm
surface = GetAsphericValue(c, k, s, sNorm, a1, a2, a3, a4)

cPow = EstimatedPower(GetAsphericValue(c, k, diameter / 2.0, 1, a1, a2, a3, a4), diameter / 2.0)

# Get base powered surface in units of mm
power = GetAsphericValue(cPow, 0, s, sNorm, 0, 0, 0, 0)

# Convert aspheric surface to units of waves
surface = surface * (10**6) / wavelength
power = power * (10**6) / wavelength

# Calculate departure from sphere and slope of departure between samples
delta = surface - power

Plot Results


In [9]:
#imshow(surface);
#title('Aspheric Surface');

In [10]:
plot(s[:,arraySizeX/2.0],surface[:,arraySizeX/2.0]);
xlim(0,arrayLength/2.0)
grid(b=True, which='major', color='k', linestyle='-')
title('Slice Through Aspheric Surface Center to Edge')
ylabel('Sag (waves)')
xlabel('Position (mm)');



In [11]:
plot(s[:,arraySizeX/2.0],(delta)[:,arraySizeX/2]);
xlim(0,arrayLength/2.0)
grid(b=True, which='major', color='k', linestyle='-')
title('Slice Through Aspheric Departure Center to Edge')
ylabel('Sag (waves)')
xlabel('Position (mm)');


Write 4Sight Compatible CSV File


In [12]:
# This strange maskFillValue was chosen to ensure that the header value and the 
# value in the array had exactly the same format.  
# A better programmer would have handled this functionally rather than choosing 
# a specific value format.
maskFillValue =  -123456.789123
surface.fill_value = maskFillValue
                                                                                             
fileHeader = 'xsize: ' + str(int(arraySizeX)) + '\nysize: ' + str(int(arraySizeY)) + \
                 '\ntype: phase\nbadpixel: ' + str(maskFillValue) + \
                     '\ntitle: Aspheric Data From iPython\nwedge: 0.50000\nxpix: ' + \
                         str(arrayLength/arraySizeX/1000.0) + \
                             '\naspect: 1.00000\ndelimiter: csv\nwavelength: ' + \
                                 str(wavelength) + '\n\n'
                      
np.savetxt(fileName, ma.filled(surface, maskFillValue), delimiter=',',fmt='%f',
               header = fileHeader, comments='')

In [13]:
# This strange maskFillValue was chosen to ensure that the header value and the 
# value in the array had exactly the same format.  
# A better programmer would have handled this functionally rather than choosing 
# a specific value format.
maskFillValue =  -123456.789123
surface.fill_value = maskFillValue
                                                                                             
fileHeader = 'xsize: ' + str(int(arraySizeX)) + '\nysize: ' + str(int(arraySizeY)) + \
                 '\ntype: phase\nbadpixel: ' + str(maskFillValue) + \
                     '\ntitle: Aspheric Data From iPython\nwedge: 0.50000\nxpix: ' + \
                         str(arrayLength/arraySizeX/1000.0) + \
                             '\naspect: 1.00000\ndelimiter: csv\nwavelength: ' + \
                                 str(wavelength) + '\n\n'
                      
np.savetxt(fileNameBFS, ma.filled(delta, maskFillValue), delimiter=',',fmt='%f',
               header = fileHeader, comments='')