Merge pull request #45 from sunpy/tweaks

Changelog plus doc tweaks/rename from pfsspy to pfss and remove fortran in plain text

.gitignore

@@ -262,6 +262,6 @@ package.json
 # Log files generated by 'vagrant up'
 *.log
 
-# repo specific
-examples/using_pfsspy/aia_map.fits
+# Repo specific
+examples/using_sunkit_magex_pfss/aia_map.fits
 docs/sg_execution_times.rst

.pre-commit-config.yaml

@@ -1,7 +1,7 @@
 repos:
     # This should be before any formatting hooks like isort
   - repo: https://github.com/astral-sh/ruff-pre-commit
-    rev: "v0.4.3"
+    rev: "v0.4.4"
     hooks:
       - id: ruff
         args: ["--fix"]

docs/changelog.rst

@@ -0,0 +1,28 @@
+.. _sunkit-magex-changelog:
+
+**************
+Full Changelog
+**************
+
+1.0.0 (2024-05-31)
+==================
+
+This is the first release and aims to keep the API the same from `pfsspy` to
+`sunkit_magex`.  The main difference is that you will need to replace the
+import like so:
+
+  .. code-block:: python
+
+    # Before
+    import pfsspy
+    from pfsspy import tracing
+
+    # After
+    from sunkit_magex import pfss as pfsspy
+    from sunkit_magex.pfss import tracing
+
+The main changes from the previous release of `pfsspy` are as follows:
+
+* The ``GongSynopticMap`` class has moved into `sunpy`, note that the new
+  class in ``sunpy`` is slightly different in that it does not modify the
+  header.

docs/conf.py

