Commit 6faebb1f authored by Médéric Boquien's avatar Médéric Boquien

Merge branch 'release/v0.12.0'

parents c0c15ec8 ff3e6298

Too many changes to show.

To preserve performance only 99 of 99+ files are displayed.
pcigale author list
pcigale authors list
====================
This document lists the various authors who wrote pcigale code with their
current email address and affiliation.
This document lists alphabetically the various authors who wrote the pcigale
code with their current email address and affiliation.
* Yannick Roehlly <yannick@iaora.eu>
* Laure Ciesla <ciesla@physics.uoc.gr>,
University of Crete, Department of Physics
* Médéric Boquien <mboquien@ast.cam.ac.uk>,
University of Cambridge, Institute of Astronomy
* Médéric Boquien <mederic.boquien@uantof.cl>,
Universidad de Antofagasta, Chile
* Denis Burgarella <denis.burgarella@lam.fr>,
Laboratoire d'Astrophysique de Marseille
Laboratoire d'Astrophysique de Marseille, France
* Laure Ciesla <ciesla@lam.fr>,
Laboratoire d'Astrophysique de Marseille, France
* Yannick Roehlly <yannick@iaora.eu>
# Change Log
## 0.12.0 (2018-02-19)
### Added
- Provide the possibility not to store a given module in cache. This can be useful on computers with a limited amount of memory. The downside is that when not caching the model generation will be slower. (Médéric Boquien)
- An option `redshift\_decimals` is now provided in `pdf\_analysis` to indicate the number of decimals to round the observed redshifts to compute the grid of models. By default the model redshifts are rounded to two decimals but this can be insufficient at low z and/or when using narrow-band filters for instance. This only applies to the grid. The physical properties are still computed for the redshift at full precision. (Médéric Boquien)
- Bands with negative fluxes are now considered valid and are fitted as any other band. (Médéric Boquien)
- Allow the models to be computed by blocks in `savefluxes`. This can be useful when computing a very large grid and/or to split the results file into various smaller files as large files can be difficult to handle. The number of blocks is set with the `blocks` parameters in the pcigale.ini. (Médéric Boquien)
- Allow the observations to be analysed by blocks of models in `pdf\_analysis`. This is useful when computing a very large grid of models that would not fit in memory. The number of blocks is set with the `blocks` parameters in the pcigale.ini. (Médéric Boquien)
- The integrated stellar luminosity is now provided as `stellar.lum`. (Médéric Boquien)
- The high resolution BC03 models have been added. They can be activated when building the database by adding `--bc03res=hr` to the build command. In that case the low resolution models are not built. (Médéric Boquien)
- Dust templates generated with THEMIS (Jones et al. 2017) have been contributed by the DustPedia team (Davis et al. 2017). Special acknowledgement to Angelos Nersesian and Frédéric Galliano for creating the dust templates and writing the code. (Dustpedia team)
- The Herschel SPIRE filters for extended sources have been added. (Médéric Boquien)
### Changed
- Make the timestamp more readable when moving the out/ directory. (Médéric Boquien)
- To accommodate the analysis of the observations by blocks, all models are now included in the estimation of the physical properties and no cut in chi² is done anymore. (Médéric Boquien)
- To accommodate the analysis of the observations by blocks, the `save_pdf` option has been eliminated. To plot PDF one needs to set `save_chi2` to True and then run `pcigale-plots pdf`. (Médéric Boquien)
- In order to capture rapid evolutionary phases, we assume that in a given period of 1 Myr, 10 small episodes of star formation occurred every 0.1 Myr, rather than one episode every 1 Myr.
- When computing the attenuation curve, the bump is now added so that its relative strength does not depend on δ. (Médéric Boquien, issue reported by Samir Salim)
- Fν was computed by calculating Fλ and then converting to Fν, which led to typical differences in fluxes of typically 1-2% and a bit more for a handful of pathological filters. Now Fν is computed directly and a bit faster. (Médéric Boquien, issue reported by Yannick Roehlly and Wouter Dobbels)
### Fixed
- Corrected a typo that prevented `restframe\_parameters` from being listed among the available modules. (Médéric Boquien)
- The filters in the residual plot of `pcigale-plots sed` are now drawn in order of increasing wavelength so that the line joining all the filters does not make loops. (Médéric Boquien)
- In the absence of a nebular component `restframe\_parameters` would crash when attempting to compute the equivalent widths of the lines listed in `EW_lines`. Now they are simply ignored. (Médéric Boquien)
- The luminosity spectrum of the best fit was saved assuming the distance corresponding to the redshift rounded to two decimals. This was an issue in particular at very low redshift as a difference of 0.005 in redshift can translate to a large difference on the luminosity distance. Now the exact luminosity distance of the object is used to compute the spectrum luminosity. (Médéric Boquien, reported by Jorge Melnick)
- When using the `parameters\_file` option, the indices of the models now correspond to the line number of the input file. (Médéric Boquien)
- When using the `parameters\_file` option, the list of modules is read from `sed\_modules` rather than being inferred from the input file. (Médéric Boquien)
- The computation of the upper limits would only work for the first few models, reverting back to regular fits for the others. (Médéric Boquien)
- A more explicit message is now given when the flux table cannot be read properly. (Médéric Boquien)
- Make sure that we do not try to fit data that have an error bar of 0 mJy. (Médéric Boquien)
- An erroneous warning was displayed when using the `restframe\_parameters` module. (Médéric Boquien)
- The formula from Sawicki et al. (2012) used to compute the χ² in the presence of upper limits was not correct. This led the χ² to depend directly on the absolute value of the upper limit. The formula has been rederived and corrected. (Médéric Boquien & Denis Burgarella)
- For some reason the wavelengths of the SCUBA 450 μm filter were a factor 10 too small. (Médéric Boquien)
- Ensure that the computation of the continuum level is correct when determining the equivalent width, in particular when the line width is very narrow. (Médéric Boquien)
- When using different line widths during a single run, ensure that the fluxes and other quantities are always computed correctly. (Médéric Boquien, special thanks to Genoveva Micheva)
- Compute Dn4000 more rigorously by integrating properly over Fν. (Médéric Boquien)
### Optimised
- The cache architecture has been simplified, making it somewhat faster. It speeds up the model generation by ~1%. (Médéric Boquien)
## 0.11.0 (2017-02-10)
### Added
- The stellar mass-weighted age is now provided. This is a much more usual measure of the age than the age of the oldest star. This is accessible with the `stellar.age_m_star` keyword in the `bc03` module with with the `stellar.age_mass` keyword in the `m2005` module. (Médéric Boquien)
- The nebular models have been expanded from log U=-3 to log U=-4. (Médéric Boquien & Akio Inoue)
- The nebular models are now sampled in steps of 0.1 dex in log U rather than 1.0 dex steps. (Médéric Boquien & Akio Inoue)
- A new set of filters from GAZPAR has been added. The pattern of the filter name is "telescope.instrument.filter", e.g. "hst.wfc3.F160W". If the telescope has one instrument, it is skipped, e.g. "galex.FUV". For now the original set of filters is still provided. (Médéric Boquien & Olivier Ilbert)
- A brand new module `restframe\_param` has been added to compute rest frame parameters: UV slope β (Calzetti et al. 1994), Dn4000 (Balogh et al. 1999), IRX, emission lines equivalent widths at any wavelength, luminosity in any filter, colour in any pair of filters. This module has to be inserted right before the redshifting module. (Médéric Boquien)
- A brand new module `restframe\_parameters` has been added to compute rest frame parameters: UV slope β (Calzetti et al. 1994), Dn4000 (Balogh et al. 1999), IRX, emission lines equivalent widths at any wavelength, luminosity in any filter, colour in any pair of filters. This module has to be inserted right before the redshifting module. (Médéric Boquien)
### Changed
- We do not output the break strength from the `bc03` module anymore as these were not computed properly. (Médéric Boquien)
- The new `restframe\_param` module replaces the unofficial `param` module, which has now been trimmed to only compute fluxes in the observed frame as the rest of its functionalities have been transferred to the much more efficient `restframe\_param` module. To reflect this, it has been renamed `fluxes`. (Médéric Boquien)
- The new `restframe\_parameters` module replaces the unofficial `param` module, which has now been trimmed to only compute fluxes in the observed frame as the rest of its functionalities have been transferred to the much more efficient `restframe\_parameters` module. To reflect this, it has been renamed `fluxes`. (Médéric Boquien)
- We now make use of new features available in Python 3.5. Previous versions are henceforth unsupported. However Python 3.6 or later is recommended for better performance. (Médéric Boquien)
### Fixed
......
......@@ -22,7 +22,7 @@ import scipy.constants as cst
from astropy.table import Table
from pcigale.data import (Database, Filter, M2005, BC03, Fritz2006,
Dale2014, DL2007, DL2014, NebularLines,
NebularContinuum, Schreiber2016)
NebularContinuum, Schreiber2016, THEMIS)
def read_bc03_ssp(filename):
......@@ -165,13 +165,13 @@ def build_filters(base):
new_filter = Filter(filter_name, filter_description, filter_table)
# We normalise the filter and compute the effective wavelength.
# If the filter is a pseudo-filter used to compute line fluxes, it
# should not be normalised.
# We normalise the filter and compute the pivot wavelength. If the
# filter is a pseudo-filter used to compute line fluxes, it should not
# be normalised.
if not filter_name.startswith('PSEUDO'):
new_filter.normalise()
else:
new_filter.effective_wavelength = np.mean(
new_filter.pivot_wavelength = np.mean(
filter_table[0][filter_table[1] > 0]
)
filters.append(new_filter)
......@@ -211,13 +211,13 @@ def build_filters_gazpar(base):
new_filter = Filter(filter_name, filter_desc, filter_table)
# We normalise the filter and compute the effective wavelength.
# If the filter is a pseudo-filter used to compute line fluxes, it
# should not be normalised.
# We normalise the filter and compute the pivot wavelength. If the
# filter is a pseudo-filter used to compute line fluxes, it should not
# be normalised.
if not filter_name.startswith('PSEUDO'):
new_filter.normalise()
else:
new_filter.effective_wavelength = np.mean(
new_filter.pivot_wavelength = np.mean(
filter_table[0][filter_table[1] > 0]
)
filters.append(new_filter)
......@@ -228,7 +228,8 @@ def build_m2005(base):
m2005_dir = os.path.join(os.path.dirname(__file__), 'maraston2005/')
# Age grid (1 Myr to 13.7 Gyr with 1 Myr step)
age_grid = np.arange(1, 13701)
time_grid = np.arange(1, 13701)
fine_time_grid = np.linspace(0.1, 13700, 137000)
# Transpose the table to have access to each value vector on the first
# axis
......@@ -258,50 +259,57 @@ def build_m2005(base):
# populations.
mass_table = mass_table[1:7, mass_table[0] == metallicity]
# Interpolate the mass table over the new age grid. We multiply per
# 1000 because the time in Maraston files is given in Gyr.
mass_table = interpolate.interp1d(mass_table[0] * 1000,
mass_table)(age_grid)
# Regrid the SSP data to the evenly spaced time grid. In doing so we
# assume 10 bursts every 0.1 Myr over a period of 1 Myr in order to
# capture short evolutionary phases.
# The time grid starts after 0.1 Myr, so we assume the value is the same
# as the first actual time step.
mass_table = interpolate.interp1d(mass_table[0] * 1e3, mass_table[1:],
assume_sorted=True)(fine_time_grid)
mass_table = np.mean(mass_table.reshape(5, -1, 10), axis=-1)
# Remove the age column from the mass table
mass_table = np.delete(mass_table, 0, 0)
# Extract the age and convert from Gyr to Myr
ssp_time = np.unique(spec_table[0]) * 1e3
spec_table = spec_table[1:]
# Remove the metallicity column from the spec table
spec_table = np.delete(spec_table, 1, 0)
# Convert the wavelength from Å to nm
spec_table[1] = spec_table[1] * 0.1
# For all ages, the lambda grid is the same.
lambda_grid = np.unique(spec_table[1])
# Creation of the age vs lambda flux table
tmp_list = []
for wavelength in lambda_grid:
[age_grid_orig, lambda_grid_orig, flux_orig] = \
spec_table[:, spec_table[1, :] == wavelength]
flux_orig = flux_orig * 10 * 1.e-7 # From erg/s^-1/Å to W/nm
age_grid_orig *= 1000 # Gyr to Myr
flux_regrid = interpolate.interp1d(age_grid_orig,
flux_orig)(age_grid)
tmp_list.append(flux_regrid)
flux_age = np.array(tmp_list)
spec_table = spec_table[1:]
# Extract the wavelength and convert from Å to nm
ssp_wave = spec_table[0][:1221] * 0.1
spec_table = spec_table[1:]
# Extra the fluxes and convert from erg/s/Å to W/nm
ssp_lumin = spec_table[0].reshape(ssp_time.size, ssp_wave.size).T
ssp_lumin *= 10 * 1e-7
# We have to do the interpolation-averaging in several blocks as it is
# a bit RAM intensive
ssp_lumin_interp = np.empty((ssp_wave.size, time_grid.size))
for i in range(0, ssp_wave.size, 100):
fill_value = (ssp_lumin[i:i+100, 0], ssp_lumin[i:i+100, -1])
ssp_interp = interpolate.interp1d(ssp_time, ssp_lumin[i:i+100, :],
fill_value=fill_value,
bounds_error=False,
assume_sorted=True)(fine_time_grid)
ssp_interp = ssp_interp.reshape(ssp_interp.shape[0], -1, 10)
ssp_lumin_interp[i:i+100, :] = np.mean(ssp_interp, axis=-1)
# To avoid the creation of waves when interpolating, we refine the grid
# beyond 10 μm following a log scale in wavelength. The interpolation
# is also done in log space as the spectrum is power-law-like
lambda_grid_resamp = np.around(np.logspace(np.log10(10000),
ssp_wave_resamp = np.around(np.logspace(np.log10(10000),
np.log10(160000), 50))
argmin = np.argmin(10000.-lambda_grid > 0)-1
flux_age_resamp = 10.**interpolate.interp1d(
np.log10(lambda_grid[argmin:]),
np.log10(flux_age[argmin:, :]),
argmin = np.argmin(10000.-ssp_wave > 0)-1
ssp_lumin_resamp = 10.**interpolate.interp1d(
np.log10(ssp_wave[argmin:]),
np.log10(ssp_lumin_interp[argmin:, :]),
assume_sorted=True,
axis=0)(np.log10(lambda_grid_resamp))
axis=0)(np.log10(ssp_wave_resamp))
lambda_grid = np.hstack([lambda_grid[:argmin+1], lambda_grid_resamp])
flux_age = np.vstack([flux_age[:argmin+1, :], flux_age_resamp])
ssp_wave = np.hstack([ssp_wave[:argmin+1], ssp_wave_resamp])
ssp_lumin = np.vstack([ssp_lumin_interp[:argmin+1, :],
ssp_lumin_resamp])
# Use Z value for metallicity, not log([Z/H])
metallicity = {-1.35: 0.001,
......@@ -309,15 +317,16 @@ def build_m2005(base):
0.0: 0.02,
0.35: 0.04}[metallicity]
base.add_m2005(M2005(imf, metallicity, age_grid, lambda_grid,
mass_table, flux_age))
base.add_m2005(M2005(imf, metallicity, time_grid, ssp_wave,
mass_table, ssp_lumin))
def build_bc2003(base):
bc03_dir = os.path.join(os.path.dirname(__file__), 'bc03//')
def build_bc2003(base, res):
bc03_dir = os.path.join(os.path.dirname(__file__), 'bc03/')
# Time grid (1 Myr to 14 Gyr with 1 Myr step)
time_grid = np.arange(1, 14000)
fine_time_grid = np.linspace(0.1, 13999, 139990)
# Metallicities associated to each key
metallicity = {
......@@ -330,12 +339,14 @@ def build_bc2003(base):
}
for key, imf in itertools.product(metallicity, ["salp", "chab"]):
base_filename = bc03_dir + "bc2003_lr_" + key + "_" + imf + "_ssp"
ssp_filename = base_filename + ".ised_ASCII"
color3_filename = base_filename + ".3color"
color4_filename = base_filename + ".4color"
ssp_filename = "{}bc2003_{}_{}_{}_ssp.ised_ASCII".format(bc03_dir, res,
key, imf)
color3_filename = "{}bc2003_lr_{}_{}_ssp.3color".format(bc03_dir, key,
imf)
color4_filename = "{}bc2003_lr_{}_{}_ssp.4color".format(bc03_dir, key,
imf)
print("Importing %s..." % base_filename)
print("Importing {}...".format(ssp_filename))
# Read the desired information from the color files
color_table = []
......@@ -344,19 +355,34 @@ def build_bc2003(base):
color_table.append(color4_table[6]) # Mstar
color_table.append(color4_table[7]) # Mgas
color_table.append(10 ** color3_table[5]) # NLy
color_table.append(color3_table[1]) # B4000
color_table.append(color3_table[2]) # B4_VN
color_table.append(color3_table[3]) # B4_SDSS
color_table.append(color3_table[4]) # B(912)
color_table = np.array(color_table)
ssp_time, ssp_wave, ssp_lumin = read_bc03_ssp(ssp_filename)
# Regrid the SSP data to the evenly spaced time grid.
color_table = interpolate.interp1d(ssp_time, color_table)(time_grid)
ssp_lumin = interpolate.interp1d(ssp_time,
ssp_lumin)(time_grid)
# Regrid the SSP data to the evenly spaced time grid. In doing so we
# assume 10 bursts every 0.1 Myr over a period of 1 Myr in order to
# capture short evolutionary phases.
# The time grid starts after 0.1 Myr, so we assume the value is the same
# as the first actual time step.
fill_value = (color_table[:, 0], color_table[:, -1])
color_table = interpolate.interp1d(ssp_time, color_table,
fill_value=fill_value,
bounds_error=False,
assume_sorted=True)(fine_time_grid)
color_table = np.mean(color_table.reshape(3, -1, 10), axis=-1)
# We have to do the interpolation-averaging in several blocks as it is
# a bit RAM intensive
ssp_lumin_interp = np.empty((ssp_wave.size, time_grid.size))
for i in range(0, ssp_wave.size, 100):
fill_value = (ssp_lumin[i:i+100, 0], ssp_lumin[i:i+100, -1])
ssp_interp = interpolate.interp1d(ssp_time, ssp_lumin[i:i+100, :],
fill_value=fill_value,
bounds_error=False,
assume_sorted=True)(fine_time_grid)
ssp_interp = ssp_interp.reshape(ssp_interp.shape[0], -1, 10)
ssp_lumin_interp[i:i+100, :] = np.mean(ssp_interp, axis=-1)
# To avoid the creation of waves when interpolating, we refine the grid
# beyond 10 μm following a log scale in wavelength. The interpolation
......@@ -366,12 +392,13 @@ def build_bc2003(base):
argmin = np.argmin(10000.-ssp_wave > 0)-1
ssp_lumin_resamp = 10.**interpolate.interp1d(
np.log10(ssp_wave[argmin:]),
np.log10(ssp_lumin[argmin:, :]),
np.log10(ssp_lumin_interp[argmin:, :]),
assume_sorted=True,
axis=0)(np.log10(ssp_wave_resamp))
ssp_wave = np.hstack([ssp_wave[:argmin+1], ssp_wave_resamp])
ssp_lumin = np.vstack([ssp_lumin[:argmin+1, :], ssp_lumin_resamp])
ssp_lumin = np.vstack([ssp_lumin_interp[:argmin+1, :],
ssp_lumin_resamp])
base.add_bc03(BC03(
imf,
......@@ -742,7 +769,77 @@ def build_schreiber2016(base):
base.add_schreiber2016(models)
def build_base():
def build_themis(base):
models = []
themis_dir = os.path.join(os.path.dirname(__file__), 'themis/')
# Mass fraction of hydrocarbon solids i.e., a-C(:H) smaller than 1.5 nm,
# also known as HAC
qhac = {"000": 0.02, "010": 0.06, "020": 0.10, "030": 0.14, "040": 0.17,
"050": 0.20, "060": 0.24, "070": 0.28, "080": 0.32, "090": 0.36,
"100": 0.40}
uminimum = ["0.100", "0.120", "0.150", "0.170", "0.200", "0.250", "0.300",
"0.350", "0.400", "0.500", "0.600", "0.700", "0.800", "1.000",
"1.200", "1.500", "1.700", "2.000", "2.500", "3.000", "3.500",
"4.000", "5.000", "6.000", "7.000", "8.000", "10.00", "12.00",
"15.00", "17.00", "20.00", "25.00", "30.00", "35.00", "40.00",
"50.00", "80.00"]
alpha = ["1.0", "1.1", "1.2", "1.3", "1.4", "1.5", "1.6", "1.7", "1.8",
"1.9", "2.0", "2.1", "2.2", "2.3", "2.4", "2.5", "2.6", "2.7",
"2.8", "2.9", "3.0"]
# Mdust/MH used to retrieve the dust mass as models as given per atom of H
MdMH = {"000": 7.4e-3, "010": 7.4e-3, "020": 7.4e-3, "030": 7.4e-3,
"040": 7.4e-3, "050": 7.4e-3, "060": 7.4e-3, "070": 7.4e-3,
"080": 7.4e-3, "090": 7.4e-3, "100": 7.4e-3}
# Here we obtain the wavelength beforehand to avoid reading it each time.
datafile = open(themis_dir + "U{}_{}_MW3.1_{}/spec_1.0.dat"
.format(uminimum[0], uminimum[0], "000"))
data = "".join(datafile.readlines()[-576:])
datafile.close()
wave = np.genfromtxt(io.BytesIO(data.encode()), usecols=(0))
# We convert wavelengths from μm to nm
wave *= 1000.
# Conversion factor from Jy cm² sr¯¹ H¯¹ to W nm¯¹ (kg of H)¯¹
conv = 4. * np.pi * 1e-30 / (cst.m_p+cst.m_e) * cst.c / (wave*wave) * 1e9
for model in sorted(qhac.keys()):
for umin in uminimum:
filename = (themis_dir + "U{}_{}_MW3.1_{}/spec_1.0.dat"
.format(umin, umin, model))
print("Importing {}...".format(filename))
with open(filename) as datafile:
data = "".join(datafile.readlines()[-576:])
lumin = np.genfromtxt(io.BytesIO(data.encode()), usecols=(2))
# Conversion from Jy cm² sr¯¹ H¯¹to W nm¯¹ (kg of dust)¯¹
lumin *= conv / MdMH[model]
models.append(THEMIS(qhac[model], umin, umin, 1.0, wave, lumin))
for al in alpha:
filename = (themis_dir + "U{}_1e7_MW3.1_{}/spec_{}.dat"
.format(umin, model, al))
print("Importing {}...".format(filename))
with open(filename) as datafile:
data = "".join(datafile.readlines()[-576:])
lumin = np.genfromtxt(io.BytesIO(data.encode()), usecols=(2))
# Conversion from Jy cm² sr¯¹ H¯¹to W nm¯¹ (kg of dust)¯¹
lumin *= conv/MdMH[model]
models.append(THEMIS(qhac[model], umin, 1e7, al, wave, lumin))
base.add_themis(models)
def build_base(bc03res='lr'):
base = Database(writable=True)
base.upgrade_base()
......@@ -759,7 +856,7 @@ def build_base():
print('#' * 78)
print("3- Importing Bruzual and Charlot 2003 SSP\n")
build_bc2003(base)
build_bc2003(base, bc03res)
print("\nDONE\n")
print('#' * 78)
......@@ -792,7 +889,12 @@ def build_base():
build_schreiber2016(base)
print("\nDONE\n")
print('#' * 78)
print("10- Importing Jones et al (2017) models)\n")
build_themis(base)
print("\nDONE\n")
print('#' * 78)
base.session.close_all()
......
This source diff could not be displayed because it is too large. You can view the blob instead.
This source diff could not be displayed because it is too large. You can view the blob instead.
This source diff could not be displayed because it is too large. You can view the blob instead.
This source diff could not be displayed because it is too large. You can view the blob instead.
This source diff could not be displayed because it is too large. You can view the blob instead.
This source diff could not be displayed because it is too large. You can view the blob instead.
This source diff could not be displayed because it is too large. You can view the blob instead.
This source diff could not be displayed because it is too large. You can view the blob instead.
This source diff could not be displayed because it is too large. You can view the blob instead.
This source diff could not be displayed because it is too large. You can view the blob instead.
This source diff could not be displayed because it is too large. You can view the blob instead.
This source diff could not be displayed because it is too large. You can view the blob instead.
# PLW_ext
# energy
# SPIRE Herschel extended sources
3438789.5 1.7189e-5
3450655.6 1.9435e-5
3462603.9 2.4090e-5
3473428.3 3.8782e-5
3485535.0 5.1049e-5
3497726.5 6.5301e-5
3510003.5 1.0005e-4
3521126.8 1.5676e-4
3533568.9 2.1130e-4
3546099.3 2.3606e-4
3558718.9 2.9043e-4
3571428.6 3.6418e-4
3582945.2 3.5368e-4
3595828.8 2.9929e-4
3608805.5 3.2039e-4
3621876.1 3.7160e-4
3635041.8 4.2752e-4
3646973.0 5.0005e-4
3660322.1 4.6930e-4
3673769.3 4.6922e-4
3687315.6 4.7055e-4
3699593.0 5.2972e-4
3713330.9 5.1960e-4
3727171.1 5.7360e-4
3741114.9 6.2072e-4
3755163.3 7.5872e-4
3767897.5 1.0568e-3
3782148.3 1.2394e-3
3796507.2 1.5078e-3
3810975.6 1.7553e-3
3825554.7 2.3163e-3
3838771.6 2.8562e-3
3853564.5 3.5237e-3
3868472.0 4.5131e-3
3883495.1 5.7435e-3
3897116.1 7.2597e-3
3912363.1 9.5634e-3
3927729.8 1.2947e-2
3943217.7 1.7169e-2
3958828.2 2.2355e-2
3972983.7 2.8932e-2
3988831.3 3.7777e-2
4004805.8 4.8490e-2
4020908.7 6.2011e-2
4037141.7 8.0548e-2
4051863.9 1.0281e-1
4068348.3 1.3181e-1
4084967.3 1.7400e-1
4101722.7 2.2578e-1
4116920.5 2.8562e-1
4133939.6 3.3612e-1
4151100.0 3.7965e-1
4168403.5 4.2041e-1
4185851.8 4.3663e-1
4201680.7 4.4928e-1
4219409.3 4.6496e-1
4237288.1 4.8225e-1
4255319.1 4.9517e-1
4271678.8 5.1090e-1
4290004.3 5.2204e-1
4308487.7 5.3735e-1
4327131.1 5.4022e-1
4345936.5 5.5033e-1
4363001.7 5.3978e-1
4382120.9 5.4205e-1
4401408.5 5.5477e-1
4420866.5 5.6928e-1
4440497.3 5.7548e-1
4458314.8 5.7328e-1
4478280.3 5.7881e-1
4498425.6 5.5404e-1
4518752.8 5.4991e-1
4537205.1 5.5117e-1
4557885.1 5.4779e-1
4578754.6 5.4823e-1
4599816.0 5.7042e-1
4621072.1 5.9189e-1
4640371.2 5.9920e-1
4662004.7 6.2448e-1
4683840.7 6.1801e-1
4705882.4 6.1402e-1
4728132.4 5.9609e-1
4748338.1 5.9431e-1
4770992.4 5.7407e-1
4793863.9 5.7422e-1
4816955.7 5.6672e-1
4837929.4 5.8706e-1
4861448.7 6.0866e-1
4885197.9 6.2440e-1
4909180.2 6.3158e-1
4933399.1 6.2554e-1
4955401.4 6.2867e-1
4980079.7 6.2821e-1
5005005.0 6.3535e-1
5030181.1 6.3911e-1
5055611.7 6.3523e-1
5078720.2 6.2544e-1
5104645.2 6.4178e-1
5130836.3 6.5917e-1
5157297.6 6.4811e-1
5181347.2 6.6172e-1
5208333.3 6.7057e-1
5235602.1 6.5462e-1
5263157.9 6.6403e-1
5291005.3 6.9016e-1
5316321.1 7.3207e-1
5344735.4 7.2943e-1
5373455.1 7.4458e-1
5402485.1 7.6016e-1
5431830.5 7.5129e-1
5458515.3 7.5642e-1
5488474.2 7.3704e-1
5518763.8 7.3503e-1
5549389.6 7.5848e-1
5577244.8 7.4298e-1
5608525.0 7.7443e-1
5640157.9 8.3959e-1
5672149.7 8.5941e-1
5704506.6 8.8025e-1
5733945.0 8.8184e-1
5767012.7 8.7001e-1
5800464.0 8.5624e-1
5834305.7 8.7827e-1
5865102.6 8.7855e-1
5899705.0 9.2660e-1
5934718.1 9.7338e-1
5970149.3 9.5633e-1
6006006.0 9.7146e-1
6038647.3 9.9422e-1
6075334.1 1.0000e+0
6112469.4 4.3820e-1
6150061.5 2.4551e-1
6188118.8 1.6527e-1
6222775.4 1.2126e-1
6261740.8 9.5505e-2
6301197.2 7.4838e-2
6341154.1 6.0805e-2
6377551.0 5.0165e-2
6418485.2 4.2153e-2
6459948.3 3.4923e-2
6501950.6 2.8790e-2
6544502.6 2.4413e-2
6583278.5 2.2513e-2
6626905.2 1.9730e-2
6671114.1 1.7073e-2
6715916.7 1.6113e-2
6761325.2 1.4671e-2