Skip to content

Commit

Permalink
Merge pull request #115 from karllark/black_all
Browse files Browse the repository at this point in the history
black all code
  • Loading branch information
karllark authored Nov 28, 2023
2 parents bb09adb + 6674e6d commit 8f5019f
Show file tree
Hide file tree
Showing 14 changed files with 96 additions and 45 deletions.
10 changes: 7 additions & 3 deletions measure_extinction/conftest.py
Original file line number Diff line number Diff line change
Expand Up @@ -8,13 +8,15 @@
from astropy.version import version as astropy_version

# For Astropy 3.0 and later, we can use the standalone pytest plugin
if astropy_version < '3.0':
if astropy_version < "3.0":
from astropy.tests.pytest_plugins import * # noqa

del pytest_report_header
ASTROPY_HEADER = True
else:
try:
from pytest_astropy_header.display import PYTEST_HEADER_MODULES, TESTED_VERSIONS

ASTROPY_HEADER = True
except ImportError:
ASTROPY_HEADER = False
Expand All @@ -28,13 +30,15 @@ def pytest_configure(config):

# Customize the following lines to add/remove entries from the list of
# packages for which version numbers are displayed when running the tests.
PYTEST_HEADER_MODULES.pop('Pandas', None)
PYTEST_HEADER_MODULES['scikit-image'] = 'skimage'
PYTEST_HEADER_MODULES.pop("Pandas", None)
PYTEST_HEADER_MODULES["scikit-image"] = "skimage"

from . import __version__

packagename = os.path.basename(os.path.dirname(__file__))
TESTED_VERSIONS[packagename] = __version__


# Uncomment the last two lines in this block to treat all DeprecationWarnings as
# exceptions. For Astropy v2.0 or later, there are 2 additional keywords,
# as follows (although default should work for most cases).
Expand Down
32 changes: 24 additions & 8 deletions measure_extinction/merge_obsspec.py
Original file line number Diff line number Diff line change
Expand Up @@ -42,7 +42,7 @@ def _wavegrid(resolution, wave_range):
np.log10(wave_range[1]) - delta_wave_log,
delta_wave_log,
)
full_wave_min = 10 ** wave_log10
full_wave_min = 10**wave_log10
full_wave_max = 10 ** (wave_log10 + delta_wave_log)

full_wave = (full_wave_min + full_wave_max) / 2.0
Expand Down Expand Up @@ -174,7 +174,9 @@ def merge_stis_obsspec(obstables, waveregion="UV", output_resolution=1000):
(indxs,) = np.where(full_npts > 0)
if len(indxs) > 0:
full_flux[indxs] /= full_unc[indxs]
full_unc[indxs] = np.sqrt(1.0 / full_unc[indxs]) # * 1e-10 # put back factor for overflow errors
full_unc[indxs] = np.sqrt(
1.0 / full_unc[indxs]
) # * 1e-10 # put back factor for overflow errors

