From ec3dc5691980cfba1b4bbad2a168f949cbde9741 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Sat, 13 Jul 2024 19:30:23 +0000 Subject: [PATCH] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- chimerapy/chimera.py | 3 +- chimerapy/chimera_copy.py | 114 ++++++++++++++++++-------------- chimerapy/chimera_legacy.py | 5 +- chimerapy/tests/test_chimera.py | 81 ++++++----------------- 4 files changed, 86 insertions(+), 117 deletions(-) diff --git a/chimerapy/chimera.py b/chimerapy/chimera.py index d5a80fe..70730e6 100644 --- a/chimerapy/chimera.py +++ b/chimerapy/chimera.py @@ -1,4 +1,5 @@ """Package for Coronal Hole Identification Algorithm""" + import glob import sys @@ -209,7 +210,7 @@ def to_helio(image: np.array): hg = hpc.transform_to(sunpy.coordinates.frames.HeliographicStonyhurst) global csys csys = wcs.WCS(hedb) - + to_helio(im171) diff --git a/chimerapy/chimera_copy.py b/chimerapy/chimera_copy.py index 68b0966..a5634fc 100644 --- a/chimerapy/chimera_copy.py +++ b/chimerapy/chimera_copy.py @@ -1,30 +1,21 @@ """Package for Coronal Hole Identification Algorithm""" -import glob -import sys import astropy.units as u import cv2 import mahotas import matplotlib.pyplot as plt import numpy as np -import scipy -from sunpy.map import Map -from sunpy.coordinates import frames -from sunpy.coordinates import propagate_with_solar_surface -import scipy.interpolate import sunpy import sunpy.map from astropy import wcs -from astropy.io import fits from astropy.modeling.models import Gaussian2D from astropy.visualization import astropy_mpl_style from skimage.util import img_as_ubyte -from astropy.wcs import WCS -from sunpy.map import all_coordinates_from_map -from sunpy.coordinates import HeliographicStonyhurst +from sunpy.coordinates import (HeliographicStonyhurst, + propagate_with_solar_surface) +from sunpy.map import Map, all_coordinates_from_map - -#--noverify +# --noverify plt.style.use(astropy_mpl_style) @@ -32,16 +23,17 @@ file_path = "./" -INPUT_FILES = {"aia171": 'http://jsoc.stanford.edu/data/aia/synoptic/2016/09/22/H1000/AIA20160922_1030_0171.fits', -"aia193": 'http://jsoc.stanford.edu/data/aia/synoptic/2016/09/22/H1000/AIA20160922_1030_0193.fits', -"aia211": 'http://jsoc.stanford.edu/data/aia/synoptic/2016/09/22/H1000/AIA20160922_1030_0211.fits', -"hmi_mag": 'http://jsoc.stanford.edu/data/hmi/fits/2016/09/22/hmi.M_720s.20160922_010000_TAI.fits', +INPUT_FILES = { + "aia171": "http://jsoc.stanford.edu/data/aia/synoptic/2016/09/22/H1000/AIA20160922_1030_0171.fits", + "aia193": "http://jsoc.stanford.edu/data/aia/synoptic/2016/09/22/H1000/AIA20160922_1030_0193.fits", + "aia211": "http://jsoc.stanford.edu/data/aia/synoptic/2016/09/22/H1000/AIA20160922_1030_0211.fits", + "hmi_mag": "http://jsoc.stanford.edu/data/hmi/fits/2016/09/22/hmi.M_720s.20160922_010000_TAI.fits", } -im171 = Map(INPUT_FILES['aia171']) -im193 = Map(INPUT_FILES['aia193']) -im211 = Map(INPUT_FILES['aia211']) -imhmi = Map(INPUT_FILES['hmi_mag']) +im171 = Map(INPUT_FILES["aia171"]) +im193 = Map(INPUT_FILES["aia193"]) +im211 = Map(INPUT_FILES["aia211"]) +imhmi = Map(INPUT_FILES["hmi_mag"]) def rescale(proj_to: sunpy.map.Map, input_map: sunpy.map.Map): @@ -56,23 +48,25 @@ def rescale(proj_to: sunpy.map.Map, input_map: sunpy.map.Map): Returns ------- array: 'np.array' - - """ + + """ with propagate_with_solar_surface(): map1 = proj_to.reproject_to(input_map.wcs) new_x_scale = map1.scale[0].to(u.arcsec / u.pixel).value new_y_scale = map1.scale[1].to(u.arcsec / u.pixel).value - map1.meta['cdelt1'] = new_x_scale - map1.meta['cdelt2'] = new_y_scale - map1.meta['cunit1'] = 'arcsec' - map1.meta['cunit2'] = 'arcsec' + map1.meta["cdelt1"] = new_x_scale + map1.meta["cdelt2"] = new_y_scale + map1.meta["cunit1"] = "arcsec" + map1.meta["cunit2"] = "arcsec" return map1 + im171 = rescale(im171, im171) im193 = rescale(im171, im193) im211 = rescale(im171, im211) imhmi = rescale(im171, imhmi) + def filter(map1: np.array, map2: np.array, map3: np.array): """ Defines headers and filters aia arrays to meet header requirements. Removes 0 values from each array. @@ -96,18 +90,18 @@ def filter(map1: np.array, map2: np.array, map3: np.array): return map1, map2, map3 + im171, im193, im211 = filter(im171, im193, im211) """defines the shape of the arrays as "s" and "rs" as the solar radius""" s = np.shape(im171.data) -#do we want the solar radius in arcsec or pixels? +# do we want the solar radius in arcsec or pixels? rs = im171.rsun_obs -rs_pixels = im171.rsun_obs/im171.scale[0] - +rs_pixels = im171.rsun_obs / im171.scale[0] def pix_arc(map: sunpy.map.Map): - ''' + """ Calculates the conversion value of pixel to arcsec Parameters @@ -117,7 +111,7 @@ def pix_arc(map: sunpy.map.Map): Returns - ''' + """ dattoarc = map.scale[0].value conver = ((s[0]) / 2) * dattoarc / map.meta["cdelt1"] - (s[1] / 2) convermul = dattoarc / map.meta["cdelt1"] @@ -149,6 +143,7 @@ def to_helio(map: sunpy.map.Map): csys = wcs.WCS(dict(map.meta)) return hpc, hg, csys + hpc, hg, csys = to_helio(im171) """Setting up arrays to be used in later processing""" @@ -166,7 +161,13 @@ def to_helio(map: sunpy.map.Map): w = np.where((xgrid - center[0]) ** 2 + (ygrid - center[1]) ** 2 > r**2) y, x = np.mgrid[0:1024, 0:1024] pix_size = (2000 * u.arcsec).value -garr = Gaussian2D(1, im171.reference_pixel.x.value, im171.reference_pixel.y.value, pix_size / im171.scale[0].value, pix_size / im171.scale[0].value)(x, y) +garr = Gaussian2D( + 1, + im171.reference_pixel.x.value, + im171.reference_pixel.y.value, + pix_size / im171.scale[0].value, + pix_size / im171.scale[0].value, +)(x, y) garr[w] = 1.0 """creates sub-arrays of props to isolate column of index 0 and column of index 1""" @@ -229,8 +230,10 @@ def to_helio(map: sunpy.map.Map): ) """define threshold values in log space""" + + def log_dat(map1: sunpy.map.Map, map2: sunpy.map.Map, map3: sunpy.map.Map): - ''' + """ Takes the log base-10 of all sunpy map data Parameters @@ -244,13 +247,14 @@ def log_dat(map1: sunpy.map.Map, map2: sunpy.map.Map, map3: sunpy.map.Map): t0: 'np.array' t1: 'np.array' t2: 'np.array' - ''' + """ with np.errstate(divide="ignore"): t0 = np.log10(map1.data) t1 = np.log10(map2.data) t2 = np.log10(map3.data) return t0, t1, t2 + t0, t1, t2 = log_dat(im171, im193, im211) @@ -276,7 +280,8 @@ def new_s(self, new_slope): t1b = Bounds(1.4, 3.0, 255) t2b = Bounds(1.2, 3.9, 255) -#set to also take in boundaries + +# set to also take in boundaries def set_contour(t0: np.array, t1: np.array, t2: np.array): """ Threshold arrays based on desired boundaries and sets contours. @@ -295,15 +300,15 @@ def set_contour(t0: np.array, t1: np.array, t2: np.array): """ if t0 is not None and t1 is not None and t2 is not None: - #set the threshold and contours for t0 + # set the threshold and contours for t0 t0[np.where(t0 < t0b.upper)] = t0b.upper t0[np.where(t0 > t0b.lower)] = t0b.lower t0 = np.array(((t0 - t0b.upper) / (t0b.lower - t0b.upper)) * t0b.slope, dtype=np.float32) - #set the threshold and contours for t1 + # set the threshold and contours for t1 t1[np.where(t1 < t1b.upper)] = t1b.upper t1[np.where(t1 > t1b.lower)] = t2b.lower t1 = np.array(((t1 - t1b.upper) / (t1b.lower - t1b.upper)) * t1b.slope, dtype=np.float32) - #set the threshold and contours for t2 + # set the threshold and contours for t2 t2[np.where(t2 < t2b.upper)] = t2b.upper t2[np.where(t2 > t2b.lower)] = t2b.lower t2 = np.array(((t2 - t2b.upper) / (t2b.lower - t2b.upper)) * t2b.slope, dtype=np.float32) @@ -314,7 +319,10 @@ def set_contour(t0: np.array, t1: np.array, t2: np.array): t0, t1, t2 = set_contour(t0, t1, t2) -def create_mask(tm1: np.array, tm2: np.array, tm3: np.array, map1: sunpy.map.Map, map2: sunpy.map.Map, map3: sunpy.map.Map): + +def create_mask( + tm1: np.array, tm2: np.array, tm3: np.array, map1: sunpy.map.Map, map2: sunpy.map.Map, map3: sunpy.map.Map +): """ Creates 3 segmented bitmasks @@ -326,7 +334,7 @@ def create_mask(tm1: np.array, tm2: np.array, tm3: np.array, map1: sunpy.map.Map map1: 'sunpy.map.Map' map2: 'sunpy.map.Map' map3: 'sunpy.map.Map' - + Returns ------- bmmix: 'np.array' @@ -343,9 +351,10 @@ def create_mask(tm1: np.array, tm2: np.array, tm3: np.array, map1: sunpy.map.Map bmmix, bmhot, bmcool = create_mask(t0, t1, t2, im171, im193, im211) -#conjunction of 3 bitmasks +# conjunction of 3 bitmasks cand = bmcool * bmmix * bmhot + def misid(can: np.array, cir: np.array, xgir: np.array, ygir: np.array, thresh_rad: int): """ Removes off-detector mis-identification @@ -362,7 +371,7 @@ def misid(can: np.array, cir: np.array, xgir: np.array, ygir: np.array, thresh_r 'np.array' """ - #make r a function argument, give name and unit + # make r a function argument, give name and unit r = thresh_rad w = np.where((xgir - center[0]) ** 2 + (ygir - center[1]) ** 2 <= thresh_rad**2) cir[w] = 1.0 @@ -396,16 +405,18 @@ def on_off(cir: np.array, can: np.array): cir[outside] = 1.0 can = can * cir plt.figure() - plt.imshow(cand, cmap='viridis') + plt.imshow(cand, cmap="viridis") plt.show return can + cand = on_off(circ, cand) + def contour_data(cand: np.array): """ Contours the identified datapoints - + Parameters ---------- cand: 'np.array' @@ -419,15 +430,17 @@ def contour_data(cand: np.array): """ cand = np.array(cand, dtype=np.uint8) cont, heir = cv2.findContours(cand, cv2.RETR_TREE, cv2.CHAIN_APPROX_SIMPLE) - #I think cont might be the x-y coordinates in pixels? + # I think cont might be the x-y coordinates in pixels? return cand, cont, heir + cand, cont, heir = contour_data(cand) + def sort(cont: tuple): """ Sorts the contours by size - + Parameters ---------- cont: 'tuple' @@ -476,7 +489,7 @@ def extent(map: sunpy.map.Map, cont: tuple): """ coord_hpc = map.world2pix(cont) - maxlat = coord_hpc.transform_to(HeliographicStonyhurst).lat.max() + maxlat = coord_hpc.transform_to(HeliographicStonyhurst).lat.max() maxlon = coord_hpc.transform_to(HeliographicStonyhurst).lon.max() minlat = coord_hpc.transform_to(HeliographicStonyhurst).lat.min() minlon = coord_hpc.transform_to(HeliographicStonyhurst).lat.min() @@ -488,7 +501,6 @@ def extent(map: sunpy.map.Map, cont: tuple): return maxlat, maxlon, minlat, minlon, centlat, centlon - def coords(i, csys, cont): """ Finds coordinates of CH boundaries @@ -921,4 +933,6 @@ def plot_mask(slate=slate): plot_mask() if __name__ == "__main__": - import_functions(INPUT_FILES['aia171'], INPUT_FILES['aia193'], INPUT_FILES['aia211'], INPUT_FILES['hmi_mag']) \ No newline at end of file + import_functions( + INPUT_FILES["aia171"], INPUT_FILES["aia193"], INPUT_FILES["aia211"], INPUT_FILES["hmi_mag"] + ) diff --git a/chimerapy/chimera_legacy.py b/chimerapy/chimera_legacy.py index 4462553..8d1a828 100644 --- a/chimerapy/chimera_legacy.py +++ b/chimerapy/chimera_legacy.py @@ -1,7 +1,4 @@ -""" - -""" - +""" """ import glob import sys diff --git a/chimerapy/tests/test_chimera.py b/chimerapy/tests/test_chimera.py index c3da28b..714c112 100644 --- a/chimerapy/tests/test_chimera.py +++ b/chimerapy/tests/test_chimera.py @@ -1,4 +1,3 @@ -import glob import os import warnings @@ -7,64 +6,24 @@ from astropy.io import fits from astropy.modeling.models import Gaussian2D -from chimerapy.chimera import ( - Bounds, - Xeb, - Xnb, - Xsb, - Xwb, - Yeb, - Ynb, - Ysb, - Ywb, - ang, - arcar, - arccent, - area, - cent, - centlat, - centlon, - chpts, - cont, - coords, - csys, - data, - datb, - datc, - datm, - dist, - eastl, - extent, - filter, - hg, - ins_prop, - mB, - mBneg, - mBpos, - npix, - pos, - remove_neg, - rescale_aia, - rescale_hmi, - set_contour, - sort, - threshold, - truarcar, - trummar, - trupixar, - westl, - width, - xpos, - ypos, -) from chimerapy.chimera import * - -INPUT_FILES = {"aia171": "https://solarmonitor.org/data/2016/09/22/fits/saia/saia_00171_fd_20160922_103010.fts.gz", -"aia193": "https://solarmonitor.org/data/2016/09/22/fits/saia/saia_00193_fd_20160922_103041.fts.gz", -"aia211": "https://solarmonitor.org/data/2016/09/22/fits/saia/saia_00211_fd_20160922_103046.fts.gz", -"hmi_mag": "https://solarmonitor.org/data/2016/09/22/fits/shmi/shmi_maglc_fd_20160922_094640.fts.gz", +from chimerapy.chimera import (Bounds, Xeb, Xnb, Xsb, Xwb, Yeb, Ynb, Ysb, Ywb, + ang, arcar, arccent, area, cent, centlat, + centlon, chpts, cont, coords, csys, data, datb, + datc, datm, dist, eastl, extent, filter, hg, + ins_prop, mB, mBneg, mBpos, npix, pos, + remove_neg, rescale_aia, rescale_hmi, + set_contour, sort, threshold, truarcar, trummar, + trupixar, westl, width, xpos, ypos) + +INPUT_FILES = { + "aia171": "https://solarmonitor.org/data/2016/09/22/fits/saia/saia_00171_fd_20160922_103010.fts.gz", + "aia193": "https://solarmonitor.org/data/2016/09/22/fits/saia/saia_00193_fd_20160922_103041.fts.gz", + "aia211": "https://solarmonitor.org/data/2016/09/22/fits/saia/saia_00211_fd_20160922_103046.fts.gz", + "hmi_mag": "https://solarmonitor.org/data/2016/09/22/fits/shmi/shmi_maglc_fd_20160922_094640.fts.gz", } + def img_present(): assert im171 != [] or im193 != [] or im211 != [] or imhmi != [], "Not all required files present" @@ -365,7 +324,7 @@ def test_fill_polygon(): iarr = np.array([[0, 0, 0, 0, 0], [0, 0, 0, 0, 0], [0, 0, 0, 1, 0], [0, 0, 0, 0, 0], [0, 0, 0, 0, 0]]) - cont = np.array([[[[1, 2], [1, 3], [2, 3], [2, 2]]]]) + cont = np.array([[[[1, 2], [1, 3], [2, 3], [2, 2]]]]) slate = np.array(iarr) @@ -384,14 +343,12 @@ def test_fill_polygon(): ] > 0 ): - mahotas.polygon.fill_polygon(polygon_vertices, slate) - print("After filling polygon:") print(slate) iarr[np.where(slate == 1)] = 0 - slate[:] = 0 + slate[:] = 0 assert np.array_equal( slate, @@ -421,8 +378,8 @@ def test_magpol(): assert npix[0][np.where(npix[0] == 0)] != 0, "Npix[0] should not be equal to zero at its zeros" assert npix[1] != 0, "Npix[1] should not be equal to 0" npixtest = [ - np.array([2, -1, 0, 3, -2, 1]), - np.array([1, -1, 0, 2, -2, 1]), + np.array([2, -1, 0, 3, -2, 1]), + np.array([1, -1, 0, 2, -2, 1]), ] wh1_expected = np.where(npixtest[1] > 0) wh1_actual = np.where(npixtest[1] > 0)