~
Calculate flux, in photons and in watts, produced by a
star of the specified visual magnitude (Mv) and temperature.

Written by  Carl Witthoft
cgw@aoainc.com
Adaptive Optics Associates, Inc
Cambridge, MA
(disclaimer: The contents of this file are the author's work only and are
not supported or warranteed by the corporation.)

In these equations, the inputs are:
	l is the wavelength in microns
	t is the temperature in Kelvins
	mv is the visual magnitude

hphot(l,t,mv) calculates the output in photons/meter^2/sec/micron;
hwatt(l,t,mv) calculates the output in watts/meter^2/micron .
These are the irradiances at the top of Earth's atmosphere. To determine
the value at ground level requires applying spectral transmissivity.

For details on the derivation of the visual magnitude scale, see
for example the RCA E-O Handbook p. 65, or any astronomy text. Two sources
differ by one magnitude when their illuminance ratio is the
fifth root of 100, or ‰ 2.512.  The 0.555 factor normalizes Mv to 555 nm,
i.e., all curves with Mv=3, regardless of temperature, cross at 555 nm.
The value 1.4388e4 = hc/k,in micron-K , is part of the exponent of the
blackbody eq'n.   Finally,  3.579e-19 is the conversion from photons to
joules at 555nm.

~

hphot(l,t,mv) = 10^(-mv/2.512)*10^(-7.45)/3.5792793e-19*(0.555/l)^4 * 
                (exp(1.4388e4/.555/t)-1)/(exp(1.4388e4/l/t)-1)

hwatt(l,t,mv) = 10^(-mv/2.512)*10^(-7.45)*(0.555/l)^5 *
                (exp(1.4388e4/.555/t)-1)/(exp(1.4388e4/l/t)-1)

--where

l[i] = i/20 dim[40]
t:= 6500:
mv:= 3:

hphot(l,t,mv):{0.000,77594.708,24541006.251,310710975.540,1164392285.067,2457359369.215,3810343953.443,4934780898.373,5716877600.202,6163014628.780,...}

l:{0.050,0.100,0.150,0.200,0.250,0.300,0.350,0.400,0.450,0.500,...}
l[20]:1.000;


Marker = {3,5}
slabel  "mv=3 and T = 9400"
slabel "mv=3 and T = 6500"
Title ="photons (upper) and watts (lower) vs T, MV, lambda (microns)"
 Xmin:= 0:
Tracecolor := 1:
-- Ylabel := "Photons/m^2/sec/µ":
plotline {l, hphot(l,t,mv)}
Tracecolor :=3:
plotline{l, hphot (l, 9400,mv)}


newaxis ; Ystrips := 2:
-- Ylabel := "Watts/m^2/µ":

Tracecolor :=1:
plotline {l,hwatt(l,t,mv)}
Tracecolor :=3:
plotline {l, hwatt(l,9400,mv)}