# Source code for irfpy.util.solarspectra

''' Solar spectrum models at the Earth orbit

.. codeauthor:: Yoshifumi Futaana

.. autosummary::

BlackBody
BlackBodyModFuv
AM0

Solar spectrum models at the Earth orbit (1 AU).
See also :mod:irfpy.util.planck.
'''
import os
import logging
_logger = logging.getLogger(__name__)
from pkg_resources import resource_filename

import numpy as np
import scipy.interpolate

from irfpy.util import planck

[docs]class BlackBody:
r''' Black body assumption of the Sun.

Sun is proxied by a black body with tempereature 5778 K.

Using the planck's formula (see :class:irfpy.util.planck.BlackBodyMksa)
the spectrum at the Sun, i.e. power by area by wavelength [W / m2 / sr/ m],
is immediately known as :math:B(\lambda, T).
This is implemented in :mod:irfpy.util.planck.

The total spectral power is the integral over the Sun surface area and
the emission angle (in solid angle).
As isotropic emission is assumed, the solid angle integral can be done

.. math::

\int{\cos\theta d\Omega} = \int_0^{2\pi}d\phi\int_0^{\pi/2}\cos\theta\sin\theta d\theta = \pi

The :math:\cos\theta appears from the inner product of the flux and
Same calculation is done.

Anyway, the total spectral power, in the unit of [W / m], is

.. math::

F(\lambda, T) = \pi B(\lambda, T) \times 4\pi Rs^2

Here Rs is the radius of the Sun.

Consider a unit area at the Earth's orbit. The spectral power through
the area is, as the emission is isotropic, calculated by multiplying the
fraction of the area to the sphare surface at 1AU, namely,

.. math::

u(\lambda, T) = F(\lambda, T) \times \frac{1}{4\pi a^2}

where a is 1AU.  The unit of u is [W / m2 / m].

In summary, the spectrum power at the Earth orbit, u [W / m2 / m] is

.. math::

u(\lambda, T) = \pi \frac{Rs^2}{a^2} B(\lambda, T)

Substituting Rs=695500 km and a=149597871 km,

.. math::

u(\lambda, T) = 6.8687\times10^5 B(\lambda, T)

'''
def __init__(self, temperature=5778.):
self.sun = planck.SunProxy(temperature=temperature)
self.scale = 6.8687e-5  # Scaling from the spectrum at the sun (W / m2 / sr / m)
# to the spectrum at the Earth (W / m2 / m).
# sr is considered into the scale.
self.name = 'BlackBody(%gK)' % temperature
''' Name of the model '''

''' Return the spectral irradiance (W / m2 / m).

Return the spectral irradiance, in the unit of W / m2 / m.

:param wavelength: Wavelength in m.
:returns: Spectral irradiance in W / m2 / m

>>> sun = BlackBody(5800)    # 5800 K blackbody approx.
>>> print('%.3e' % irr500nm)
1.846e+09
'''
spectral_radiance_at_sun = self.sun.Blambda(wavelength)  # in W/m2 sr m

[docs]    def photon_flux(self, wavelength):
r''' Return the photon's differential flux (# / m2 / s / m).

The spectral irradiance, u [W / m2 / m = J / m2 / s / m],
can be converted to the photon flux.

A photon with wavelength, :math:\lambda, will have the energy
:math:E = hc/\lambda.
Thus, the number of photons through a unit area per unit time is
:math:F = u / E.

>>> sun = BlackBody(5800)
>>> fph = sun.photon_flux(500e-9)   # in /m2 /s /m
>>> fph = fph / 1e9   # in m2 /s /nm
>>> print('%.3e' % fph)
4.648e+18
'''
ene = planck.wavlen2energy(wavelength)
f = u / ene
return f

[docs]class BlackBodyModFuv(BlackBody):
''' Black body assumption with constant FUV (<~100 nm) modification
'''
def __init__(self, *args, **kwds):
BlackBody.__init__(self, *args, **kwds)
self.name = self.name.replace('BlackBody', 'BlackBody+constFUV')
''' Name of the model '''

'''

>>> sun = BlackBodyModFuv(5800)
1.846e+09
'''
irrarr = np.array(irr)
irrarr = np.where(np.logical_and(irrarr < 1e5, wavelength < 500e-9), 1e5, irrarr)

#        if irr < 1e5 and wavelength < 500e-9:
#            return 1e5
#        else:
#            return irr

return irrarr

[docs]class AM0:
''' Air mass 0 model (at extraterrestorial) under ASTM E-490.

>>> am0 = AM0()
>>> print('%.3e' % irr)
1.914e+09
>>> fph = am0.photon_flux(np.array([500e-9, 1000e-9]))
>>> print('%.3e %.3e' % (fph[0] / 1e9, fph[1] / 1e9))   # in m2 /nm /s
4.816e+18 3.765e+18
'''
logger = logging.getLogger(__name__ + '.AM0')

def __init__(self):
fn = resource_filename(__name__,
os.path.join('data', 'astm_e490_am0.txt.gz'))
self.logger.debug('Filename %s, exists: %s' % (fn,
os.path.exists(fn)))
self.logger.debug('Data size = %s' % str(dat.shape))
wavlen = dat[:, 0]  # Data file in nm
irrad = dat[:, 1]   # Data file in W/m2/nm

'''

self.name = 'ASTM/E490/AM0'

''' Return the irradiance with interface of MKSA
'''
wavlen_nm = np.array(wavelength) * 1e9   # m -> nm.

[docs]    def photon_flux(self, wavelength):
ene = planck.wavlen2energy(wavelength)
f = u / ene
return f

import unittest
import doctest

[docs]def doctests():
return unittest.TestSuite((
doctest.DocTestSuite(),
))
if __name__ == '__main__':
unittest.main(defaultTest='doctests')