# Settings
set xlabel "Apparent Magnitude"
set ylabel "Exposure time, in seconds"
set title "Exposure times with a background magnitude of 22"
set xrange [5:15]
set yrange [0.001:1000]
set logscale y
set nokey
set mxtics 4
set term postscript

A = 0.4064   # Aperture in meters
L = 500      # 500nm bandwidth
QE = 0.39    # Quantum efficiency
E = 0.7      # Optical efficiency
# x is the X variable, being apparent magnitude
# t is the Y variable, being exposure time
Vsky = 22    # Magnitude of the sky
S = 240      # Apparent diameter of object in arcseconds
Theta = 0.23 # Number of arcseconds per pixel on the CCD
Ndark = 0.5  # Number of electrons per pixel per second from dark current
Nread = 12   # Number of electrons per pixel from read noise
SN_100 = 100 # Signal-to-noise ratio of 100
SN_10 = 10   # Signal-to-noise ratio of 10
SN_5 = 5     # Signal-to-noise ratio of 5

# Number of electrons from the star
Nstar(x) = 5.5E7 * A * L * QE * (1 - E) / 2.5 ** x
# Number of electrons from the sky
Nsky = 5.5E7 * A * L * QE * (1 - E) / 2.5 ** Vsky
# Number of pixels collecting light from the object.
npix = 1.5 * (S / Theta) ** 2

# Quadratic variables
A(x) = Nstar(x) ** 2
B_100(x) = -(SN_100) ** 2 * (Nstar(x) + npix * (Nsky + Ndark))
C_100 = -(SN_100) ** 2 * (npix * Nread ** 2)
B_10(x) = -(SN_10) ** 2 * (Nstar(x) + npix * (Nsky + Ndark))
C_10 = -(SN_10) ** 2 * (npix * Nread ** 2)
B_5(x) = -(SN_5) ** 2 * (Nstar(x) + npix * (Nsky + Ndark))
C_5 = -(SN_5) ** 2 * (npix * Nread ** 2)

# Exposure time
t_100(x) = (-B_100(x) + (B_100(x) ** 2 - 4 * A(x) * C_100) ** (1/2)) / (2 * A(x))
t_10(x) = (-B_10(x) + (B_10(x) ** 2 - 4 * A(x) * C_10) ** (1/2)) / (2 * A(x))
t_5(x) = (-B_5(x) + (B_5(x) ** 2 - 4 * A(x) * C_5) ** (1/2)) / (2 * A(x))

set label "SNR of 100" at 9.5, 0.8 rotate by 43
set label "SNR of 10" at 12, 0.8 rotate by 43
set label "SNR of 5" at 12.75, 0.8 rotate by 43
set label "60 seconds" at 6, 85

plot t_100(x) w l, t_10(x) w l lt 1, t_5(x) w l lt 1, 60 w l lt 2