Plotting the uv-tracks for CS 791 A2

uv_tracks.py

#!/usr/bin/env python3

import numpy as np
import matplotlib.pyplot as plt


def plot_uv(antennas):
    """Plots the uv-tracks for antennas."""

    # get baseline for antenna's 1-4
    baseline = antennas[3] - antennas[2] - antennas[1] - antennas[0]

    # calculate distance
    sqx = np.square(baseline[0])
    sqy = np.square(baseline[1])
    sqz = np.square(baseline[2])
    dist = np.sqrt(sqx + sqy + sqz)

    # convert delta_0 to radians
    delta_0 = (-74. -39./60 -37.481/3600) / 180 * np.pi

    # calculate azimuth and elevation angles
    azimuth = np.arctan(baseline[0] / baseline[1])
    elevation = np.arcsin(baseline[2] / dist)

    # calculate latitude
    L = (-30. -43./60 -17.34/3600) / 180 * np.pi

    # calculate X, Y and Z from lecture slides
    X = dist * (np.cos(L) * np.sin(elevation) - np.sin(L) * np.cos(elevation) * np.cos(azimuth))
    Y = dist * np.cos(elevation) * np.sin(azimuth)
    Z = dist * (np.sin(L) * np.sin(elevation) + np.cos(L) * np.cos(elevation) * np.cos(azimuth))

    # pack into matrix
    XYZ = np.matrix([[X], [Y], [Z]])

    # calculate wavelength using speed of light
    c = 3e8
    frequency = 1.4e9
    wavelength = c / frequency

    # calculate misc
    sq_lambda = np.sqrt(np.square(X) + np.square(Y)) * (1 / wavelength)
    sin_sq_lambda = abs(np.sin(delta_0)) * sq_lambda
    cos_z = np.cos(delta_0) * Z * (1 / wavelength)

    H = np.linspace(-14, 14, 100)
    H = H * 15 / 180 * np.pi

    tracks = np.matrix([[np.sin(H), np.cos(H), 0], [-np.sin(delta_0) * np.cos(H), np.sin(delta_0) * np.sin(H), np.cos(delta_0)]])

    uv = tracks * XYZ / wavelength
    u = uv[0, 0]
    v = uv[1, 0]

    # plot tracks
    plt.plot(u, v, label="$b_{14}^{xyz}$")
    plt.plot(-1 * u, -1 * v, label="$b_{41}^{xyz}$")
    plt.xlabel("u (in rad$^{-1}$)")
    plt.ylabel("v (in rad$^{-1}$)")
    plt.title("$uv$-tracks generated over 24 hours")
    plt.legend()
    plt.grid('on')
    plt.show()


if __name__ == "__main__":
    # load in antennas
    antennas = np.array([[ 25.095, -9.095,  0.045],
                         [ 90.284, 26.380, -0.226],
                         [  3.985, 26.839,  0.000],
                         [-21.605, 25.494,  0.019]])

    plot_uv(antennas)