From ebcdbba66b88552f446c9aa6f645b513d826c453 Mon Sep 17 00:00:00 2001 From: Alexandros Maragkoudakis Date: Thu, 29 Feb 2024 10:20:48 -0800 Subject: [PATCH] update plotting, ci fixes (#41) Nice job! --- pypahdb/decomposer.py | 59 +++++++++++++++++-------------- pypahdb/observation.py | 17 ++++----- pypahdb/tests/test_decomposer.py | 4 +-- pypahdb/tests/test_observation.py | 10 +++--- 4 files changed, 49 insertions(+), 41 deletions(-) diff --git a/pypahdb/decomposer.py b/pypahdb/decomposer.py index 6ea1cf8..fb0dc2e 100644 --- a/pypahdb/decomposer.py +++ b/pypahdb/decomposer.py @@ -60,7 +60,7 @@ def save_pdf(self, filename, header="", domaps=True, doplots=True): d = pdf.infodict() d["Title"] = "pyPAHdb Results Summary" d["Author"] = ( - "Dr. C. Boersma, Dr. M.J. Shannon, and Dr. A. " "Maragkoudakis" + "Dr. C. Boersma, Dr. M.J. Shannon, and Dr. A. Maragkoudakis" ) d["Producer"] = "NASA Ames Research Center" d["Creator"] = "pypahdb v{}(Python {}.{}.{})".format( @@ -158,7 +158,7 @@ def _fits_to_disk(hdr, filename): "Software used to create this file", ) hdr["AUTHOR"] = ( - "Dr. C. Boersma, Dr. M.J. Shannon, and Dr. A. " "Maragkoudakis", + "Dr. C. Boersma, Dr. M.J. Shannon, and Dr. A. Maragkoudakis", "Authors of the software", ) cards = [ @@ -382,11 +382,6 @@ def plot_fit(self, i=0, j=0): """ - # Enable quantity_support. - from astropy.visualization import quantity_support - - quantity_support() - # Create figure on shared axes. fig = plt.figure() gs = gridspec.GridSpec(4, 1, height_ratios=[2, 1, 2, 2], figure=fig) @@ -397,10 +392,18 @@ def plot_fit(self, i=0, j=0): ax1 = fig.add_subplot(gs[1], sharex=ax0) ax2 = fig.add_subplot(gs[2], sharex=ax0) ax3 = fig.add_subplot(gs[3], sharex=ax0) - - # Convenience defintions. + for ax in [ax0, ax2, ax3]: + fmt = ax.yaxis.get_major_formatter() + fmt.set_scientific(True) + fmt.set_powerlimits((-3, 1)) + for ax in [ax0, ax1, ax2]: + plt.setp(ax.get_xticklabels(), visible=False) + + # Convenience definitions. abscissa = self.spectrum.spectral_axis charge = self.charge + # Check if size of datapoints are too large and change marker size. + ms = 5 if len(abscissa) < 1000 else 2 # ax0: Best fit. data = self.spectrum.flux.T[:, i, j] @@ -412,8 +415,8 @@ def plot_fit(self, i=0, j=0): abscissa, data, yerr=unc, - marker="x", - ms=5, + marker=".", + ms=ms, mew=0.5, lw=0, color="black", @@ -422,12 +425,13 @@ def plot_fit(self, i=0, j=0): label="input", zorder=0, ) - ax0.plot(abscissa, model, label="fit", color="red", lw=1.5) + ax0.plot(abscissa, model, label="fit", color="tab:red", lw=1.5) error_str = "$error$=%-4.2f" % (self.error[i][j]) - ax0.text(0.025, 0.9, error_str, ha="left", va="center", transform=ax0.transAxes) + ax0.text(0.025, 0.88, error_str, ha="left", va="center", transform=ax0.transAxes) + ax0.set_ylabel(f'{self.spectrum.meta["colnames"][1]} [{self.spectrum.flux.unit}]') # ax1: Residual. - ax1.plot(abscissa, data - model, lw=1, label="residual", color="black") + ax1.plot(abscissa, data - model, lw=1, label="residual", color="gray") ax1.axhline(y=0, color="0.5", ls="--", dashes=(12, 16), zorder=-10, lw=0.5) # ax2: Size breakdown. @@ -435,8 +439,8 @@ def plot_fit(self, i=0, j=0): abscissa, data, yerr=unc, - marker="x", - ms=5, + marker=".", + ms=ms, mew=0.5, lw=0, color="black", @@ -444,23 +448,24 @@ def plot_fit(self, i=0, j=0): capsize=2, zorder=0, ) - ax2.plot(abscissa, model, color="red", lw=1.5) + ax2.plot(abscissa, model, color="tab:red", lw=1.5) ax2.plot( - abscissa, self.size["large"][:, i, j], label="large", lw=1, color="purple" + abscissa, self.size["large"][:, i, j], label="large", lw=1, color="tab:green" ) ax2.plot( - abscissa, self.size["small"][:, i, j], label="small", lw=1, color="crimson" + abscissa, self.size["small"][:, i, j], label="small", lw=1, color="tab:blue" ) size_str = "$f_{large}$=%3.1f" % (self.large_fraction[i][j]) - ax2.text(0.025, 0.9, size_str, ha="left", va="center", transform=ax2.transAxes) + ax2.text(0.025, 0.88, size_str, ha="left", va="center", transform=ax2.transAxes) + ax2.set_ylabel(f'{self.spectrum.meta["colnames"][1]} [{self.spectrum.flux.unit}]') # ax3: Charge breakdown. ax3.errorbar( abscissa, data, yerr=unc, - marker="x", - ms=5, + marker=".", + ms=ms, mew=0.5, lw=0, color="black", @@ -470,16 +475,18 @@ def plot_fit(self, i=0, j=0): ) ax3.plot(abscissa, model, color="red", lw=1.5) ax3.plot( - abscissa, charge["anion"][:, i, j], label="anion", lw=1, color="orange" + abscissa, charge["anion"][:, i, j], label="anion", lw=1, color="tab:orange" ) ax3.plot( - abscissa, charge["neutral"][:, i, j], label="neutral", lw=1, color="green" + abscissa, charge["neutral"][:, i, j], label="neutral", lw=1, color="tab:cyan" ) ax3.plot( - abscissa, charge["cation"][:, i, j], label="cation", lw=1, color="blue" + abscissa, charge["cation"][:, i, j], label="cation", lw=1, color="tab:purple" ) ion_str = "$f_{ionized}$=%3.1f" % (self.ionized_fraction[i][j]) - ax3.text(0.025, 0.9, ion_str, ha="left", va="center", transform=ax3.transAxes) + ax3.text(0.025, 0.88, ion_str, ha="left", va="center", transform=ax3.transAxes) + ax3.set_xlabel(f'{self.spectrum.meta["colnames"][0]} [{self.spectrum.spectral_axis.unit}]') + ax3.set_ylabel(f'{self.spectrum.meta["colnames"][1]} [{self.spectrum.flux.unit}]') # Set tick parameters and add legends to axes. for ax in (ax0, ax1, ax2, ax3): diff --git a/pypahdb/observation.py b/pypahdb/observation.py index 92d911f..8b330b5 100644 --- a/pypahdb/observation.py +++ b/pypahdb/observation.py @@ -43,6 +43,7 @@ def __init__(self, file_path): warnings.simplefilter("ignore", category=VerifyWarning) self.spectrum = Spectrum1D.read(self.file_path) + self.spectrum.meta['colnames'] = 3 * ["",] # Always work as if spectrum is a cube. if len(self.spectrum.flux.shape) == 1: @@ -92,6 +93,7 @@ def __init__(self, file_path): flux = h.data.T * u.Unit(h.header["BUNIT"]) wave = hdu[h0].data[h1] * u.Unit(hdu[h0].columns[h1].unit) self.spectrum = Spectrum1D(flux, spectral_axis=wave) + self.spectrum.meta['colnames'] = 3 * ["",] return None @@ -125,31 +127,30 @@ def __init__(self, file_path): data = ascii.read(self.file_path) # Always work as if spectrum is a cube. flux = np.reshape( - data["FLUX"].quantity, + data[data.colnames[1]].quantity, ( 1, 1, ) - + data["FLUX"].quantity.shape, + + data[data.colnames[1]].quantity.shape, ) # Create Spectrum1D object. - for name in data.colnames: - data.rename_column(name, name.upper()) - wave = data["WAVELENGTH"].quantity + wave = data[data.colnames[0]].quantity unc = None - if "FLUX_UNCERTAINTY" in data.colnames: + if len(data.colnames) > 2: unc = StdDevUncertainty( np.reshape( - data["FLUX_UNCERTAINTY"].quantity, + data[data.colnames[2]].quantity, ( 1, 1, ) - + (data["FLUX_UNCERTAINTY"].quantity.shape), + + (data[data.colnames[2]].quantity.shape), ) ) self.spectrum = Spectrum1D(flux, spectral_axis=wave, uncertainty=unc) + self.spectrum.meta['colnames'] = data.colnames hdr = "" for card in data.meta["keywords"].keys(): diff --git a/pypahdb/tests/test_decomposer.py b/pypahdb/tests/test_decomposer.py index 5b9834e..c69dd28 100644 --- a/pypahdb/tests/test_decomposer.py +++ b/pypahdb/tests/test_decomposer.py @@ -20,10 +20,10 @@ class DecomposerTestCase(unittest.TestCase): """Unit tests for `decomposer.py`.""" def setUp(self): - import pkg_resources + import importlib_resources import tempfile file_name = 'resources/sample_data_NGC7023.tbl' - file_path = pkg_resources.resource_filename('pypahdb', file_name) + file_path = importlib_resources.files('pypahdb') / file_name self.observation = Observation(file_path) self.decomposer = Decomposer(self.observation.spectrum) self.tmpdir = tempfile.gettempdir() diff --git a/pypahdb/tests/test_observation.py b/pypahdb/tests/test_observation.py index 3264faa..6587375 100644 --- a/pypahdb/tests/test_observation.py +++ b/pypahdb/tests/test_observation.py @@ -6,7 +6,7 @@ """ import unittest -import pkg_resources +import importlib_resources from pypahdb.observation import Observation @@ -17,21 +17,21 @@ class SpectrumTestCase(unittest.TestCase): def test_read_spectrum1d(self): """Can we create an instance of Observation from a Spectrum1D file?""" file_name = 'resources/sample_data_jwst.fits' - file_path = pkg_resources.resource_filename('pypahdb', file_name) + file_path = importlib_resources.files('pypahdb') / file_name assert isinstance(Observation(file_path), Observation) def test_read_ascii(self): """Can we create an instance of Observation from an ASCII file?""" file_name = 'resources/sample_data_NGC7023.tbl' - file_path = pkg_resources.resource_filename('pypahdb', file_name) + file_path = importlib_resources.files('pypahdb') / file_name assert isinstance(Observation(file_path), Observation) def test_read_fits(self): """Can we create an instance of Observation from a FITS file?""" file_name = 'resources/sample_data_NGC7023.fits' - file_path = pkg_resources.resource_filename('pypahdb', file_name) + file_path = importlib_resources.files('pypahdb') / file_name assert isinstance(Observation(file_path), Observation) @@ -44,7 +44,7 @@ def test_file_not_found(self): def test_file_malformed(self): """Can we detect when a file is malformed?""" file_name = 'resources/sample_malformed.fits' - file_path = pkg_resources.resource_filename('pypahdb', file_name) + file_path = importlib_resources.files('pypahdb') / file_name self.assertRaises(OSError, Observation, file_path)