import numpy as np
from pathlib import Path
import pandas as pd
import xarray as xr
from scipy.interpolate import RectBivariateSpline
import logging
[docs]def blur(mu, variance):
""" The gaussian weights used in the convolution. """
return .001 / np.sqrt(2 * np.pi * variance) * np.exp(-.5 * mu ** 2 / variance)
[docs]class BlurredData(RectBivariateSpline):
""" Class defining the Solar Spectrum and useful functions. """
window_size = 3.0 # Must make an integer when multiplied by 2000. Should not be less then 5 sigma
num_window = int(window_size * 2000)
min_wavelength = 290.0
max_wavelength = 370.0
def __init__(self, spectra, sigma, measurement_sigma=0.0):
mu = np.linspace(-BlurredData.window_size, BlurredData.window_size, BlurredData.num_window + 1)
convolution_array = np.array([blur(mu, sigma ** 2.0 - measurement_sigma ** 2.0) for sigma in sigma])
blurred2d = np.array(
[np.convolve(row, spectra.interp(
wavelength=np.arange(
BlurredData.min_wavelength - BlurredData.window_size,
BlurredData.max_wavelength + BlurredData.window_size, .001)).values)
for row in convolution_array])
super().__init__(sigma, np.arange(BlurredData.min_wavelength, BlurredData.max_wavelength, .01),
blurred2d[:, BlurredData.num_window:-BlurredData.num_window:10])
[docs]class SAO2010(BlurredData):
""" Implements the Solar Spectrum Class with the SAO2010 data. """
def __init__(self, sigma=np.arange(.35, 0.61, 0.01)):
"""Constructor for SAO2010."""
sao = np.loadtxt(Path(__file__).parent.joinpath('SAO2010', 'sao2010.solref.converted'), skiprows=5)
solar = np.array([(nm, sun / nm ** 4) for (nm, sun) in zip(sao[:, 0], sao[:, 2]) if 285 <= nm < 375])
solar = xr.DataArray(solar[:, 1], dims={'wavelength'}, coords={'wavelength': solar[:, 0]})
super().__init__(solar, sigma)
logging.info('Blurred Solar Spectra (SAO2010) loaded.')
[docs]class O3BDM(BlurredData):
""" Implements the Brion–Daumont–Malicet (BDM) O3 Cross Section. """
def __init__(self, sigma=np.arange(.35, 0.61, 0.01), temperature=243):
"""Constructor for BDM."""
df = pd.concat(
[pd.read_csv(str(file), header=None, skiprows=1, index_col=0, names=['wavelength', int(file.name[11:14])],
sep=r'\s+') for file in Path(__file__).parent.joinpath('BDM').glob('*.dat')], axis=1)
df.index = df.index / 10 # Convert wavelength to nm
o3 = xr.DataArray.from_series(df[temperature])
fwhm = 0.02 # Measurement specific FWHM that needs to be factored in before additional blurring.
measurement_sigma = (fwhm / np.sqrt(8 * np.log(2)))
super().__init__(o3, sigma, measurement_sigma)
logging.info('Blurred O3 Cross-Section (BDM) loaded.')