otable = Table()
otable["WAVELENGTH"] = Column(full_wave, unit=u.angstrom)
Expand Down Expand Up @@ -218,7 +220,7 @@ def merge_spex_obsspec(obstable, mask=[], output_resolution=2000):
# take out data points with NaN fluxes
npts[np.isnan(fluxes)] = 0
# quadratically add 1 percent uncertainty to account for unknown uncertainties
uncs = np.sqrt(uncs ** 2 + (0.01 * fluxes) ** 2)
uncs = np.sqrt(uncs**2 + (0.01 * fluxes) ** 2)
# take out data points with low SNR
npts[np.less(fluxes / uncs, 10, where=~np.isnan(fluxes / uncs))] = 0
# take out wavelength regions affected by the atmosphere
Expand Down Expand Up @@ -307,8 +309,16 @@ def merge_gen_obsspec(obstables, wave_range, output_resolution=100):
full_npts = np.zeros((n_waves), dtype=int)
for ctable in obstables:
cwaves = ctable["WAVELENGTH"].to(u.angstrom).value
cfluxes = ctable["FLUX"].to(fluxunit, equivalencies=u.spectral_density(ctable["WAVELENGTH"])).value
cuncs = ctable["ERROR"].to(fluxunit, equivalencies=u.spectral_density(ctable["WAVELENGTH"])).value
cfluxes = (
ctable["FLUX"]
.to(fluxunit, equivalencies=u.spectral_density(ctable["WAVELENGTH"]))
.value
)
cuncs = (
ctable["ERROR"]
.to(fluxunit, equivalencies=u.spectral_density(ctable["WAVELENGTH"]))
.value
)
cnpts = ctable["NPTS"].value
for k in range(n_waves):
(indxs,) = np.where(
Expand Down Expand Up @@ -356,7 +366,9 @@ def merge_irs_obsspec(obstables, output_resolution=150):
merged spectra
"""
wave_range = [5.0, 40.0] * u.micron
otable = merge_gen_obsspec(obstables, wave_range, output_resolution=output_resolution)
otable = merge_gen_obsspec(
obstables, wave_range, output_resolution=output_resolution
)
return otable


Expand All @@ -381,7 +393,9 @@ def merge_nircam_ss_obsspec(obstables, output_resolution=1600):
merged spectra
"""
wave_range = [2.35, 5.55] * u.micron
otable = merge_gen_obsspec(obstables, wave_range, output_resolution=output_resolution)
otable = merge_gen_obsspec(
obstables, wave_range, output_resolution=output_resolution
)
return otable


Expand All @@ -406,5 +420,7 @@ def merge_miri_ifu_obsspec(obstables, output_resolution=3000):
merged spectra
"""
wave_range = [4.8, 29.0] * u.micron
otable = merge_gen_obsspec(obstables, wave_range, output_resolution=output_resolution)
otable = merge_gen_obsspec(
obstables, wave_range, output_resolution=output_resolution
)
return otable
17 changes: 11 additions & 6 deletions measure_extinction/modeldata.py
Original file line number Diff line number Diff line change
Expand Up @@ -312,11 +312,11 @@ def dust_extinguished_sed(self, params, sed, velocity=0.0):
-0.426 + 1.0044 * Rv,
-0.050 + 1.0016 * Rv,
0.701 + 1.0016 * Rv,
1.208 + 1.0032 * Rv - 0.00033 * (Rv ** 2),
1.208 + 1.0032 * Rv - 0.00033 * (Rv**2),
]
)
# updated NIR curve from F04, note R dependence
nir_axebv_y = (0.63 * Rv - 0.84) * nir_axav_x ** 1.84
nir_axebv_y = (0.63 * Rv - 0.84) * nir_axav_x**1.84

optnir_axebv_y = np.concatenate([nir_axebv_y, opt_axebv_y])

Expand Down Expand Up @@ -390,7 +390,9 @@ def hi_abs_sed(self, params, hi_velocities, sed):
hi_sed = {}
for cspec in self.fluxes.keys():
hi_sed[cspec] = np.copy(sed[cspec])
indxs, = np.where(np.absolute((self.waves[cspec] - h_lines[0]) <= h_width))
(indxs,) = np.where(
np.absolute((self.waves[cspec] - h_lines[0]) <= h_width)
)
if len(indxs) > 0:
for i, cvel in enumerate(hi_velocities):
# compute the Ly-alpha abs: from Bohlin et al. (197?)
Expand Down Expand Up @@ -425,7 +427,7 @@ def SED_to_StarData(self, sed):
# populate the BAND info
sd.data["BAND"] = BandData("BAND")
sd.data["BAND"].fluxes = sed["BAND"] * (
u.erg / ((u.cm ** 2) * u.s * u.angstrom)
u.erg / ((u.cm**2) * u.s * u.angstrom)
)
for k, cband in enumerate(self.band_names):
sd.data["BAND"].band_fluxes[cband] = (sed["BAND"][k], 0.0)
Expand All @@ -435,13 +437,16 @@ def SED_to_StarData(self, sed):
# populate the spectral info
sd.data[cspec] = SpecData(cspec)
sd.data[cspec].fluxes = sed[cspec] * (
u.erg / ((u.cm ** 2) * u.s * u.angstrom)
u.erg / ((u.cm**2) * u.s * u.angstrom)
)

sd.data[cspec].waves = self.waves[cspec]
sd.data[cspec].n_waves = len(sd.data[cspec].waves)
sd.data[cspec].uncs = 0.0 * sd.data[cspec].fluxes
sd.data[cspec].npts = np.full((sd.data[cspec].n_waves), 1.0)
sd.data[cspec].wave_range = [min(sd.data[cspec].waves), max(sd.data[cspec].waves)]
sd.data[cspec].wave_range = [
min(sd.data[cspec].waves),
max(sd.data[cspec].waves),
]

return sd
4 changes: 3 additions & 1 deletion measure_extinction/plotting/plot_spec.py
Original file line number Diff line number Diff line change
Expand Up @@ -213,7 +213,9 @@ def plot_multi_spectra(
gkeys = list(starobs.data.keys())
gkeys.remove("BAND")
for src in gkeys:
starobs.data[src].rebin_constres(starobs.data[src].wave_range, rebin_res)
starobs.data[src].rebin_constres(
starobs.data[src].wave_range, rebin_res
)

# spread out the spectra if requested
# add extra whitespace when the luminosity class changes from main sequence to giant
Expand Down
24 changes: 13 additions & 11 deletions measure_extinction/stardata.py
Original file line number Diff line number Diff line change
Expand Up @@ -27,7 +27,7 @@
# const = 1e-26*1e7*1e-4*1e-4 = 1e-27
Jy_to_cgs_const = 1e-27 * const.c.to("micron/s").value

fluxunit = u.erg / ((u.cm ** 2) * u.s * u.angstrom)
fluxunit = u.erg / ((u.cm**2) * u.s * u.angstrom)


class BandData:
Expand Down Expand Up @@ -484,8 +484,8 @@ def get_band_fluxes(self):
# add the units
self.waves = self.waves * u.micron
self.wave_range = self.wave_range * u.micron
self.fluxes = self.fluxes * (u.erg / ((u.cm ** 2) * u.s * u.angstrom))
self.uncs = self.uncs * (u.erg / ((u.cm ** 2) * u.s * u.angstrom))
self.fluxes = self.fluxes * (u.erg / ((u.cm**2) * u.s * u.angstrom))
self.uncs = self.uncs * (u.erg / ((u.cm**2) * u.s * u.angstrom))

def get_band_mags_from_fluxes(self):
"""
Expand Down Expand Up @@ -606,7 +606,7 @@ def read_spectra(self, line, path="./"):
# open and read the spectrum
# ignore units warnings as non-standard units are explicitly handled a few lines later
with warnings.catch_warnings():
warnings.simplefilter('ignore', UnitsWarning)
warnings.simplefilter("ignore", UnitsWarning)
tdata = Table.read(full_filename)

self.waves = tdata["WAVELENGTH"].quantity
Expand All @@ -626,8 +626,8 @@ def read_spectra(self, line, path="./"):
if self.waves.unit == "MICRON":
self.waves = self.waves.value * u.micron
if self.fluxes.unit == "ERG/CM2/S/A":
self.fluxes = self.fluxes.value * (u.erg / ((u.cm ** 2) * u.s * u.angstrom))
self.uncs = self.uncs.value * (u.erg / ((u.cm ** 2) * u.s * u.angstrom))
self.fluxes = self.fluxes.value * (u.erg / ((u.cm**2) * u.s * u.angstrom))
self.uncs = self.uncs.value * (u.erg / ((u.cm**2) * u.s * u.angstrom))

# compute the min/max wavelengths
self.wave_range = (
Expand Down Expand Up @@ -714,8 +714,8 @@ def read_stis(self, line, path="./"):
self.read_spectra(line, path)

# add units
self.fluxes = self.fluxes.value * (u.erg / ((u.cm ** 2) * u.s * u.angstrom))
self.uncs = self.uncs.value * (u.erg / ((u.cm ** 2) * u.s * u.angstrom))
self.fluxes = self.fluxes.value * (u.erg / ((u.cm**2) * u.s * u.angstrom))
self.uncs = self.uncs.value * (u.erg / ((u.cm**2) * u.s * u.angstrom))

def read_spex(self, line, path="./", use_corfac=True, corfac=None):
"""
Expand Down Expand Up @@ -757,8 +757,8 @@ def read_spex(self, line, path="./", use_corfac=True, corfac=None):
self.uncs *= corfac

# add units
self.fluxes = self.fluxes.value * (u.erg / ((u.cm ** 2) * u.s * u.angstrom))
self.uncs = self.uncs.value * (u.erg / ((u.cm ** 2) * u.s * u.angstrom))
self.fluxes = self.fluxes.value * (u.erg / ((u.cm**2) * u.s * u.angstrom))
self.uncs = self.uncs.value * (u.erg / ((u.cm**2) * u.s * u.angstrom))

def read_irs(self, line, path="./", use_corfac=True, corfac=None):
"""
Expand Down Expand Up @@ -837,7 +837,9 @@ def read_miri_ifu(self, line, path="./"):
# add units
self.fluxes = self.fluxes.value * u.Jy
self.uncs = self.uncs.value * u.Jy
self.fluxes = self.fluxes.to(fluxunit, equivalencies=u.spectral_density(self.waves))
self.fluxes = self.fluxes.to(
fluxunit, equivalencies=u.spectral_density(self.waves)
)
self.uncs = self.uncs.to(fluxunit, equivalencies=u.spectral_density(self.waves))

def rebin_constres(self, waverange, resolution):
Expand Down
2 changes: 1 addition & 1 deletion measure_extinction/tests/test_rebin.py
Original file line number Diff line number Diff line change
Expand Up @@ -22,7 +22,7 @@ def test_spec_rebin_constres():
spec.fluxes = np.empty((a.size + b.size), dtype=a.dtype)
spec.fluxes[0::2] = a
spec.fluxes[1::2] = b
spec.fluxes *= u.erg / ((u.cm ** 2) * u.s * u.angstrom)
spec.fluxes *= u.erg / ((u.cm**2) * u.s * u.angstrom)
spec.uncs = spec.fluxes * 0.05
spec.npts = np.full((n_test), 1)

Expand Down
2 changes: 1 addition & 1 deletion measure_extinction/tests/test_stardata.py
Original file line number Diff line number Diff line change
Expand Up @@ -30,7 +30,7 @@ def test_units_stardata():
assert isinstance(star.data[cursrc].fluxes, u.Quantity)
assert isinstance(star.data[cursrc].uncs, u.Quantity)

fluxunit = u.erg / ((u.cm ** 2) * u.s * u.angstrom)
fluxunit = u.erg / ((u.cm**2) * u.s * u.angstrom)
# check that the wavelengths can be converted to microns and the
# flux units can be converted to spectral density
for cursrc in star.data.keys():
Expand Down
11 changes: 9 additions & 2 deletions measure_extinction/utils/calc_ext.py
Original file line number Diff line number Diff line change
Expand Up @@ -9,8 +9,15 @@
from measure_extinction.extdata import ExtData, AverageExtData


def calc_extinction(redstarname, compstarname, path, savepath="./", deredden=False,
elvebv=False, alav=False):
def calc_extinction(
redstarname,
compstarname,
path,
savepath="./",
deredden=False,
elvebv=False,
alav=False,
):
# read in the observed data for both stars
redstarobs = StarData("%s.dat" % redstarname.lower(), path=path)
compstarobs = StarData(
Expand Down
4 changes: 3 additions & 1 deletion measure_extinction/utils/make_all_tlusty_obsdata.py
Original file line number Diff line number Diff line change
Expand Up @@ -62,7 +62,9 @@ def decode_params(filename):


if __name__ == "__main__":
tlusty_models = glob.glob("/home/kgordon/Python/extstar_data/Models/Tlusty_2023/*v10.spec.gz")
tlusty_models = glob.glob(
"/home/kgordon/Python/extstar_data/Models/Tlusty_2023/*v10.spec.gz"
)

for cfname in tlusty_models:
# parse the filename to get the model parameters
Expand Down
12 changes: 6 additions & 6 deletions measure_extinction/utils/make_obsdata_from_model.py
Original file line number Diff line number Diff line change
Expand Up @@ -290,7 +290,7 @@ def make_obsdata_from_model(

# rebin to R=10000 for speed
# use a wavelength range that spans FUSE to Spitzer IRS
rbres = 10000.
rbres = 10000.0
wave_rebin, flux_rebin, npts_rebin = rebin_spectrum(
mwave.value, mflux.value, rbres, [912.0, 310000.0]
)
Expand All @@ -299,9 +299,7 @@ def make_obsdata_from_model(
otable = QTable()
otable["WAVELENGTH"] = Column(wave_rebin, unit=u.angstrom)
otable["FLUX"] = Column(flux_rebin, unit=fluxunit)
otable["SIGMA"] = Column(
flux_rebin * 0.01, unit=fluxunit
)
otable["SIGMA"] = Column(flux_rebin * 0.01, unit=fluxunit)
otable["NPTS"] = Column(npts_rebin)
otable.write(
"%s/Models/%s_full.fits" % (output_path, output_filebase), overwrite=True
Expand Down Expand Up @@ -429,7 +427,7 @@ def make_obsdata_from_model(
if show_plot:
fig, ax = plt.subplots(figsize=(10, 8))
# indxs, = np.where(npts_rebin > 0)
ax.plot(wave_rebin * 1e-4, flux_rebin * 2., "b-")
ax.plot(wave_rebin * 1e-4, flux_rebin * 2.0, "b-")
ax.plot(bandinfo.waves, bandinfo.fluxes, "ro")

(indxs,) = np.where(rb_stis_uv["NPTS"] > 0)
Expand Down Expand Up @@ -461,7 +459,9 @@ def make_obsdata_from_model(


if __name__ == "__main__":
mname = "/home/kgordon/Python/extstar_data/Models/Tlusty_2023/z050t23000g250v2.spec.gz"
mname = (
"/home/kgordon/Python/extstar_data/Models/Tlusty_2023/z050t23000g250v2.spec.gz"
)
model_params = {}
model_params["origin"] = "tlusty"
model_params["Teff"] = 23000.0
Expand Down
1 change: 1 addition & 0 deletions measure_extinction/utils/merge_iue_spec.py
Original file line number Diff line number Diff line change
@@ -1,4 +1,5 @@
import argparse

# import numpy as np
import pkg_resources

Expand Down
18 changes: 15 additions & 3 deletions measure_extinction/utils/merge_stis_spec.py
Original file line number Diff line number Diff line change
Expand Up @@ -117,12 +117,24 @@ def read_stis_archive_format(filename):
fig, ax = plt.subplots(nrows=1, ncols=1, figsize=(10, 5.5))

gvals = stable[0]["NPTS"] > 0
ax.plot(stable[0]["WAVELENGTH"][gvals], stable[0]["FLUX"][gvals], "k-", alpha=0.5, label="orig")
ax.plot(
stable[0]["WAVELENGTH"][gvals],
stable[0]["FLUX"][gvals],
"k-",
alpha=0.5,
label="orig",
)
gvals = rb_stis["NPTS"] > 0
ax.plot(rb_stis["WAVELENGTH"][gvals], rb_stis["FLUX"][gvals], "b-", alpha=0.5, label="merged")
ax.plot(
rb_stis["WAVELENGTH"][gvals],
rb_stis["FLUX"][gvals],
"b-",
alpha=0.5,
label="merged",
)

# set min/max ignoring Ly-alpha as it often has a strong core emission
gvals = (rb_stis["WAVELENGTH"] > 1300.) & (rb_stis["NPTS"] > 0)
gvals = (rb_stis["WAVELENGTH"] > 1300.0) & (rb_stis["NPTS"] > 0)
miny = np.nanmin(rb_stis["FLUX"][gvals])
maxy = np.nanmax(rb_stis["FLUX"][gvals])
delt = maxy - miny
Expand Down
2 changes: 1 addition & 1 deletion measure_extinction/utils/plot_bandpasses.py
Original file line number Diff line number Diff line change
Expand Up @@ -11,7 +11,7 @@ def plot_bandpasses(bands):
for band in bands:
bp = SpectralElement.from_file("%s%s.dat" % (band_path, band))
if "MIPS" in band:
wavelengths = bp.waveset * 10 ** 4
wavelengths = bp.waveset * 10**4
else:
wavelengths = bp.waveset
plt.plot(wavelengths, bp(bp.waveset), label=band)
Expand Down
2 changes: 1 addition & 1 deletion measure_extinction/utils/plot_mspec.py
Original file line number Diff line number Diff line change
Expand Up @@ -56,7 +56,7 @@ def plot_mspec_parser():
for k, cstarname in enumerate(starnames):
fstarname, file_path = get_full_starfile(cstarname)
starobs = StarData(fstarname, path=file_path)
starobs.plot(ax, norm_wave_range=[0.2, 0.3] * u.micron, yoffset=2 ** k)
starobs.plot(ax, norm_wave_range=[0.2, 0.3] * u.micron, yoffset=2**k)

# finish configuring the plot
ax.set_yscale("log")
Expand Down

0 comments on commit 8f5019f

Please sign in to comment.