From 22a9e4fe231420e0d42cfc554e9fdeabb7af9c4f Mon Sep 17 00:00:00 2001 From: wuziqi Date: Fri, 15 Sep 2023 22:46:02 +0800 Subject: [PATCH 01/10] Add feature "plot current sheet" --- examples/{field_lines.py => current_sheet.py} | 20 ++++------ sunkit_pyvista/plotter.py | 37 +++++++++++++++++++ 2 files changed, 45 insertions(+), 12 deletions(-) rename examples/{field_lines.py => current_sheet.py} (89%) diff --git a/examples/field_lines.py b/examples/current_sheet.py similarity index 89% rename from examples/field_lines.py rename to examples/current_sheet.py index 3e0be06..30a2bdf 100644 --- a/examples/field_lines.py +++ b/examples/current_sheet.py @@ -14,6 +14,7 @@ import sunpy.map from astropy.constants import R_sun from astropy.coordinates import SkyCoord + from matplotlib import colors from sunkit_magex import pfss from sunkit_magex.pfss import tracing @@ -44,13 +45,16 @@ # We load a gong_map from sunkit-magex gong_fname = get_gong_map() gong_map = sunpy.map.Map(gong_fname) +gong_map.plot() + +plotter = SunpyPlotter() # Now we plot the Gong Map to fill in the farside. plotter.plot_map(gong_map, cmap="hmimag", clip_interval=[1, 99] * u.percent) # Create points spaced between lat={-90, 90} degrees -lat = np.linspace(-np.pi / 2, np.pi / 2, 32, endpoint=False) +lat = np.linspace(-np.pi / 2, np.pi / 2, 16, endpoint=False) # Create 32 points spaced between long={0, 360} degrees -lon = np.linspace(0, 2 * np.pi, 32, endpoint=False) +lon = np.linspace(0, np.pi*2, 16, endpoint=False) # Make a 2D grid from these 1D points lat, lon = np.meshgrid(lat, lon, indexing="ij") # Create lon, lat and radial coordinate values by using a sunkit-magex @@ -80,15 +84,7 @@ def my_fline_color_func(field_line): # Plotting the field lines plotter.plot_field_lines(field_lines, color_func=my_fline_color_func) - -# Set the camera coordinate to view the plot correctly -camera_coord = SkyCoord( - 0 * u.deg, - 0 * u.deg, - 6 * R_sun, - frame=frames.HeliographicStonyhurst, - obstime=LOW_RES_AIA_193.date, -) -plotter.set_camera_coordinate(camera_coord) +# Plotting the current sheet +plotter.plot_current_sheet(output_) plotter.show() diff --git a/sunkit_pyvista/plotter.py b/sunkit_pyvista/plotter.py index 3416028..9d91c04 100644 --- a/sunkit_pyvista/plotter.py +++ b/sunkit_pyvista/plotter.py @@ -12,6 +12,7 @@ from sunpy.coordinates import HeliocentricInertial, Helioprojective from sunpy.coordinates.utils import get_rectangle_coordinates from sunpy.map.maputils import all_corner_coords_from_map +from pfsspy.coords import strum2cart from sunkit_pyvista.utils import get_limb_coordinates @@ -509,6 +510,42 @@ def color_func(field_line): self._add_mesh_to_dict(block_name="field_lines", mesh=spline) + def plot_current_sheet(self, pfss_out, **kwargs): + """ + Plot current sheet (Br=0) from `pfsspy`. + + Parameters + ---------- + pfss_out: `pfsspy.output.Output` + Magnetic field calculated by pfss. + **kwargs + Keyword arguments are handed to `pyvista.Plotter.add_mesh`. + ------- + + """ + lon0 = ((pfss_out.input_map.meta['crval1']) * u.deg).to(u.rad).value + + sc_vect = pfss_out.grid.sc + pc_vect = np.insert(pfss_out.grid.pc, 0, pfss_out.grid.pc[-1]) + rg_vect = pfss_out.grid.rg + + c = SkyCoord(lon0*u.rad, 0.*u.rad, frame=HeliocentricInertial) + c = c.transform_to(self.coordinate_frame) + + bc_r = np.insert(pfss_out.bc[0], 0, pfss_out.bc[0][-1], axis=0) + + [s_arr, p_arr, r_arr] = np.meshgrid(sc_vect, pc_vect, rg_vect) + x_arr, y_arr, z_arr = strum2cart(r_arr, s_arr, p_arr + c.lon.to(u.rad).value) + + pfss_out_pv = pv.StructuredGrid(x_arr, y_arr, z_arr) + pfss_out_pv['Br'] = bc_r.ravel('F') + + isos_br = pfss_out_pv.contour(isosurfaces=1, rng=[0, 0]) + + self.plotter.add_mesh(isos_br,**kwargs) + self._add_mesh_to_dict(block_name='current_sheet', mesh=isos_br) + self.plotter.show_grid() + def save(self, filepath, *, overwrite=False): """ Save all the meshes. From 6d2ebf9c8df3e565abdefcfee2feb361ffefc83c Mon Sep 17 00:00:00 2001 From: wuziqi Date: Wed, 27 Sep 2023 13:26:38 +0800 Subject: [PATCH 02/10] Add changelog file --- changelog/143.feature.rst | 1 + 1 file changed, 1 insertion(+) create mode 100644 changelog/143.feature.rst diff --git a/changelog/143.feature.rst b/changelog/143.feature.rst new file mode 100644 index 0000000..8853f5b --- /dev/null +++ b/changelog/143.feature.rst @@ -0,0 +1 @@ +Add feature 'plot_current_sheet' in sunkit-pyvista. \ No newline at end of file From 73be5da8c1eab5daaa914e1874e76f0e4ebe864c Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Wed, 27 Sep 2023 03:40:08 +0000 Subject: [PATCH 03/10] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- examples/current_sheet.py | 3 +-- sunkit_pyvista/plotter.py | 13 ++++++------- 2 files changed, 7 insertions(+), 9 deletions(-) diff --git a/examples/current_sheet.py b/examples/current_sheet.py index 30a2bdf..ce677af 100644 --- a/examples/current_sheet.py +++ b/examples/current_sheet.py @@ -14,7 +14,6 @@ import sunpy.map from astropy.constants import R_sun from astropy.coordinates import SkyCoord - from matplotlib import colors from sunkit_magex import pfss from sunkit_magex.pfss import tracing @@ -54,7 +53,7 @@ # Create points spaced between lat={-90, 90} degrees lat = np.linspace(-np.pi / 2, np.pi / 2, 16, endpoint=False) # Create 32 points spaced between long={0, 360} degrees -lon = np.linspace(0, np.pi*2, 16, endpoint=False) +lon = np.linspace(0, np.pi * 2, 16, endpoint=False) # Make a 2D grid from these 1D points lat, lon = np.meshgrid(lat, lon, indexing="ij") # Create lon, lat and radial coordinate values by using a sunkit-magex diff --git a/sunkit_pyvista/plotter.py b/sunkit_pyvista/plotter.py index 9d91c04..fb0955d 100644 --- a/sunkit_pyvista/plotter.py +++ b/sunkit_pyvista/plotter.py @@ -9,10 +9,10 @@ from astropy.visualization import AsymmetricPercentileInterval from astropy.visualization.wcsaxes import Quadrangle from matplotlib import colors +from pfsspy.coords import strum2cart from sunpy.coordinates import HeliocentricInertial, Helioprojective from sunpy.coordinates.utils import get_rectangle_coordinates from sunpy.map.maputils import all_corner_coords_from_map -from pfsspy.coords import strum2cart from sunkit_pyvista.utils import get_limb_coordinates @@ -521,15 +521,14 @@ def plot_current_sheet(self, pfss_out, **kwargs): **kwargs Keyword arguments are handed to `pyvista.Plotter.add_mesh`. ------- - """ - lon0 = ((pfss_out.input_map.meta['crval1']) * u.deg).to(u.rad).value + lon0 = ((pfss_out.input_map.meta["crval1"]) * u.deg).to(u.rad).value sc_vect = pfss_out.grid.sc pc_vect = np.insert(pfss_out.grid.pc, 0, pfss_out.grid.pc[-1]) rg_vect = pfss_out.grid.rg - c = SkyCoord(lon0*u.rad, 0.*u.rad, frame=HeliocentricInertial) + c = SkyCoord(lon0 * u.rad, 0.0 * u.rad, frame=HeliocentricInertial) c = c.transform_to(self.coordinate_frame) bc_r = np.insert(pfss_out.bc[0], 0, pfss_out.bc[0][-1], axis=0) @@ -538,12 +537,12 @@ def plot_current_sheet(self, pfss_out, **kwargs): x_arr, y_arr, z_arr = strum2cart(r_arr, s_arr, p_arr + c.lon.to(u.rad).value) pfss_out_pv = pv.StructuredGrid(x_arr, y_arr, z_arr) - pfss_out_pv['Br'] = bc_r.ravel('F') + pfss_out_pv["Br"] = bc_r.ravel("F") isos_br = pfss_out_pv.contour(isosurfaces=1, rng=[0, 0]) - self.plotter.add_mesh(isos_br,**kwargs) - self._add_mesh_to_dict(block_name='current_sheet', mesh=isos_br) + self.plotter.add_mesh(isos_br, **kwargs) + self._add_mesh_to_dict(block_name="current_sheet", mesh=isos_br) self.plotter.show_grid() def save(self, filepath, *, overwrite=False): From 3f83599d773a36fd9c43644bd097fb18c21dce9f Mon Sep 17 00:00:00 2001 From: wuziqi Date: Wed, 27 Sep 2023 13:29:54 +0800 Subject: [PATCH 04/10] add blank line --- changelog/143.feature.rst | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/changelog/143.feature.rst b/changelog/143.feature.rst index 8853f5b..11345ab 100644 --- a/changelog/143.feature.rst +++ b/changelog/143.feature.rst @@ -1 +1 @@ -Add feature 'plot_current_sheet' in sunkit-pyvista. \ No newline at end of file +Add feature 'plot_current_sheet' in sunkit-pyvista. From f91f1f9e06ddbb8ed78e58959c3b92bc03ef95a2 Mon Sep 17 00:00:00 2001 From: wuziqi Date: Tue, 3 Oct 2023 20:36:34 +0800 Subject: [PATCH 05/10] add test current sheet figure. But fail to run any tests in the folder. message: "_jb_pytest_runner.py: error: unrecognized arguments: --doctest-rst" --- sunkit_pyvista/tests/test_plotting.py | 36 +++++++++++++++++++++++++++ 1 file changed, 36 insertions(+) diff --git a/sunkit_pyvista/tests/test_plotting.py b/sunkit_pyvista/tests/test_plotting.py index 9a3e6dd..0d8c68d 100644 --- a/sunkit_pyvista/tests/test_plotting.py +++ b/sunkit_pyvista/tests/test_plotting.py @@ -72,3 +72,39 @@ def test_plot_map_with_functionality( plotter.load(filepath) plotter.show(cpos=(0, 1, 0), before_close_callback=verify_cache_image) + + +def test_field_lines_figure(aia171_test_map, plotter, verify_cache_image): + gong_fname = get_gong_map() + gong_map = smap.Map(gong_fname) + nrho = 35 + rss = 2.5 + lat = np.linspace(-np.pi / 2, np.pi / 2, 8, endpoint=False) + lon = np.linspace(0, 2 * np.pi, 8, endpoint=False) + lat, lon = np.meshgrid(lat, lon, indexing="ij") + lat, lon = lat.ravel() * u.rad, lon.ravel() * u.rad + radius = 1.2 + tracer = tracing.PythonTracer() + input_ = pfsspy.Input(gong_map, nrho, rss) + output_ = pfsspy.pfss(input_) + seeds = SkyCoord(lon, lat, radius * const.R_sun, frame=gong_map.coordinate_frame) + field_lines = tracer.trace(seeds, output_) + + def color_function(field_line): + norm = colors.LogNorm(vmin=1, vmax=1000) + cmap = plt.get_cmap("magma") + return cmap(norm(np.abs(field_line.expansion_factor))) + + plotter.plot_map(aia171_test_map) + plotter.plot_field_lines(field_lines, color_func=color_function) + plotter.show(cpos=(0, 1, 0), before_close_callback=verify_cache_image) + +def test_current_sheet_figure(plotter,verify_cache_image): + gong_fname = get_gong_map() + gong_map = smap.Map(gong_fname) + nrho = 35 + rss = 2.5 + input_ = pfsspy.Input(gong_map, nrho, rss) + output_ = pfsspy.pfss(input_) + plotter.plot_current_sheet(output_) + plotter.show(before_close_callback=verify_cache_image) From 0f3a875242c661d475a9f3c12e439fcd21b48c53 Mon Sep 17 00:00:00 2001 From: wuziqi Date: Tue, 3 Oct 2023 20:46:30 +0800 Subject: [PATCH 06/10] Reformat --- sunkit_pyvista/tests/test_plotting.py | 11 ++++++----- 1 file changed, 6 insertions(+), 5 deletions(-) diff --git a/sunkit_pyvista/tests/test_plotting.py b/sunkit_pyvista/tests/test_plotting.py index 0d8c68d..c111760 100644 --- a/sunkit_pyvista/tests/test_plotting.py +++ b/sunkit_pyvista/tests/test_plotting.py @@ -27,10 +27,10 @@ def plotter(): def test_plot_map_with_functionality( - aia171_test_map, - plotter, - verify_cache_image, - tmp_path, + aia171_test_map, + plotter, + verify_cache_image, + tmp_path, ): plotter.plot_map(aia171_test_map, clip_interval=(0, 99) * u.percent) plotter.plot_solar_axis() @@ -99,7 +99,8 @@ def color_function(field_line): plotter.plot_field_lines(field_lines, color_func=color_function) plotter.show(cpos=(0, 1, 0), before_close_callback=verify_cache_image) -def test_current_sheet_figure(plotter,verify_cache_image): + +def test_current_sheet_figure(plotter, verify_cache_image): gong_fname = get_gong_map() gong_map = smap.Map(gong_fname) nrho = 35 From 93517b9294342148acd88a08bd17ba75c4c0b426 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Tue, 3 Oct 2023 12:45:37 +0000 Subject: [PATCH 07/10] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci From 27e99513fbc26e77ba633861833d9fde80cc5006 Mon Sep 17 00:00:00 2001 From: wuziqi Date: Tue, 3 Oct 2023 21:12:18 +0800 Subject: [PATCH 08/10] Change the order of tests --- sunkit_pyvista/tests/test_plotting.py | 22 +++++++++++----------- 1 file changed, 11 insertions(+), 11 deletions(-) diff --git a/sunkit_pyvista/tests/test_plotting.py b/sunkit_pyvista/tests/test_plotting.py index c111760..50828e8 100644 --- a/sunkit_pyvista/tests/test_plotting.py +++ b/sunkit_pyvista/tests/test_plotting.py @@ -74,6 +74,17 @@ def test_plot_map_with_functionality( plotter.show(cpos=(0, 1, 0), before_close_callback=verify_cache_image) +def test_current_sheet_figure(plotter, verify_cache_image): + gong_fname = get_gong_map() + gong_map = smap.Map(gong_fname) + nrho = 35 + rss = 2.5 + input_ = pfsspy.Input(gong_map, nrho, rss) + output_ = pfsspy.pfss(input_) + plotter.plot_current_sheet(output_) + plotter.show(before_close_callback=verify_cache_image) + + def test_field_lines_figure(aia171_test_map, plotter, verify_cache_image): gong_fname = get_gong_map() gong_map = smap.Map(gong_fname) @@ -98,14 +109,3 @@ def color_function(field_line): plotter.plot_map(aia171_test_map) plotter.plot_field_lines(field_lines, color_func=color_function) plotter.show(cpos=(0, 1, 0), before_close_callback=verify_cache_image) - - -def test_current_sheet_figure(plotter, verify_cache_image): - gong_fname = get_gong_map() - gong_map = smap.Map(gong_fname) - nrho = 35 - rss = 2.5 - input_ = pfsspy.Input(gong_map, nrho, rss) - output_ = pfsspy.pfss(input_) - plotter.plot_current_sheet(output_) - plotter.show(before_close_callback=verify_cache_image) From 881fbe9383e8c25fa2c37a7846a33050c6f9876c Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Tue, 3 Oct 2023 12:55:56 +0000 Subject: [PATCH 09/10] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- sunkit_pyvista/tests/test_plotting.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/sunkit_pyvista/tests/test_plotting.py b/sunkit_pyvista/tests/test_plotting.py index 50828e8..dbc15ce 100644 --- a/sunkit_pyvista/tests/test_plotting.py +++ b/sunkit_pyvista/tests/test_plotting.py @@ -27,10 +27,10 @@ def plotter(): def test_plot_map_with_functionality( - aia171_test_map, - plotter, - verify_cache_image, - tmp_path, + aia171_test_map, + plotter, + verify_cache_image, + tmp_path, ): plotter.plot_map(aia171_test_map, clip_interval=(0, 99) * u.percent) plotter.plot_solar_axis() From 6f3e509ac47bb3a91a5bad9cb6f56dedee9d9d29 Mon Sep 17 00:00:00 2001 From: Nabil Freij Date: Wed, 12 Jun 2024 18:15:13 -0700 Subject: [PATCH 10/10] rebase plus update example --- changelog/143.feature.rst | 2 +- ...sheet.py => field_lines_current_sheets.py} | 62 +++++++++++-------- sunkit_pyvista/plotter.py | 28 +++------ sunkit_pyvista/tests/test_magex.py | 28 ++++++++- sunkit_pyvista/tests/test_plotting.py | 43 +------------ 5 files changed, 73 insertions(+), 90 deletions(-) rename examples/{current_sheet.py => field_lines_current_sheets.py} (74%) diff --git a/changelog/143.feature.rst b/changelog/143.feature.rst index 11345ab..969b6c7 100644 --- a/changelog/143.feature.rst +++ b/changelog/143.feature.rst @@ -1 +1 @@ -Add feature 'plot_current_sheet' in sunkit-pyvista. +Added a method to plot a current sheet via :meth:`~sunkit_pyvista.plotter.SunpyPlotter.plot_current_sheet`. diff --git a/examples/current_sheet.py b/examples/field_lines_current_sheets.py similarity index 74% rename from examples/current_sheet.py rename to examples/field_lines_current_sheets.py index ce677af..104f2ba 100644 --- a/examples/current_sheet.py +++ b/examples/field_lines_current_sheets.py @@ -1,11 +1,14 @@ """ -====================================== -Plotting Field Lines from sunkit-magex -====================================== +========================================================= +Plotting field lines and current sheets from sunkit-magex +========================================================= ``sunkit-pyvista`` can be used to plot field lines from ``sunkit-magex``. This example requires the `streamtracer `__ package. + +This example uses specific sample data from ``sunkit-magex`` to demonstrate the plotting of field lines and current sheets. +If you want to adapt this, you have to make sure that the input magnetic field data has the same date as the AIA image. """ import astropy.units as u @@ -24,36 +27,16 @@ from sunkit_pyvista.sample import LOW_RES_AIA_193 ############################################################################### -# We will be using an AIA 193 image from the sunpy sample data as the base image. - -# sphinx_gallery_defer_figures - -# Start by creating a plotter -plotter = SunpyPlotter() - -# Plot a map -plotter.plot_map(LOW_RES_AIA_193, clip_interval=[1, 99] * u.percent, assume_spherical_screen=False) -# Add an arrow to show the solar rotation axis -plotter.plot_solar_axis() - -############################################################################### -# We now need to do the magnetic field extrapolation using ``sunkit-magex``. - -# sphinx_gallery_defer_figures +# We will start by doing the magnetic field extrapolation using ``sunkit-magex``. # We load a gong_map from sunkit-magex gong_fname = get_gong_map() gong_map = sunpy.map.Map(gong_fname) -gong_map.plot() - -plotter = SunpyPlotter() -# Now we plot the Gong Map to fill in the farside. -plotter.plot_map(gong_map, cmap="hmimag", clip_interval=[1, 99] * u.percent) -# Create points spaced between lat={-90, 90} degrees -lat = np.linspace(-np.pi / 2, np.pi / 2, 16, endpoint=False) +# Create points spaced between lat={-89, 89} degrees +lat = np.linspace(-np.pi / 2, np.pi / 2, 32, endpoint=False) # Create 32 points spaced between long={0, 360} degrees -lon = np.linspace(0, np.pi * 2, 16, endpoint=False) +lon = np.linspace(1, 2 * np.pi, 32, endpoint=False) # Make a 2D grid from these 1D points lat, lon = np.meshgrid(lat, lon, indexing="ij") # Create lon, lat and radial coordinate values by using a sunkit-magex @@ -68,6 +51,21 @@ tracer = tracing.FortranTracer() field_lines = tracer.trace(seeds, output_) +############################################################################### +# We will be using an AIA 193 image from the sunpy sample data as the base image. + +# sphinx_gallery_defer_figures + +# Start by creating a plotter +plotter = SunpyPlotter() + +# Plot the AIA map +plotter.plot_map(LOW_RES_AIA_193, clip_interval=[1, 99] * u.percent, assume_spherical_screen=False) +# Add an arrow to show the solar rotation axis +plotter.plot_solar_axis() +# Now we plot the Gong Map to fill in the farside. +plotter.plot_map(gong_map, cmap="hmimag", clip_interval=[1, 99] * u.percent) + ############################################################################### # We can also specify a color function while plotting the field lines. # This function takes a single field line, and returns a color either @@ -86,4 +84,14 @@ def my_fline_color_func(field_line): # Plotting the current sheet plotter.plot_current_sheet(output_) +# Set the camera coordinate to view the plot correctly +camera_coord = SkyCoord( + 0 * u.deg, + 0 * u.deg, + 6 * R_sun, + frame=frames.HeliographicStonyhurst, + obstime=LOW_RES_AIA_193.date, +) +plotter.set_camera_coordinate(camera_coord) + plotter.show() diff --git a/sunkit_pyvista/plotter.py b/sunkit_pyvista/plotter.py index fb0955d..26f5b7d 100644 --- a/sunkit_pyvista/plotter.py +++ b/sunkit_pyvista/plotter.py @@ -9,7 +9,7 @@ from astropy.visualization import AsymmetricPercentileInterval from astropy.visualization.wcsaxes import Quadrangle from matplotlib import colors -from pfsspy.coords import strum2cart +from sunkit_magex.pfss.coords import strum2cart from sunpy.coordinates import HeliocentricInertial, Helioprojective from sunpy.coordinates.utils import get_rectangle_coordinates from sunpy.map.maputils import all_corner_coords_from_map @@ -512,38 +512,28 @@ def color_func(field_line): def plot_current_sheet(self, pfss_out, **kwargs): """ - Plot current sheet (Br=0) from `pfsspy`. + Plot current sheet "(Br=0)". Parameters ---------- - pfss_out: `pfsspy.output.Output` - Magnetic field calculated by pfss. - **kwargs - Keyword arguments are handed to `pyvista.Plotter.add_mesh`. - ------- + pfss_out : `sunkit_magex.pfss.output.Output` + Magnetic field calculated by `sunkit_magex`. + **kwargs : + Keyword arguments are passed to `pyvista.Plotter.add_mesh`. """ - lon0 = ((pfss_out.input_map.meta["crval1"]) * u.deg).to(u.rad).value - sc_vect = pfss_out.grid.sc pc_vect = np.insert(pfss_out.grid.pc, 0, pfss_out.grid.pc[-1]) rg_vect = pfss_out.grid.rg - - c = SkyCoord(lon0 * u.rad, 0.0 * u.rad, frame=HeliocentricInertial) - c = c.transform_to(self.coordinate_frame) - + coord = SkyCoord((pfss_out.input_map.meta["crval1"] * u.deg).to(u.rad), 0.0 * u.rad, frame=HeliocentricInertial) + coord = coord.transform_to(self.coordinate_frame) bc_r = np.insert(pfss_out.bc[0], 0, pfss_out.bc[0][-1], axis=0) - [s_arr, p_arr, r_arr] = np.meshgrid(sc_vect, pc_vect, rg_vect) - x_arr, y_arr, z_arr = strum2cart(r_arr, s_arr, p_arr + c.lon.to(u.rad).value) - + x_arr, y_arr, z_arr = strum2cart(r_arr, s_arr, p_arr + coord.lon.to(u.rad).value) pfss_out_pv = pv.StructuredGrid(x_arr, y_arr, z_arr) pfss_out_pv["Br"] = bc_r.ravel("F") - isos_br = pfss_out_pv.contour(isosurfaces=1, rng=[0, 0]) - self.plotter.add_mesh(isos_br, **kwargs) self._add_mesh_to_dict(block_name="current_sheet", mesh=isos_br) - self.plotter.show_grid() def save(self, filepath, *, overwrite=False): """ diff --git a/sunkit_pyvista/tests/test_magex.py b/sunkit_pyvista/tests/test_magex.py index 1519253..f113cf1 100644 --- a/sunkit_pyvista/tests/test_magex.py +++ b/sunkit_pyvista/tests/test_magex.py @@ -1,3 +1,7 @@ +""" +This file contains tests for any methods that use sunkit-magex. +""" + import astropy.constants as const import astropy.units as u import matplotlib.pyplot as plt @@ -10,6 +14,8 @@ from sunkit_pyvista import SunpyPlotter +pv.OFF_SCREEN = True + def test_field_lines_figure(aia171_test_map, plotter, verify_cache_image): pfss = pytest.importorskip("sunkit_magex.pfss") @@ -26,7 +32,7 @@ def test_field_lines_figure(aia171_test_map, plotter, verify_cache_image): lat, lon = np.meshgrid(lat, lon, indexing="ij") lat, lon = lat.ravel() * u.rad, lon.ravel() * u.rad radius = 1.2 - tracer = tracing.PythonTracer() + tracer = tracing.FortranTracer() input_ = pfss.Input(gong_map, nrho, rss) output_ = pfss.pfss(input_) seeds = SkyCoord(lon, lat, radius * const.R_sun, frame=gong_map.coordinate_frame) @@ -42,7 +48,7 @@ def color_function(field_line): plotter.show(cpos=(0, 1, 0), before_close_callback=verify_cache_image) -def test_field_lines_and_color_func(plotter): +def test_field_lines_and_color_func(plotter, verify_cache_image): pfss = pytest.importorskip("sunkit_magex.pfss") from sunkit_magex.pfss import tracing @@ -57,7 +63,7 @@ def test_field_lines_and_color_func(plotter): lat, lon = np.meshgrid(lat, lon, indexing="ij") lat, lon = lat.ravel() * u.rad, lon.ravel() * u.rad radius = 1.2 - tracer = tracing.PythonTracer() + tracer = tracing.FortranTracer() input_ = pfss.Input(gong_map, nrho, rss) output_ = pfss.pfss(input_) seeds = SkyCoord(lon, lat, radius * const.R_sun, frame=gong_map.coordinate_frame) @@ -72,3 +78,19 @@ def color_func(field_line): plotter = SunpyPlotter() plotter.plot_field_lines(field_lines, color_func=color_func) + plotter.show(before_close_callback=verify_cache_image) + + +def test_current_sheet_figure(plotter, verify_cache_image): + pfss = pytest.importorskip("sunkit_magex.pfss") + + from sunkit_magex.pfss.sample_data import get_gong_map + + gong_fname = get_gong_map() + gong_map = smap.Map(gong_fname) + nrho = 35 + rss = 2.5 + input_ = pfss.Input(gong_map, nrho, rss) + output_ = pfss.pfss(input_) + plotter.plot_current_sheet(output_) + plotter.show(before_close_callback=verify_cache_image) diff --git a/sunkit_pyvista/tests/test_plotting.py b/sunkit_pyvista/tests/test_plotting.py index dbc15ce..f1f2d8e 100644 --- a/sunkit_pyvista/tests/test_plotting.py +++ b/sunkit_pyvista/tests/test_plotting.py @@ -1,11 +1,11 @@ """ -This file contains figure comparison tests. +This file contains tests for the main plotting routines. """ import astropy.constants as const import astropy.units as u import pytest -import pyvista +import pyvista as pv import sunpy.data.test as test import sunpy.map as smap from astropy.coordinates import SkyCoord @@ -13,7 +13,7 @@ from sunkit_pyvista import SunpyPlotter -pyvista.OFF_SCREEN = True +pv.OFF_SCREEN = True @pytest.fixture() @@ -72,40 +72,3 @@ def test_plot_map_with_functionality( plotter.load(filepath) plotter.show(cpos=(0, 1, 0), before_close_callback=verify_cache_image) - - -def test_current_sheet_figure(plotter, verify_cache_image): - gong_fname = get_gong_map() - gong_map = smap.Map(gong_fname) - nrho = 35 - rss = 2.5 - input_ = pfsspy.Input(gong_map, nrho, rss) - output_ = pfsspy.pfss(input_) - plotter.plot_current_sheet(output_) - plotter.show(before_close_callback=verify_cache_image) - - -def test_field_lines_figure(aia171_test_map, plotter, verify_cache_image): - gong_fname = get_gong_map() - gong_map = smap.Map(gong_fname) - nrho = 35 - rss = 2.5 - lat = np.linspace(-np.pi / 2, np.pi / 2, 8, endpoint=False) - lon = np.linspace(0, 2 * np.pi, 8, endpoint=False) - lat, lon = np.meshgrid(lat, lon, indexing="ij") - lat, lon = lat.ravel() * u.rad, lon.ravel() * u.rad - radius = 1.2 - tracer = tracing.PythonTracer() - input_ = pfsspy.Input(gong_map, nrho, rss) - output_ = pfsspy.pfss(input_) - seeds = SkyCoord(lon, lat, radius * const.R_sun, frame=gong_map.coordinate_frame) - field_lines = tracer.trace(seeds, output_) - - def color_function(field_line): - norm = colors.LogNorm(vmin=1, vmax=1000) - cmap = plt.get_cmap("magma") - return cmap(norm(np.abs(field_line.expansion_factor))) - - plotter.plot_map(aia171_test_map) - plotter.plot_field_lines(field_lines, color_func=color_function) - plotter.show(cpos=(0, 1, 0), before_close_callback=verify_cache_image)