From 5fae8baac5194c225bf3c6b68cbd71588b73bbac Mon Sep 17 00:00:00 2001 From: gillyspace27 Date: Sun, 22 Sep 2024 00:27:10 -0600 Subject: [PATCH 1/9] begins to generalize past AIA-only data --- sunkit_image/radial.py | 121 +++++++++++++++++++++++------------------ 1 file changed, 68 insertions(+), 53 deletions(-) diff --git a/sunkit_image/radial.py b/sunkit_image/radial.py index 9a460ae5..d375a5a9 100644 --- a/sunkit_image/radial.py +++ b/sunkit_image/radial.py @@ -3,6 +3,7 @@ radius. """ +import matplotlib.pyplot as plt import numpy as np from tqdm import tqdm @@ -10,6 +11,7 @@ import sunpy.map from sunpy.coordinates import frames +from sunpy.map import Map from sunkit_image.utils import ( apply_upsilon, @@ -207,6 +209,8 @@ def intensity_enhance( enhancement[map_r < normalization_radius] = 1 # Return a map with the intensity enhanced above the normalization radius + # and the same meta data as the input map. + new_map = sunpy.map.Map(smap.data * enhancement, smap.meta) new_map.plot_settings["norm"] = None return new_map @@ -221,7 +225,6 @@ def nrgf( width_function=np.std, width_function_kwargs=None, application_radius=1 * u.R_sun, - progress=True, ): """ Implementation of the normalizing radial gradient filter (NRGF). @@ -265,9 +268,6 @@ def nrgf( application_radius : `astropy.units.Quantity`, optional The NRGF is applied to emission at radii above the application_radius. Defaults to 1 solar radii. - progress : `bool`, optional - Display a progressbar on the main loop. - Defaults to True. Returns ------- @@ -317,7 +317,7 @@ def nrgf( data = np.zeros_like(smap.data) # Calculate the filter value for each radial bin. - for i in tqdm(range(radial_bin_edges.shape[1]), desc="NRGF: ", disable=not progress): + for i in tqdm(range(radial_bin_edges.shape[1]), desc="NRGF: "): here = np.logical_and(map_r >= radial_bin_edges[0, i], map_r < radial_bin_edges[1, i]) here = np.logical_and(here, map_r > application_radius) data[here] = smap.data[here] - radial_intensity[i] @@ -351,7 +351,7 @@ def set_attenuation_coefficients(order, range_mean=None, range_std=None, cutoff= .. note:: - The returned maps have their ``plot_settings`` changed to remove the extra normalization step. + The returned maps have their ``plot_settings`` changed to remove the extra normalization step. Parameters ---------- @@ -402,7 +402,6 @@ def fnrgf( width_function=np.std, application_radius=1 * u.R_sun, number_angular_segments=130, - progress=True, ): """ Implementation of the fourier normalizing radial gradient filter (FNRGF). @@ -450,9 +449,6 @@ def fnrgf( number_angular_segments : `int` Number of angular segments in a circular annulus. Defaults to 130. - progress : `bool`, optional - Display a progressbar on the main loop. - Defaults to True. Returns ------- @@ -502,7 +498,7 @@ def fnrgf( data = np.zeros_like(smap.data) # Iterate over each circular ring - for i in tqdm(range(nbins), desc="FNRGF: ", disable=not progress): + for i in tqdm(range(nbins), desc="FNRGF: "): # Finding the pixels which belong to a certain circular ring annulus = np.logical_and(map_r >= radial_bin_edges[0, i], map_r < radial_bin_edges[1, i]) annulus = np.logical_and(annulus, map_r > application_radius) @@ -634,13 +630,14 @@ def _percentile_ranks_numpy_inplace(arr): @u.quantity_input(application_radius=u.R_sun, vignette=u.R_sun) def rhef( smap, + *, radial_bin_edges=None, application_radius=0 * u.R_sun, upsilon=0.35, method="numpy", - *, - vignette=1.5 * u.R_sun, + vignette=None, # 1.5 * u.R_sun, progress=True, + animate=False, ): """ Implementation of the Radial Histogram Equalizing Filter (RHEF). @@ -658,33 +655,28 @@ def rhef( Parameters ---------- smap : `sunpy.map.Map` - The sunpy map to enhance. - radial_bin_edges : `astropy.units.Quantity` - A two-dimensional array of bin edges of size ``[2, nbins]`` where ``nbins`` is - the number of bins. + The SunPy map to enhance using the RHEF algorithm. + radial_bin_edges : `astropy.units.Quantity`, optional + A two-dimensional array of bin edges of size ``[2, nbins]`` where ``nbins`` is the number of bins. + These define the radial segments where filtering is applied. If None, radial bins will be generated automatically. application_radius : `astropy.units.Quantity`, optional - The RHEF is applied to emission at radii above the application_radius. + The radius above which to apply the RHEF. Only regions with radii above this value will be filtered. Defaults to 0 solar radii. - upsilon : None, float, or tuple of `float`, optional - A double-sided gamma function applied to the equalized histograms. - See Equation (4.15) in the thesis. - Defaults to 0.35. - method : str - vignette: `astropy.units.Quantity`, optional - Set pixels above this radius to black. - Defaults to ``1.5*u.R_sun``. - If you want to disable this, pass in None. - Set pixels above this radius to black. - Defaults to None which is no vignette. - One suggested value is ``1.5*u.R_sun``. - progress: bool, optional - Display a progressbar on the main loop. - Defaults to True. + upsilon : float or None, optional + A double-sided gamma function to apply to modify the equalized histograms. Defaults to 0.35. + method : str, optional + Method used to rank the pixels for equalization. Defaults to 'inplace'. + vignette : `astropy.units.Quantity`, optional + Radius beyond which pixels will be set to black. Defaults to ``1.5*u.R_sun``. Set to None to disable vignetting. + progress : bool, optional + If True, display a progress bar during the filtering process. Defaults to True. + animate : bool, optional + If True, display an animated plot during the filtering process. Defaults to False. Returns ------- `sunpy.map.Map` - A SunPy map that has had the RHEF applied to it. + A SunPy map with the Radial Histogram Equalizing Filter applied to it. References ---------- @@ -695,46 +687,69 @@ def rhef( https://www.proquest.com/docview/2759080511 """ - # Get the radii for every pixel - map_r = find_pixel_radii(smap).to(u.R_sun) + # Get the radii for every pixel, ensuring units are correct (in terms of pixels or solar radii) + map_r = find_pixel_radii(smap) + # Automatically generate radial bin edges if none are provided if radial_bin_edges is None: - radial_bin_edges = equally_spaced_bins(0, 2, smap.data.shape[0] // 2) - radial_bin_edges *= u.R_sun + radial_bin_edges = equally_spaced_bins(0, 10 * np.max(map_r.value), smap.data.shape[0] // 2) * u.R_sun - # Make sure bins are in the map. + # Ensure radial_bin_edges are within the bounds of the map_r values if radial_bin_edges[1, -1] > np.max(map_r): - radial_bin_edges = equally_spaced_bins( - inner_value=radial_bin_edges[0, 0], - outer_value=np.max(map_r), - nbins=radial_bin_edges.shape[1], + radial_bin_edges = ( + equally_spaced_bins( + inner_value=radial_bin_edges[0, 0].to(u.R_sun).value, + outer_value=np.max(map_r.to(u.R_sun)).value, + nbins=radial_bin_edges.shape[1], + ) + * u.R_sun ) - # Allocate storage for the filtered data data = np.zeros_like(smap.data) - meta = smap.meta + # Enable progress bar for large numbers of bins if radial_bin_edges.shape[1] > 2000: progress = True - # Calculate the filter values for each radial bin. + # Create figure and axis for animation if enabled + if animate: + fig, ax = plt.subplots(1, 1) + im = ax.imshow(data, origin="lower", cmap="gray", vmin=0, vmax=1) + plt.show(block=False) + + # Select the ranking method + ranking_func = _select_rank_method(method) + + # Loop over each radial bin to apply the filter for i in tqdm(range(radial_bin_edges.shape[1]), desc="RHEF: ", disable=not progress): - # Identify the appropriate radial slice - here = np.logical_and(map_r >= radial_bin_edges[0, i], map_r < radial_bin_edges[1, i]) + # Identify pixels within the current radial bin + here = np.logical_and( + map_r >= radial_bin_edges[0, i].to(u.R_sun), map_r < radial_bin_edges[1, i].to(u.R_sun) + ) if application_radius is not None and application_radius > 0: here = np.logical_and(here, map_r >= application_radius) - - # Perform the filtering operation - ranking_func = _select_rank_method(method) + # Apply ranking function data[here] = ranking_func(smap.data[here]) if upsilon is not None: data[here] = apply_upsilon(data[here], upsilon) - new_map = sunpy.map.Map(data, meta, autoalign=True) + # If animation is enabled, update the image data in the figure + if animate: + im.set_data(data) + plt.pause(0.0001) # Pause to refresh the plot + if animate: + plt.show(block=True) + + # Final update to the map + new_map = Map(data, smap.meta, autoalign=True) + # Apply the vignette if specified if vignette is not None: - new_map = blackout_pixels_above_radius(new_map, vignette) + new_map = blackout_pixels_above_radius(new_map, vignette.to(u.R_sun)) + # Adjust plot settings to remove extra normalization # This must be done whenever one is adjusting the overall statistical distribution of values + new_map.plot_settings["norm"] = None + # Return the new SunPy map with RHEF applied return new_map From 3ac0e98731a11fa830ba1ce4ced9c8cd94bdb1e2 Mon Sep 17 00:00:00 2001 From: gillyspace27 Date: Mon, 23 Sep 2024 14:21:54 -0600 Subject: [PATCH 2/9] Fix call signature in example --- examples/radial_histogram_equalization.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/examples/radial_histogram_equalization.py b/examples/radial_histogram_equalization.py index 0334a66e..d6fa5544 100644 --- a/examples/radial_histogram_equalization.py +++ b/examples/radial_histogram_equalization.py @@ -27,7 +27,7 @@ radial_bin_edges = equally_spaced_bins(0, 2, aia_map.data.shape[0] // 2) radial_bin_edges *= u.R_sun -rhef_map = radial.rhef(aia_map, radial_bin_edges, progress=False) +rhef_map = radial.rhef(aia_map, radial_bin_edges=radial_bin_edges, progress=False) fig = plt.figure(figsize=(15, 10)) From 8f0dffbca83594f87f312362d38c52827284df2d Mon Sep 17 00:00:00 2001 From: gillyspace27 Date: Mon, 23 Sep 2024 14:23:24 -0600 Subject: [PATCH 3/9] apply changes from comments --- sunkit_image/radial.py | 25 +++++-------------------- 1 file changed, 5 insertions(+), 20 deletions(-) diff --git a/sunkit_image/radial.py b/sunkit_image/radial.py index d375a5a9..d9436af6 100644 --- a/sunkit_image/radial.py +++ b/sunkit_image/radial.py @@ -351,7 +351,7 @@ def set_attenuation_coefficients(order, range_mean=None, range_std=None, cutoff= .. note:: - The returned maps have their ``plot_settings`` changed to remove the extra normalization step. + The returned maps have their ``plot_settings`` changed to remove the extra normalization step. Parameters ---------- @@ -637,7 +637,6 @@ def rhef( method="numpy", vignette=None, # 1.5 * u.R_sun, progress=True, - animate=False, ): """ Implementation of the Radial Histogram Equalizing Filter (RHEF). @@ -670,8 +669,6 @@ def rhef( Radius beyond which pixels will be set to black. Defaults to ``1.5*u.R_sun``. Set to None to disable vignetting. progress : bool, optional If True, display a progress bar during the filtering process. Defaults to True. - animate : bool, optional - If True, display an animated plot during the filtering process. Defaults to False. Returns ------- @@ -711,12 +708,6 @@ def rhef( if radial_bin_edges.shape[1] > 2000: progress = True - # Create figure and axis for animation if enabled - if animate: - fig, ax = plt.subplots(1, 1) - im = ax.imshow(data, origin="lower", cmap="gray", vmin=0, vmax=1) - plt.show(block=False) - # Select the ranking method ranking_func = _select_rank_method(method) @@ -733,21 +724,15 @@ def rhef( if upsilon is not None: data[here] = apply_upsilon(data[here], upsilon) - # If animation is enabled, update the image data in the figure - if animate: - im.set_data(data) - plt.pause(0.0001) # Pause to refresh the plot - if animate: - plt.show(block=True) - # Final update to the map - new_map = Map(data, smap.meta, autoalign=True) - # Apply the vignette if specified + new_map = Map(data, smap.meta) + if vignette is not None: new_map = blackout_pixels_above_radius(new_map, vignette.to(u.R_sun)) # Adjust plot settings to remove extra normalization - # This must be done whenever one is adjusting the overall statistical distribution of values + # This must be done whenever one is adjusting + # the overall statistical distribution of values new_map.plot_settings["norm"] = None From e998b6c513d3fc99f3f3db7536b815b849fe433b Mon Sep 17 00:00:00 2001 From: gillyspace27 Date: Mon, 23 Sep 2024 14:23:53 -0600 Subject: [PATCH 4/9] Switch from 0-masking to NAN-masking invalid regions --- sunkit_image/utils/utils.py | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/sunkit_image/utils/utils.py b/sunkit_image/utils/utils.py index 52696379..4a763788 100644 --- a/sunkit_image/utils/utils.py +++ b/sunkit_image/utils/utils.py @@ -194,7 +194,9 @@ def reform2d(array, factor=1): msg = "Input array must be 2d!" raise ValueError(msg) if factor > 1: - congridx = RectBivariateSpline(np.arange(0, array.shape[0]), np.arange(0, array.shape[1]), array, kx=1, ky=1) + congridx = RectBivariateSpline( + np.arange(0, array.shape[0]), np.arange(0, array.shape[1]), array, kx=1, ky=1 + ) return congridx(np.arange(0, array.shape[0], 1 / factor), np.arange(0, array.shape[1], 1 / factor)) return array @@ -392,7 +394,7 @@ def blackout_pixels_above_radius(smap, radius_limit=1.5 * u.R_sun): mask = map_r > radius_limit # Apply the mask to the map data - masked_data = np.where(mask, 0, smap.data) + masked_data = np.where(mask, np.NAN, smap.data) # Create a new map with the masked data return sunpy.map.Map(masked_data, smap.meta) From 41166e913d20bacd29b559076f210e5b5bde9515 Mon Sep 17 00:00:00 2001 From: gillyspace27 Date: Mon, 23 Sep 2024 15:32:20 -0600 Subject: [PATCH 5/9] Respond to PR Review Comments --- sunkit_image/radial.py | 79 +++++++++++++++---------------------- sunkit_image/utils/utils.py | 2 +- 2 files changed, 33 insertions(+), 48 deletions(-) diff --git a/sunkit_image/radial.py b/sunkit_image/radial.py index d9436af6..cb41e965 100644 --- a/sunkit_image/radial.py +++ b/sunkit_image/radial.py @@ -106,7 +106,7 @@ def _normalize_fit_radial_intensity(radii, polynomial, normalization_radius): def intensity_enhance( smap, - radial_bin_edges, + radial_bin_edges=None, scale=None, summarize_bin_edges="center", summary=np.mean, @@ -171,8 +171,8 @@ def intensity_enhance( `sunpy.map.Map` A SunPy map that has the emission above the normalization radius enhanced. """ - # Get the radii for every pixel - map_r = find_pixel_radii(smap).to(u.R_sun) + # Handle the bin edges and radius array + radial_bin_edges, map_r = find_radial_bin_edges(smap, radial_bin_edges) # Get the radial intensity distribution radial_intensity = get_radial_intensity_summary( @@ -218,7 +218,7 @@ def intensity_enhance( def nrgf( smap, - radial_bin_edges, + radial_bin_edges=None, scale=None, intensity_summary=np.nanmean, intensity_summary_kwargs=None, @@ -285,15 +285,9 @@ def nrgf( width_function_kwargs = {} if intensity_summary_kwargs is None: intensity_summary_kwargs = {} - map_r = find_pixel_radii(smap).to(u.R_sun) - # To make sure bins are in the map. - if radial_bin_edges[1, -1] > np.max(map_r): - radial_bin_edges = equally_spaced_bins( - inner_value=radial_bin_edges[0, 0], - outer_value=np.max(map_r), - nbins=radial_bin_edges.shape[1], - ) + # Handle the bin edges and radius array + radial_bin_edges, map_r = find_radial_bin_edges(smap, radial_bin_edges) # Radial intensity radial_intensity = get_radial_intensity_summary( @@ -470,16 +464,8 @@ def fnrgf( msg = "Minimum value of order is 1" raise ValueError(msg) - # Get the radii for every pixel - map_r = find_pixel_radii(smap).to(u.R_sun) - - # To make sure bins are in the map. - if radial_bin_edges[1, -1] > np.max(map_r): - radial_bin_edges = equally_spaced_bins( - inner_value=radial_bin_edges[0, 0], - outer_value=np.max(map_r), - nbins=radial_bin_edges.shape[1], - ) + # Handle the bin edges and radius + radial_bin_edges, map_r = find_radial_bin_edges(smap, radial_bin_edges) # Get the Helioprojective coordinates of each pixel x, y = np.meshgrid(*[np.arange(v.value) for v in smap.dimensions]) * u.pix @@ -627,6 +613,26 @@ def _percentile_ranks_numpy_inplace(arr): return ranking_func +def find_radial_bin_edges(smap, radial_bin_edges=None): + # Get the radii for every pixel, ensuring units are correct (in terms of pixels or solar radii) + map_r = find_pixel_radii(smap) + # Automatically generate radial bin edges if none are provided + if radial_bin_edges is None: + radial_bin_edges = equally_spaced_bins(0, 10 * np.max(map_r.value), smap.data.shape[0] // 2) * u.R_sun + + # Ensure radial_bin_edges are within the bounds of the map_r values + if radial_bin_edges[1, -1] > np.max(map_r): + radial_bin_edges = ( + equally_spaced_bins( + inner_value=radial_bin_edges[0, 0].to(u.R_sun).value, + outer_value=np.max(map_r.to(u.R_sun)).value, + nbins=radial_bin_edges.shape[1] // 2, + ) + * u.R_sun + ) + return radial_bin_edges, map_r + + @u.quantity_input(application_radius=u.R_sun, vignette=u.R_sun) def rhef( smap, @@ -635,7 +641,7 @@ def rhef( application_radius=0 * u.R_sun, upsilon=0.35, method="numpy", - vignette=None, # 1.5 * u.R_sun, + vignette=None, progress=True, ): """ @@ -664,9 +670,9 @@ def rhef( upsilon : float or None, optional A double-sided gamma function to apply to modify the equalized histograms. Defaults to 0.35. method : str, optional - Method used to rank the pixels for equalization. Defaults to 'inplace'. + Method used to rank the pixels for equalization. Defaults to 'inplace', with 'scipy' and 'numpy' as other options. vignette : `astropy.units.Quantity`, optional - Radius beyond which pixels will be set to black. Defaults to ``1.5*u.R_sun``. Set to None to disable vignetting. + Radius beyond which pixels will be set to NAN. Defaults to None, must be in units that are compatible with "R_sun" as the value will be transformed. progress : bool, optional If True, display a progress bar during the filtering process. Defaults to True. @@ -684,30 +690,10 @@ def rhef( https://www.proquest.com/docview/2759080511 """ - # Get the radii for every pixel, ensuring units are correct (in terms of pixels or solar radii) - map_r = find_pixel_radii(smap) - - # Automatically generate radial bin edges if none are provided - if radial_bin_edges is None: - radial_bin_edges = equally_spaced_bins(0, 10 * np.max(map_r.value), smap.data.shape[0] // 2) * u.R_sun + radial_bin_edges, map_r = find_radial_bin_edges(smap, radial_bin_edges) - # Ensure radial_bin_edges are within the bounds of the map_r values - if radial_bin_edges[1, -1] > np.max(map_r): - radial_bin_edges = ( - equally_spaced_bins( - inner_value=radial_bin_edges[0, 0].to(u.R_sun).value, - outer_value=np.max(map_r.to(u.R_sun)).value, - nbins=radial_bin_edges.shape[1], - ) - * u.R_sun - ) - # Allocate storage for the filtered data data = np.zeros_like(smap.data) - # Enable progress bar for large numbers of bins - if radial_bin_edges.shape[1] > 2000: - progress = True - # Select the ranking method ranking_func = _select_rank_method(method) @@ -724,7 +710,6 @@ def rhef( if upsilon is not None: data[here] = apply_upsilon(data[here], upsilon) - # Final update to the map new_map = Map(data, smap.meta) if vignette is not None: diff --git a/sunkit_image/utils/utils.py b/sunkit_image/utils/utils.py index 4a763788..ae7db383 100644 --- a/sunkit_image/utils/utils.py +++ b/sunkit_image/utils/utils.py @@ -394,7 +394,7 @@ def blackout_pixels_above_radius(smap, radius_limit=1.5 * u.R_sun): mask = map_r > radius_limit # Apply the mask to the map data - masked_data = np.where(mask, np.NAN, smap.data) + masked_data = np.where(mask, np.nan, smap.data) # Create a new map with the masked data return sunpy.map.Map(masked_data, smap.meta) From 2f4c5ea9abf2a135d3fa1d23d7751b0b8eb105c9 Mon Sep 17 00:00:00 2001 From: gillyspace27 Date: Mon, 23 Sep 2024 16:00:44 -0600 Subject: [PATCH 6/9] remove matplotlib --- sunkit_image/radial.py | 1 - 1 file changed, 1 deletion(-) diff --git a/sunkit_image/radial.py b/sunkit_image/radial.py index cb41e965..3b45ef3c 100644 --- a/sunkit_image/radial.py +++ b/sunkit_image/radial.py @@ -3,7 +3,6 @@ radius. """ -import matplotlib.pyplot as plt import numpy as np from tqdm import tqdm From 4533c896be2e2d43175632676bce0342d896c7ae Mon Sep 17 00:00:00 2001 From: gillyspace27 Date: Mon, 23 Sep 2024 21:17:49 -0600 Subject: [PATCH 7/9] fix double calc of radial bin edges, replacement of progressbars --- sunkit_image/radial.py | 16 +++++++++++----- sunkit_image/utils/utils.py | 6 ++++-- 2 files changed, 15 insertions(+), 7 deletions(-) diff --git a/sunkit_image/radial.py b/sunkit_image/radial.py index 3b45ef3c..6ea21688 100644 --- a/sunkit_image/radial.py +++ b/sunkit_image/radial.py @@ -10,7 +10,7 @@ import sunpy.map from sunpy.coordinates import frames -from sunpy.map import Map + from sunkit_image.utils import ( apply_upsilon, @@ -224,6 +224,7 @@ def nrgf( width_function=np.std, width_function_kwargs=None, application_radius=1 * u.R_sun, + progress=True, ): """ Implementation of the normalizing radial gradient filter (NRGF). @@ -267,6 +268,8 @@ def nrgf( application_radius : `astropy.units.Quantity`, optional The NRGF is applied to emission at radii above the application_radius. Defaults to 1 solar radii. + progress : ``bool``, optional + Show a progressbar while computing Returns ------- @@ -310,7 +313,7 @@ def nrgf( data = np.zeros_like(smap.data) # Calculate the filter value for each radial bin. - for i in tqdm(range(radial_bin_edges.shape[1]), desc="NRGF: "): + for i in tqdm(range(radial_bin_edges.shape[1]), desc="NRGF: ", disable=not progress): here = np.logical_and(map_r >= radial_bin_edges[0, i], map_r < radial_bin_edges[1, i]) here = np.logical_and(here, map_r > application_radius) data[here] = smap.data[here] - radial_intensity[i] @@ -395,6 +398,7 @@ def fnrgf( width_function=np.std, application_radius=1 * u.R_sun, number_angular_segments=130, + progress=True, ): """ Implementation of the fourier normalizing radial gradient filter (FNRGF). @@ -442,6 +446,8 @@ def fnrgf( number_angular_segments : `int` Number of angular segments in a circular annulus. Defaults to 130. + progress : ``bool``, optional + Show a progressbar while computing Returns ------- @@ -483,7 +489,7 @@ def fnrgf( data = np.zeros_like(smap.data) # Iterate over each circular ring - for i in tqdm(range(nbins), desc="FNRGF: "): + for i in tqdm(range(nbins), desc="FNRGF: ", disable=not progress): # Finding the pixels which belong to a certain circular ring annulus = np.logical_and(map_r >= radial_bin_edges[0, i], map_r < radial_bin_edges[1, i]) annulus = np.logical_and(annulus, map_r > application_radius) @@ -617,7 +623,7 @@ def find_radial_bin_edges(smap, radial_bin_edges=None): map_r = find_pixel_radii(smap) # Automatically generate radial bin edges if none are provided if radial_bin_edges is None: - radial_bin_edges = equally_spaced_bins(0, 10 * np.max(map_r.value), smap.data.shape[0] // 2) * u.R_sun + radial_bin_edges = equally_spaced_bins(0, np.max(map_r.value), smap.data.shape[0] // 2) * u.R_sun # Ensure radial_bin_edges are within the bounds of the map_r values if radial_bin_edges[1, -1] > np.max(map_r): @@ -709,7 +715,7 @@ def rhef( if upsilon is not None: data[here] = apply_upsilon(data[here], upsilon) - new_map = Map(data, smap.meta) + new_map = sunpy.map.Map(data, smap.meta) if vignette is not None: new_map = blackout_pixels_above_radius(new_map, vignette.to(u.R_sun)) diff --git a/sunkit_image/utils/utils.py b/sunkit_image/utils/utils.py index ae7db383..bc506675 100644 --- a/sunkit_image/utils/utils.py +++ b/sunkit_image/utils/utils.py @@ -371,7 +371,7 @@ def apply_upsilon(data, upsilon=(0.5, 0.5)): return out_curve -def blackout_pixels_above_radius(smap, radius_limit=1.5 * u.R_sun): +def blackout_pixels_above_radius(smap, radius_limit=1.5 * u.R_sun, fill=np.nan): """ Black out any pixels above a certain radius in a SunPy map. @@ -381,6 +381,8 @@ def blackout_pixels_above_radius(smap, radius_limit=1.5 * u.R_sun): The input sunpy map. radius_limit : `astropy.units.Quantity` The radius limit above which to black out pixels. + fill : `any` + The value to use above the radius_limit. Returns ------- @@ -394,7 +396,7 @@ def blackout_pixels_above_radius(smap, radius_limit=1.5 * u.R_sun): mask = map_r > radius_limit # Apply the mask to the map data - masked_data = np.where(mask, np.nan, smap.data) + masked_data = np.where(mask, fill, smap.data) # Create a new map with the masked data return sunpy.map.Map(masked_data, smap.meta) From f361f6fac1430258668dc62c833e63be2aa89b03 Mon Sep 17 00:00:00 2001 From: gillyspace27 Date: Mon, 23 Sep 2024 21:25:27 -0600 Subject: [PATCH 8/9] fix typo --- changelog/205.doc.rst | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/changelog/205.doc.rst b/changelog/205.doc.rst index f46453c4..3e85ba65 100644 --- a/changelog/205.doc.rst +++ b/changelog/205.doc.rst @@ -1,2 +1,2 @@ Added a new example for `sunkit_image.radial.rhef` at :ref:`sphx_glr_generated_gallery_radial_histogram_equalization.py` -Added RFEF to :ref:`sphx_glr_generated_gallery_radial_gradient_filters.py` as a comparison to the other filters. +Added RHEF to :ref:`sphx_glr_generated_gallery_radial_gradient_filters.py` as a comparison to the other filters. From 0b558645ca37978ef4df7abfeb6c834125f41a39 Mon Sep 17 00:00:00 2001 From: gillyspace27 Date: Mon, 23 Sep 2024 22:29:18 -0600 Subject: [PATCH 9/9] default to nan outside of bounds, or allow user to select value --- sunkit_image/radial.py | 19 +++++++++++++++---- 1 file changed, 15 insertions(+), 4 deletions(-) diff --git a/sunkit_image/radial.py b/sunkit_image/radial.py index 6ea21688..1cd3103d 100644 --- a/sunkit_image/radial.py +++ b/sunkit_image/radial.py @@ -225,6 +225,7 @@ def nrgf( width_function_kwargs=None, application_radius=1 * u.R_sun, progress=True, + fill=np.nan, ): """ Implementation of the normalizing radial gradient filter (NRGF). @@ -270,6 +271,9 @@ def nrgf( Defaults to 1 solar radii. progress : ``bool``, optional Show a progressbar while computing + fill : ``any``, optional + The value to be placed outside of the bounds of the algorithm + Defaults to NAN. Returns ------- @@ -310,7 +314,7 @@ def nrgf( ) # Storage for the filtered data - data = np.zeros_like(smap.data) + data = np.ones_like(smap.data) * fill # Calculate the filter value for each radial bin. for i in tqdm(range(radial_bin_edges.shape[1]), desc="NRGF: ", disable=not progress): @@ -399,6 +403,7 @@ def fnrgf( application_radius=1 * u.R_sun, number_angular_segments=130, progress=True, + fill=np.nan, ): """ Implementation of the fourier normalizing radial gradient filter (FNRGF). @@ -448,6 +453,9 @@ def fnrgf( Defaults to 130. progress : ``bool``, optional Show a progressbar while computing + fill : ``any``, optional + The value to be placed outside of the bounds of the algorithm + Defaults to NAN. Returns ------- @@ -486,7 +494,7 @@ def fnrgf( nbins = radial_bin_edges.shape[1] # Storage for the filtered data - data = np.zeros_like(smap.data) + data = np.ones_like(smap.data) * fill # Iterate over each circular ring for i in tqdm(range(nbins), desc="FNRGF: ", disable=not progress): @@ -648,6 +656,7 @@ def rhef( method="numpy", vignette=None, progress=True, + fill=np.nan, ): """ Implementation of the Radial Histogram Equalizing Filter (RHEF). @@ -680,7 +689,9 @@ def rhef( Radius beyond which pixels will be set to NAN. Defaults to None, must be in units that are compatible with "R_sun" as the value will be transformed. progress : bool, optional If True, display a progress bar during the filtering process. Defaults to True. - + fill : ``any``, optional + The value to be placed outside of the bounds of the algorithm + Defaults to NAN. Returns ------- `sunpy.map.Map` @@ -697,7 +708,7 @@ def rhef( radial_bin_edges, map_r = find_radial_bin_edges(smap, radial_bin_edges) - data = np.zeros_like(smap.data) + data = np.ones_like(smap.data) * fill # Select the ranking method ranking_func = _select_rank_method(method)