@@ -99,7 +99,7 @@
     "plot_gallery": True,
     "remove_config_comments": True,
     "subsection_order": ExplicitOrder([
-        "../examples/using_pfsspy",
+        "../examples/using_sunkit_magex_pfss",
         "../examples/finding_data",
         "../examples/utils",
         "../examples/internals",

docs/index.rst

@@ -15,3 +15,4 @@ This is the documentation for sunkit-magex.
    performance
    numerical_methods_pfss/index.rst
    synoptic_fits
+   changelog

docs/performance.rst

@@ -11,8 +11,8 @@ To enable this simply `install numba`_  and use `sunkit_magex.pfss` as normal.
 Streamline tracing
 ==================
 
-`sunkit_magex.pfss` has two streamline tracers: a pure python `sunkit_magex.pfss.tracing.PythonTracer` and a FORTRAN `sunkit_magex.pfss.tracing.FortranTracer`.
-The FORTRAN version is significantly faster, using the `streamtracer`_ package.
+`sunkit_magex.pfss` has two streamline tracers: a pure python `sunkit_magex.pfss.tracing.PythonTracer` and a complied tracer `sunkit_magex.pfss.tracing.FortranTracer`.
+The compiled version is significantly faster, using the `streamtracer`_ package.
 
 .. _numba: https://numba.pydata.org
 .. _install numba: http://numba.pydata.org/numba-doc/latest/user/installing.html

docs/reference/index.rst

@@ -1,6 +1,7 @@
-API reference
+API Reference
 =============
 
 .. toctree::
 
-    pfsspy
+    magex
+    pfss

docs/reference/magex.rst

@@ -0,0 +1 @@
+.. automodapi:: sunkit_magex

examples/internals/tracer_performance.py

@@ -68,7 +68,7 @@ def dipole_Br(r, theta):
 
 fig, ax = plt.subplots()
 ax.scatter(nseeds[1:len(times[0])], times[0][1:], label='python')
-ax.scatter(nseeds[1:], times[1][1:], label='fortran')
+ax.scatter(nseeds[1:], times[1][1:], label='compiled')
 
 pydt = (times[0][4] - times[0][3]) / (nseeds[4] - nseeds[3])
 ax.plot([1, 1e5], [pydt, 1e5 * pydt])

examples/testing/helpers.py

@@ -43,11 +43,11 @@ def pffspy_output(nphi, ns, nrho, rss, l, m):
     return sunkit_magex.pfss.pfss(pfss_input)
 
 
-def brss_pfsspy(nphi, ns, nrho, rss, l, m):
+def brss_pfss(nphi, ns, nrho, rss, l, m):
     # Return the radial component of the source surface mangetic field
     # for given input parameters
-    pfsspy_out = pffspy_output(nphi, ns, nrho, rss, l, m)
-    return pfsspy_out.bc[0][:, :, -1].T.astype(float)
+    pfss_out = pffspy_output(nphi, ns, nrho, rss, l, m)
+    return pfss_out.bc[0][:, :, -1].T.astype(float)
 
 
 def brss_analytic(nphi, ns, rss, l, m):
@@ -79,7 +79,7 @@ def open_flux_numeric(l, m, zss, nrho):
     """
     nphi = 360
     ns = 180
-    br = brss_pfsspy(nphi, ns, nrho, zss, l, m)
+    br = brss_pfss(nphi, ns, nrho, zss, l, m)
     return np.sum(np.abs(br)) * (4 * np.pi) / nphi / ns
 
 

examples/testing/plot_error_map.py

@@ -27,7 +27,7 @@
 ns = 180
 nr = 40
 rss = 2
-pfsspy_out = pffspy_output(nphi, ns, nr, rss, l, m)
+pfss_out = pffspy_output(nphi, ns, nr, rss, l, m)
 rss = rss * const.R_sun
 
 ###############################################################################
@@ -42,14 +42,14 @@
 theta, phi = np.meshgrid(theta, phi, indexing='ij')
 theta, phi = theta * u.rad, phi * u.deg
 seeds = SkyCoord(radius=rss, lat=theta.ravel(), lon=phi.ravel(),
-                 frame=pfsspy_out.coordinate_frame)
+                 frame=pfss_out.coordinate_frame)
 
 step_size = 1
 dthetas = []
 print(f'Tracing {step_size}...')
 # Trace
 tracer = pfss.tracing.FortranTracer(step_size=step_size)
-flines = tracer.trace(seeds, pfsspy_out)
+flines = tracer.trace(seeds, pfss_out)
 # Set a mask of open field lines
 mask = flines.connectivities.astype(bool).reshape(theta.shape)
 

examples/testing/plot_harmonic_comparison.py

@@ -7,7 +7,7 @@
 import matplotlib.pyplot as plt
 import matplotlib.ticker as mticker
 
-from examples.testing.helpers import LMAxes, brss_analytic, brss_pfsspy
+from examples.testing.helpers import LMAxes, brss_analytic, brss_pfss
 
 ###############################################################################
 # Compare the the `sunkit_magex.pfss` solution to the analytic solutions.
@@ -25,10 +25,10 @@
         print(f'l={l}, m={m}')
         ax = axs[l, m]
 
-        br_pfsspy = brss_pfsspy(nphi, ns, nrho, rss, l, m)
+        br_pfss = brss_pfss(nphi, ns, nrho, rss, l, m)
         br_actual = brss_analytic(nphi, ns, rss, l, m)
 
-        ax.plot(br_pfsspy[:, 15], label='sunkit_magex.pfss')
+        ax.plot(br_pfss[:, 15], label='sunkit_magex.pfss')
         ax.plot(br_actual[:, 15], label='analytic')
         if l == 1 and m == 0:
             ax.xaxis.set_major_formatter(mticker.StrMethodFormatter('{x}°'))

examples/testing/tracer_step_size.py

@@ -22,7 +22,7 @@
 
 ###############################################################################
 # Calculate PFSS solution
-pfsspy_out = pffspy_output(nphi, ns, nr, rss, l, m)
+pfss_out = pffspy_output(nphi, ns, nr, rss, l, m)
 
 ###############################################################################
 # Trace an array of field lines from the source surface down to the solar
@@ -37,7 +37,7 @@
 theta, phi = theta * u.rad, phi * u.deg
 rss = rss * const.R_sun
 seeds = SkyCoord(radius=rss, lat=theta.ravel(), lon=phi.ravel(),
-                 frame=pfsspy_out.coordinate_frame)
+                 frame=pfss_out.coordinate_frame)
 
 step_sizes = [32, 16, 8, 4, 2, 1, 0.5]
 dthetas = []
@@ -46,7 +46,7 @@
     print(f'Tracing {step_size}...')
     # Trace
     tracer = tracing.FortranTracer(step_size=step_size)
-    flines = tracer.trace(seeds, pfsspy_out)
+    flines = tracer.trace(seeds, pfss_out)
     # Set a mask of open field lines
     mask = flines.connectivities.astype(bool).reshape(theta.shape)
 

pyproject.toml

@@ -54,6 +54,8 @@ performance = [
 analytic = [
   "sympy",
 ]
+all = ["sunkit-magex[performance,analytic]"]
+dev = ["sunkit-magex[all,tests,docs]"]
 
 [project.urls]
 Homepage = "https://sunpy.org"

pytest.ini

@@ -44,3 +44,5 @@ filterwarnings =
     # Can be removed when https://github.com/dateutil/dateutil/issues/1314 is resolved
     # deprecated in Python 3.12, needs a release of dateutil 2.8.3 or higher
     ignore:datetime.datetime.utcfromtimestamp():DeprecationWarning
+    # When we depend on sunpy 6.0.0, this can go
+    ignore::sunpy.util.exceptions.SunpyMetadataWarning

sunkit_magex/pfss/map.py

@@ -1,8 +1,6 @@
 """
 Custom `sunpy.map.GenericMap` sub-classes for different magnetogram sources.
 """
-import astropy.units as u
-
 import sunpy.map
 
 __all__ = ['ADAPTMap']
@@ -24,24 +22,3 @@ def __init__(self, data, header, **kwargs):
     def is_datasource_for(cls, data, header, **kwargs):
         """Determines if header corresponds to an ADAPT map."""
         return header.get('model') == 'ADAPT'
-
-
-def _observer_coord_meta(observer_coord):
-    """
-    Convert an observer coordinate into FITS metadata.
-    """
-    new_obs_frame = sunpy.coordinates.HeliographicStonyhurst(
-        obstime=observer_coord.obstime)
-    observer_coord = observer_coord.transform_to(new_obs_frame)
-
-    new_meta = {'hglt_obs': observer_coord.lat.to_value(u.deg)}
-    new_meta['hgln_obs'] = observer_coord.lon.to_value(u.deg)
-    new_meta['dsun_obs'] = observer_coord.radius.to_value(u.m)
-    return new_meta
-
-
-def _earth_obs_coord_meta(obstime):
-    """
-    Return metadata for an Earth observer coordinate.
-    """
-    return _observer_coord_meta(sunpy.coordinates.get_earth(obstime))

sunkit_magex/pfss/tests/test_map.py

@@ -1,29 +1,12 @@
-
-import astropy.units as u
+import astropy.io
 
 import sunpy.map
 
 import sunkit_magex.pfss.map
 
 
-def test_gong_source(gong_map):
-    m = sunpy.map.Map(gong_map)
-    # Check round-trip is robust against sunpy changes to the meta
-    m = sunpy.map.Map(m.data, m.meta)
-    assert m.date.isot == '2020-09-01T13:04:00.000'
-    # Construct a WCS to check no warnings are thrown
-    m.wcs
-    # Check observer coordinate is populated
-    observer = m.coordinate_frame.observer
-    assert observer.obstime.isot == m.date.isot
-    assert observer.lon == 0 * u.deg
-    assert u.allclose(observer.lat, 7.20584924 * u.deg)
-    assert u.allclose(observer.radius, 1.50953137e+11 * u.m)
-
-
 def test_adapt_map(adapt_test_file):
-    import astropy.io
-    adapt_fits = astropy.io.fits.open(adapt_test_file)
-    for map_slice in adapt_fits[0].data:
-        m = sunpy.map.Map((map_slice, adapt_fits[0].header))
-        assert isinstance(m, sunkit_magex.pfss.map.ADAPTMap)
+    with astropy.io.fits.open(adapt_test_file) as adapt_fits:
+        for map_slice in adapt_fits[0].data:
+            m = sunpy.map.Map((map_slice, adapt_fits[0].header))
+            assert isinstance(m, sunkit_magex.pfss.map.ADAPTMap)

sunkit_magex/pfss/tests/test_tracers.py

@@ -11,7 +11,7 @@
 
 @pytest.fixture(params=[tracing.PythonTracer(),
                         tracing.FortranTracer()],
-                ids=['python', 'fortran'])
+                ids=['python', 'compiled'])
 def flines(dipole_result, request):
     tracer = request.param
     _, out = dipole_result

sunkit_magex/pfss/tracing.py

@@ -75,7 +75,7 @@ def coords_to_xyz(seeds, output):
 
 class FortranTracer(Tracer):
     r"""
-    Tracer using Fortran code.
+    Tracer using compiled code via streamtracer.
 
     Parameters
     ----------
@@ -103,7 +103,7 @@ def __init__(self, max_steps='auto', step_size=1):
             from streamtracer import StreamTracer
         except ModuleNotFoundError as e:
             raise RuntimeError(
-                'Using FortranTracer requires the streamtracer module, '
+                'Using this tracer requires the streamtracer module, '
                 'but streamtracer could not be loaded') from e
         self.max_steps = max_steps
         self.step_size = step_size

sunkit_magex/pfss/utils.py

@@ -9,11 +9,33 @@
 from astropy import units as u
 from astropy.wcs import WCS
 
+import sunpy.coordinates
 import sunpy.map
 import sunpy.time
 
 __all__ = ['fix_hmi_meta', 'load_adapt', 'carr_cea_wcs_header', 'is_cea_map', 'is_car_map', 'is_full_sun_synoptic_map', 'car_to_cea', 'roll_map']
 
+def _observer_coord_meta(observer_coord):
+    """
+    Convert an observer coordinate into FITS metadata.
+    """
+    new_obs_frame = sunpy.coordinates.HeliographicStonyhurst(
+        obstime=observer_coord.obstime)
+    observer_coord = observer_coord.transform_to(new_obs_frame)
+
+    new_meta = {'hglt_obs': observer_coord.lat.to_value(u.deg)}
+    new_meta['hgln_obs'] = observer_coord.lon.to_value(u.deg)
+    new_meta['dsun_obs'] = observer_coord.radius.to_value(u.m)
+    return new_meta
+
+
+def _earth_obs_coord_meta(obstime):
+    """
+    Return metadata for an Earth observer coordinate.
+    """
+    return _observer_coord_meta(sunpy.coordinates.get_earth(obstime))
+
+
 def fix_hmi_meta(hmi_map):
     """
     Fix non-compliant FITS metadata in HMI maps.
@@ -53,9 +75,8 @@ def fix_hmi_meta(hmi_map):
 
     # Fix observer coordinate
     if 'hglt_obs' not in hmi_map.meta:
-        from sunkit_magex import pfss
 
-        hmi_map.meta.update(pfss.map._earth_obs_coord_meta(hmi_map.meta['date-obs']))
+        hmi_map.meta.update(_earth_obs_coord_meta(hmi_map.meta['date-obs']))
 
 
 def load_adapt(adapt_path):
@@ -184,7 +205,7 @@ def is_full_sun_synoptic_map(m, error=False):
         return checks[projection](m, error)
     else:
         raise NotImplementedError(
-            f'is_full_sun_synoptic_map is only implemented for {list(checks.keys())} projections.'
+            f'is_full_sun_synoptic_map is only implemented for {list(checks.keys())} projections and not {projection}'
         )