In [79]:
# Uses a rotation matrix to plot the focal plane in mm
def plot_focal_plane():
'''
Plots the focal plane with the 5000 fiber positioners and 10 GFAs
using the initial four corners of the 0th indexed GFA and rotating the
points by 36 degrees counter-clockwise. Uses the reference projection of
the active area for each GFA.
'''
%pylab
import desimodel.io
# Sets the title of the graph
title('Focal Plane Overhead View')
# Plots the location of each of the fiber positioners
fp = desimodel.io.load_fiberpos()
plot(fp['X'],fp['Y'], 'g.')
#x = [318.529, 330.901, 348.947, 336.574]
#y = [225.702, 234.691, 209.830, 200.841]
# Experiments with "Reference projection of active area" coordinates
#x = [318.703, 331.075, 349.121, 336.748]
#y = [225.816, 234.805, 209.944, 200.955]
""" Uses the x and y from the petal indexed at 9 so the first petal
added to the table is indexed at 0
[[ 313.24842144 -233.32358331]
[ 325.62062672 -242.31230077]
[ 307.55293135 -267.15753816]
[ 295.18041705 -258.16786964]]"""
x = [313.24842144, 325.62062672, 307.55293135, 295.18041705]
y = [-233.32358331, -242.31230077, -267.15753816, -258.16786964]
# Creates a rotation matrix for 36 degrees counter-clockwise
rotatemat = numpy.zeros(shape=(2,2))
rotatemat[0] = [cos(36*pi/180), -sin(36*pi/180)]
rotatemat[1] = [sin(36*pi/180), cos(36*pi/180)]
draw_gfas(x, y, rotatemat)
# Function that draws the GFAs on the focal plane
def draw_gfas(x, y, rotatemat):
"""
Draws the 10 GFAs given the initial x and y coordinates of four corners of a GFA
Parameters
----------
x : Array of four x initial coordinates of the GFA
y : Array of four y initial coordinates of the GFA
"""
coord = numpy.zeros(shape=(2,1))
gfacoord = numpy.zeros(shape=(4, 2))
oldxcoord = x
oldycoord = y
for j in range(10):
for i in range(4):
coord[0] = oldxcoord[i]
coord[1] = oldycoord[i]
newcoord = matmul(rotatemat, coord)
oldxcoord[i] = newcoord[0]
oldycoord[i] = newcoord[1]
gfacoord[i] = [newcoord[0], newcoord[1]]
plot(newcoord[0], newcoord[1], 'k.')
draw_single_gfa(gfacoord)
def draw_single_gfa(gfacoord):
"""
Draws a single GFA given a 4X2 array of coordinates for the four corners of a GFA
Parameters
----------
gfacoord: 4X2 array of x and y coordinates with each row representing a corner of the GFA
"""
# Prints all of the GFA coordinates for debugging
#print(gfaCoord)
gfapolygon = Polygon(gfacoord)
plt.gca().add_patch(gfapolygon)
In [80]:
plot_focal_plane()
In [81]:
# Import statements to run code for debugging
import desimodel.io
import scipy.interpolate
def build_gfa_table(outfile = 'gfa.ecsv'):
'''
Builds the GFA table given the data from DESI-0530-v13 Excel spreadsheet
and writes a .ecsv file using the astropy table library. The data is pulled from
the "GFALocation" tab on the spreadsheet and from rows 16-23 and columns A-I.
Parameters
----------
outfile: a default parameter that represents the desired filename which is returned by
this function. The filename defaults to "gfa.ecsv" if no parameters are given.
'''
# Uses the reference projection of active area to create data table of GFAs
from astropy.table import Table
# Initial x and y coordinates for the GFAs
""" Uses the x and y from the petal indexed at 9 so the first petal
added to the table is indexed at 0
[[-125.10482863 -370.01790486]
[-129.83038525 -384.56223777]
[-159.04283509 -375.05643893]
[-154.31646944 -360.51151824]]
"""
# Data obtained from DESI-0530-v13 Excel spreadsheet
x = [-125.10482863, -129.83038525, -159.04283509, -154.31646944]
y = [-370.01790486, -384.56223777, -375.05643893, -360.51151824]
z = [-17.053, -18.487, -18.631, -17.198]
rotatemat = numpy.zeros(shape=(2,2))
rotatemat[0] = [cos(36*pi/180), -sin(36*pi/180)]
rotatemat[1] = [sin(36*pi/180), cos(36*pi/180)]
# Note: the corners are 0 indexed
gfatable = Table(names = ('PETAL', 'CORNER', 'X', 'Y', 'Z', 'Q', 'RADIUS_DEG', 'RADIUS_MM'),
dtype = ('int', 'int', 'float', 'float', 'float', 'float', 'float', 'float'))
# Sets the units for the GFA table
gfatable['X'].unit = 'mm'
gfatable['Y'].unit = 'mm'
gfatable['Z'].unit = 'mm'
gfatable['Q'].unit = 'degrees'
gfatable['RADIUS_DEG'] = 'degrees'
gfatable['RADIUS_MM'] = 'mm'
find_gfa_coordinates(x, y, z, gfatable, rotatemat)
# Saves the table of data as an ecsv file
gfatable.write(outfile, format='ascii.ecsv')
# Function that obtains the x and y coordinates for each corner of the GFAs
def find_gfa_coordinates(x, y, z, gfatable, rotatemat):
'''
Finds all the GFA coordinates by rotating the initial coordinates and adding
the respective coordinates to the gfaTable
Parameters
----------
x : Array of four x initial coordinates of the GFA
y : Array of four y initial coordinates of the GFA
z: Array of four z initial coordinates of the GFA
gfaTable: Astropy Table object which stores the petal number, corner number,
and x, y, and z coordinates in mm within each row
rotateMat: Rotation matrix to rotate the coordinates 36 degrees counterclockwise
'''
coord = numpy.zeros(shape=(2,1))
gfacoord = numpy.zeros(shape=(4, 2))
oldxcoord = x
oldycoord = y
for j in range(10):
for i in range(4):
coord[0] = oldxcoord[i]
coord[1] = oldycoord[i]
newcoord = matmul(rotatemat, coord)
oldxcoord[i] = newcoord[0]
oldycoord[i] = newcoord[1]
gfacoord[i] = [newcoord[0], newcoord[1]]
theta = cartesian_to_polar_angle(newcoord[0], newcoord[1])
# radius is the radius in mm
radius = np.sqrt(newcoord[0]**2 + newcoord[1]**2)
# degree is the radius in degrees
degree = get_radius_deg(newcoord[0], newcoord[1])
# Could be building the table in O(N^2), which is notably inefficient
gfatable.add_row([j, i, newcoord[0], newcoord[1], z[i], theta, degree, radius])
In [82]:
# Code from Sky Coordinates.ipynb for building the gfa table with polar coordinates
def get_radius_mm(theta):
"""
Returns the radius in mm given an radius in degrees using the platescale data
relative to the center of the focal plane as (0,0)
Parameters
----------
theta: A float that represents the angle from the center of the focal plane
"""
platescale = desimodel.io.load_platescale()
# Uses a quadratic one-dimensional interpolation to approximate the radius in degrees versus radius in mm
fn = scipy.interpolate.interp1d(platescale['theta'], platescale['radius'], kind = 'quadratic')
radius = float(fn(theta))
return radius
def get_radius_deg(x, y):
"""
Returns the radius in degrees given x, y coordinates using the platescale data
Parameters
----------
x: The x coordinate in mm of a location on the focal plane
y: The y coordinate in mm of a location on the focal plane
"""
radius = np.sqrt(x**2 + y**2)
platescale = desimodel.io.load_platescale()
fn = scipy.interpolate.interp1d(platescale['radius'], platescale['theta'], kind = 'quadratic')
degree = float(fn(radius))
return degree
def cartesian_to_polar_angle(x, y):
"""
Given cartesian coordinates, this function returns the polar angle in radians
for use in polar coordinates
Parameters
----------
x: The x coordinate in mm of a location on the focal plane
y: The y coordinate in mm of a location on the focal plane
"""
# Prevents division by 0
if(x == 0):
if(y > 0):
return 90
else:
return 270
# Case for quadrant 1
elif(x > 0 and y >= 0):
return arctan(y / x) * 180 / pi
# Case for quadrants 2 and 3
elif(x < 0):
return (arctan(y / x) * 180 / pi) + 180
# Case for quadrant 4
elif(x > 0 and y <= 0):
return (arctan(y / x) * 180 / pi) + 360
else:
raise ValueError
In [83]:
build_gfa_table()