ADVANCED WARNING FOR SOLAR FLARES FROM PHOTOSPHERE IMAGE ANALYSIS

Abstract

A method for quantifying disorder and extracting a corresponding numerical value of an order parameter from contrast analysis applied to optical images acquired of the solar photosphere. Temporal variation of the order parameter may be utilized to predict events such as solar flares, which have the ability to disrupt both communication systems and satellite orbits. The degree of order of the photosphere may be monitored to predict solar flares and other solar events. The method may utilize a spin-based (Ising/Potts) model of disorder.

Claims

1. A method of predicting solar events, the method comprising: repeatedly extracting numerical values corresponding to an order parameter from a plurality of sets of solar atmosphere data measured over time; predicting a solar event if a change in the extracted numerical values satisfies prediction criteria.

2. The method of claim 1, wherein: the prediction criteria comprises a predefined increase in numerical values corresponding to the order parameter.

3. The method of claim 2, wherein: the prediction criteria comprises a predefined increase in numerical values corresponding to the order parameter over a predefined interval of time.

4. The method of claim 3, wherein: the plurality of sets of solar atmosphere data comprise images of the Sun; and the numerical values corresponding to order parameter are extracted from a series of images of the Sun.

5. The method of claim 4, wherein: at least some of the images of the Sun include bright regions and dark regions; and including: extracting the numerical values corresponding to order parameter comprises thresholding the images of the Sun to determine an area of the bright regions; and including: selecting regions of interest of the images of the Sun; and determining a squared order parameter value (S.sup.2) by calculating a ratio of the area of the bright regions to a total area of the regions of interest of each image of the Sun.

6. The method of claim 5, wherein: thresholding the images of the Sun includes converting the images to greyscale images; and including: calculating a pixel intensity histogram for each of the greyscale images; fitting the pixel intensity histogram with first and second curves that sum to an overall curve fit of the pixel intensity histogram, wherein the first and second curves represent pixel intensity distributions from disordered and ordered regions, respectively, and wherein a higher one of the first and second curves has a peak that is higher than the other of the first and second curves; determining an intersection of the first and second curves; determining a distance comprising a number of standard deviations that the intersection is away from the peak of the higher one of the first and second curves; setting a threshold to a value of a highest center point value of the higher one of the first and second curves less the floor of that number of standard deviations; performing a binary threshold on each greyscale image to form a binary image; and calculating an order parameter squared (S.sup.2) value for each greyscale image by counting the number of bright pixels in the binary image and dividing the number of bright pixels in the binary image by the total number of pixels contained in each greyscale image.

7. The method of claim 6, including: selecting a region of interest prior to calculating a pixel intensity histogram; and wherein: the pixel intensity histogram is calculated for the region of interest; and the S.sup.2 value is calculated for the region of interest.

8. The method of claim 3, wherein: the predefined increase in numerical values corresponding to order parameter comprises an increase of at least 5% over a period of 30 minutes or less.

9. The method of claim 3, wherein: the predefined increase in numerical values corresponding to order parameter comprises an increase of at least 0.04.

10. The method of claim 1, wherein: the numerical values corresponding to order parameter comprise a squared order parameter (S.sup.2).

11. The method of claim 10, including: fitting a curve to the S.sup.2 values, wherein the curve comprises S.sup.2 values over time; and wherein: the prediction criteria comprises a slope in the curve that meets or exceeds a predefined slope.

12. The method of claim 1, wherein: the photosphere data comprises frames from a solar photosphere movie.

13. The method of claim 1, wherein: the image data comprises frames from a movie comprising a sequence of images of a specific layer of the solar atmosphere.

14. The method of claim 1, including: utilizing a computer to repeatedly extract numerical values corresponding to an order parameter from a plurality of sets of solar atmosphere data measured over time; utilizing a computer to predict a solar event if a change in the extracted numerical values satisfies prediction criteria.

15. A computer-implemented method of predicting solar events, the method comprising: utilizing a computer to repeatedly extract numerical values corresponding to an order parameter from sets of digital solar atmosphere data measured over time; utilizing a computer to predict a solar event if a change in the extracted numerical values satisfies prediction criteria.

16. The method of claim 15, wherein: the prediction criteria comprises a predefined increase in numerical values corresponding to the order parameter.

17. The method of claim 16, wherein: the prediction criteria comprises a predefined increase in numerical values corresponding to the order parameter over a predefined interval of time.

18. The method of claim 3, wherein: the sets of solar atmosphere data comprise digital images of the Sun; and the numerical values corresponding to order parameter are extracted from a series of digital images of the Sun utilizing a computer.

19. The method of claim 4, wherein: at least some of the digital images of the Sun include bright regions and dark regions; and including: extracting the numerical values corresponding to order parameter comprises utilizing a computer to threshold the digital images of the Sun to determine an area of the bright regions; and including: selecting regions of interest of the digital images of the Sun; and utilizing a computer to determine a squared order parameter value (S.sup.2) by calculating, using a computer, a ratio of the area of the bright regions to a total area of the regions of interest of each digital image of the Sun.

20. The method of claim 19, wherein: thresholding the digital images of the Sun includes utilizing a computer to convert the digital images to greyscale digital images; and including: utilizing a computer to calculate a pixel intensity histogram for each of the greyscale images; utilizing a computer to fit the pixel intensity histogram with first and second curves that sum to an overall curve fit of the pixel intensity histogram, wherein the first and second curves represent pixel intensity distributions from disordered and ordered regions, respectively, and wherein a higher one of the first and second curves has a peak that is higher than the other of the first and second curves; utilizing a computer to determine an intersection of the first and second curves; utilizing a computer to determine a distance comprising a number of standard deviations that the intersection is away from the peak of the higher one of the first and second curves; utilizing a computer to set a threshold to a value of a highest center point value of the higher one of the first and second curves less the floor of that number of standard deviations; utilizing a computer to perform a binary threshold on each greyscale digital image to form a binary digital image; and calculating an order parameter squared (S.sup.2) value for each greyscale digital image by counting the number of bright pixels in the binary digital image and dividing the number of bright pixels by the total number of pixels contained in each greyscale digital image.

Description

BRIEF DESCRIPTION OF THE DRAWINGS

[0011] The patent or application file contains at least one drawing executed in color. Copies of this patent or patent application publication with color drawing(s) will be provided by the Office upon request and payment of the necessary fee.

[0012] FIG. 1 is a graph showing the squared order parameter (S.sup.2) of the photosphere extracted from images captured by DKIST at a wavelength of 705 nm;

[0013] FIG. 2 is a graph of the squared order parameter (S.sup.2) as a function of time (years), wherein the data was obtained by the Atmospheric Imaging Assembly (AIA) instrument in NASA's Solar Dynamics Observatory (SDO) at a wavelength of 17.1 nm;

[0014] FIG. 3A is a representative still frame from a solar photosphere movie acquired at 705 nm by the Inouye telescope;

[0015] FIG. 3B is a pixel intensity histogram, showing the fitting of the overall intensity distribution and the decomposition into two curves, corresponding to the disordered and ordered components of the structure, respectively;

[0016] FIG. 3C is a thresholded image based on the peak height of the ordered curve in FIG. 3B, yielding a squared order parameter (S.sup.2) value of 0.3150;

[0017] FIG. 4 is a graph showing squared order parameter (S.sup.2) as a function of time in the vicinity of a documented A3.4 class solar flare;

[0018] FIG. 5A is a graph showing the squared order parameter (S.sup.2) as a function of time in the vicinity of the 2017 X9.3 solar flare;

[0019] FIG. 5B is a graph showing the squared order parameter (S.sup.2) for C1.6, C1.0, and X8.12 solar flares;

[0020] FIG. 6A is an S.sup.2 mapping for the photosphere just before the 2017 X9.3 solar flare at 3:06 am UTC, wherein the image to the top left is a full disk image at 17.1 nm for the corresponding timeframe obtained by the Atmospheric Imaging Assembly (AIA) in NASA's SDO;

[0021] FIG. 6B is a mapping of the squared order parameter (S.sup.2) for the photosphere just before the 2017 X9.3 solar flare at 3:07 am UTC, wherein the image to the top left is a full disk image at 17.1 nm for the corresponding timeframe obtained by the Atmospheric Imaging Assembly (AIA) in NASA's SDO;

[0022] FIG. 6C is a mapping of the squared order parameter (S.sup.2) for the photosphere just before the 2017 X9.3 solar flare at 3:08 am UTC, wherein the image to the top left is a full disk image at 17.1 nm for the corresponding timeframe obtained by the Atmospheric Imaging Assembly (AIA) in NASA's SDO;

[0023] FIG. 7 is an Ising model plot of the solar polar field strength as a function of squared order parameter (S.sup.2), wherein the red, green, and yellow data represent contiguous time periods, each approximately one year long, wherein the inset is an enlarged view of the data; and

[0024] FIG. 8 is a flowchart showing a procedure for calculating the threshold value of an image.

DETAILED DESCRIPTION

[0025] It is to be understood that the processes described herein may assume various alternative step sequences, except where expressly specified to the contrary. It is also to be understood that the specific data and processes illustrated in the attached drawings, and described in the following specification, are simply exemplary embodiments of the inventive concepts defined in the appended claims. Hence, specific dimensions and other physical characteristics relating to the embodiments disclosed herein are not to be considered as limiting, unless the claims expressly state otherwise.

[0026] Known spin-based models of disorder (e.g., Ising/Potts models) and its corresponding quantitative measure have been utilized in various contexts. The present disclosure involves quantifying the order parameter of the photosphere of the Sun. The degree of order (e.g., the squared order parameter (S.sup.2) of the solar atmosphere) may be monitored to predict solar flares and other significant solar events.

[0027] In modeling disorder in the Sun's photosphere, hydrogen (H) and helium (He) atoms present in the photosphere represent the two “spins” in the Ising model. Other elements may be present in the photosphere, and they may be included in a model according to the present disclosure by adding each element as an additional spin. In general, this does not result in changes or errors in the extraction of the order parameter (S.sup.2) from the data.

[0028] As discussed in more detail below, the squared order parameter (S.sup.2) can be extracted from solar atmosphere images by thresholding image data into dark and bright regions (areas), followed by dividing the number of bright pixels by the total number of pixels within a region of interest. As also discussed in more detail below, a disorder analysis, according to the present disclosure, may be accomplished by software that is capable of performing near real-time, semi-automated image analysis from streaming or “live” sources. In one example (FIG. 1), this methodology of quantifying the order parameter was applied to solar photosphere images released from the Daniel K. Inouye Solar Telescope (DKIST). The images were captured at a wavelength of 705 nm. The resulting analysis detects (i.e., shows) oscillations in order parameter (S) as a function of time. In this example, the period of oscillations is approximately 200 seconds, which is somewhat similar to the oscillation period of approximately 300 seconds obtained from published observations of the power output of the Sun. In FIG. 1, an additional phase-delayed oscillatory component of smaller peak magnitude is also in evidence, and it may be a signature of two different-sized regions of cells oscillating out of phase within the image area. As discussed in more detail below, FIG. 1 shows the of the order parameter (S), which is calculated (determined) from the value of S.sup.2 that is directly obtained from image analysis.

[0029] With further reference to FIG. 2, the analysis has also been applied to a historical 10-year set of photosphere images from the Solar Dynamics Observatory operated by NASA. As shown in FIG. 2, the data contains multiple high-frequency components. The analysis also revealed an underlying multi-year period oscillation, which is highlighted in red in FIG. 2. Fourier analysis and/or curve fitting may be utilized to determine the period of oscillations. In the example of FIG. 2, a combination of Fourier analysis and curve fitting was used to determine that the multi-year oscillation has a period of 9.14 years, which is similar to the current estimate of 11 years for the period of the solar cycle.

[0030] The portion of FIG. 2 with higher spike density (2011-2015) corresponds to active years of the solar cycle where many solar flares have been observed. The data in FIG. 2 shows that the order parameter (S) and the directly related quantity of squared order parameter (S.sup.2) tend to spike approximately an hour before an observed solar flare. The data in FIG. 2 has a sample rate of one hour, however, so an additional analysis was performed on selected easily identifiable solar flares on a finer timescale of 2-minute intervals as shown in FIG. 4. FIG. 4 is a chart showing squared order parameter (S.sup.2) as a function of time in the vicinity of a document A3.4-class solar flare on May 17, 2018. The squared order parameter (S.sup.2) values of FIG. 4 were extracted from full disk images of the Sun taken at 17.1 nm wavelength by the Atmospheric Imaging Assembly (AIA) on NASA's Solar Dynamics Observatory (SDO). The data shows that approximately an hour before the solar flare there is a sharp jump in the order parameter, followed by a spike in squared order parameter (S.sup.2) at the time of the start of the solar flare. Thus, squared order parameter (S.sup.2) can be extracted from images of the solar atmosphere and monitored to provide an advance warning of approximately one hour prior to an impending solar flare of related solar event.

[0031] This process may be utilized to provide/generate an early warning to protect critical communication systems and power grids. Warnings are currently performed in the case of solar flares through “weather” modeling, which may not be sufficiently accurate. A disorder-based analysis, according to the present disclosure, may provide approximately one hour of warning, which may be sufficient to issue a protective directive to earth-orbiting satellites, other space probes in the solar system, and earth-based power grids, thereby reducing the disruption that otherwise occurs during a solar flare event. The disorder analysis may also be utilized to predict types of solar events or eruptions such as Coronal Mass Ejections (“CME”), which take longer to travel to earth, but may also be extremely disruptive.

[0032] FIG. 3A is a representative still frame from a photosphere movie acquired at 705 nm and released by DKIST. FIG. 3B is a pixel intensity histogram, showing the fitting of the overall intensity distribution and decomposition into two curves, corresponding to the disordered and ordered components of the structure of FIG. 3A, respectively. FIG. 3C is a thresholded image based on the peak height of the ordered curve in FIG. 3B, yielding a squared order parameter (S.sup.2) value of 0.3150. The thresholding process and pixel intensity histogram are discussed in more detail below in connection with FIG. 8.

[0033] In general, an image will include bright regions and dark regions, and the squared order parameter (S.sup.2) value of a sample is equal to the percentage of bright regions to the total image area. The bright and dark areas corresponding to the ordered and disordered regions, respectively, can be determined by thresholding the image near the average pixel intensity of the bright regions. This pixel value can be found by fitting the pixel intensity histogram with two curves representing the pixel intensity distribution from the disordered and ordered regions. The image threshold is generally chosen at the peak of the curve in FIG. 3B corresponding to the ordered regions, although it may be selected at an integer multiple of the standard deviation, σ, away from the peak depending on the relative location of the intersection of the two curves.

[0034] FIGS. 5A and 5B show detection of the onset of conditions that led to the eruption of an A3.4-class solar flare almost an hour in advance, using data extracted from the Atmospheric Imaging Assembly (AIA) on NASA's SDO. Far in advance of any intensity variation, the ordering of the solar atmosphere in the vicinity of the flare undergoes a significant change, indicating a physical reorganization that may enable the buildup of energy subsequently released in the flare. A primary interest in flare prediction is to protect systems on earth or in earth orbit so that it is worthwhile to consider what, if any, order parameter changes accompanied the X9.3 flare (FIG. 5A) in September 2017. Similar to the A-class flare in FIG. 4, there is an onset condition in the squared order parameter (S.sup.2) value of the full disk image of the Sun approximately an hour in advance of the X9.3 flare, as illustrated in FIG. 5A. This same approximate one-hour advance onset condition occurs as well for the C1.6, C1.0, and X8.2 class flares in FIG. 5B. Thus, the present disclosure provides a straightforward approach to solar flare prediction over a wide range of flare intensities based on order parameter extraction from telescope images.

[0035] In the examples discussed above, changes in squared order parameter (S.sup.2) sufficient to predict a solar flare occurred approximately one hour before the solar flares. However, it may be possible to predict solar flares more than an hour prior to the solar flare utilizing measurements and calculations of the order parameter according to other aspects of the present disclosure. Also, the information extracted from squared order parameter (S.sup.2) analysis of the full disk and smaller area images of the Sun can also provide information that could be used to enhance modeling of solar dynamics. For example, a “heatmap” type plot of squared order parameter (S.sup.2) can be generated from an image, including a full disk image, by dividing the thresholded image into a two-dimensional mesh of equal area triangles. The squared order parameter (S.sup.2) value of each region can then be calculated by taking the fraction of bright pixels to the total number of pixels within the segment to produce a map of squared order parameter (S.sup.2) across the image.

[0036] FIGS. 6A-6C show progression of squared order parameter (S.sup.2) mapping over a 5-minute interval 4 hours prior to the X9.3 flare identified in FIG. 5A. Such mappings aid visualization of the spatial redistribution of disorder within the system, which in FIGS. 6A-6C appears to be most active in the right near-equatorial side of the full disk, where a disorder-dominated ring surrounds a highly-ordered region which is the eventual eruption location of the X9.3 flare. The temporal variation of disorder and interactions between regions correspond to the dynamic behavior of kinetic and potential energy within the system. It can be shown that the kinetic energy associated with ordering is

[00001]$\frac{1}{2}{m\left[\frac{\mathrm{dS}}{dt}\right]}^2$

and the potential energy associated with ordering is ½kS.sup.2, where m and k are constants related to the characteristics of disorder within the system. Adding in a damping force that is proportional to order parameter (S), a simple damped harmonic oscillator Lagrangian of the form =(T−U)e.sup.αt can be constructed in terms of order parameter (S), specifically:

[00002]$\begin{array}{cc}ℒ=\frac{1}{2}\left[{m\left[\frac{dS}{dt}\right]}^2-k{S}^2\right]e^{\alpha t}& \left(1\right)\end{array}$

[0037] The oscillations seen in the order parameter (S) plot in FIG. 1 may be understood as undamped simple harmonic motion of order parameter (S) in the solar atmosphere. This result provides evidence that the state of ordering, represented through either the order parameter (S) or the squared order parameter (S.sup.2), represents a fundamental parameter of solar dynamics whose experimentally measured values can be used to potentially validate models and inform their development. Different layers of the solar atmosphere can be examined by, for example, acquiring images at specific wavelengths.

[0038] The ability to extract a value for the order parameter of the solar atmosphere may provide several advantages. Previous work has demonstrated that the order parameter can be related to specific system properties when those properties are dominated by pair interactions (see, e.g., Makin, R. A. et al., “Alloy-free band gap tuning across the visible spectrum,” Physical Review Letters 122, 256403, 2019; and Makin, R. A. et al., Revisiting semiconductor band gaps through structural motifs: An Ising model perspective,” Physical Review B 102, 115202, 2020). In the case of semiconductors, one such property is the band gap energy of the material, which exhibits a linear relationship with squared order parameter (S.sup.2). Using a cluster expansion up to pair-wise terms along with a spin-based representation of the system, such as the Ising-model, a system level property P can be expressed in terms of squared order parameter (S.sup.2) as:

P(x, S)=S.sup.2[P(0.5,1)−P(x, 0)]+P(x, 0)   (2)

[0039] In the context of the solar atmosphere, the mean polar field strength appears to exhibit a linear relationship with squared order parameter (S.sup.2) as predicted by Eq. 2. For example, FIG. 7 shows the mean polar field strength plotted as a function of squared order parameter (S.sup.2) over an entire decade of measurements (2010 to 2020) obtained from the Wilcox Solar Observatory and NASA's SDO. FIG. 7 is an Ising model plot of the solar polar field strength as a function of the squared order parameter (S.sup.2), determined from full disk images of the sun taken at a wavelength of 17.1 nm. The red, green, and yellow data represent contiguous time periods, each approximately one year long. The inset is a magnified view of the data to better highlight the linear trends. The data of the FIG. 7 was obtained from the Wilcox Solar Observatory and NASA's SDO.

[0040] For contiguous time periods of approximately one year in length, the measured mean polar field strength values all lie on the same S.sup.2 line—highlighted for three such time periods in FIG. 7. In terms of Eq. 2, this would imply that the composition of the system being imaged (x) is nominally the same for points that lie on the same line. These results suggest that the composition of the atmosphere varies with time. It can be seen that measurements of the order parameter of the solar atmosphere provide a variety of insights into the characteristics and behaviors of the Sun beyond simply the degree of ordering in the atmosphere.

[0041] With reference to FIG. 8, a process 10 for calculating the threshold value of a digital image 10 is shown. Process 10 may be implemented utilizing a computer that is programmed to one or more of the steps of process 10. Process 10 starts at 12, and proceeds to step 14, which includes converting the image to greyscale. A region of interest may be selected as shown in step 16. In general, the region of interest may be selected before or after converting an image to greyscale. The region of interest may be selected by an individual inspecting one or more solar atmosphere images, or the region of interest may be selected by a computer algorithm. For example, the algorithm may select a region of interest corresponding to a portion of the Sun's surface whereby the same portion of the Sun's surface is utilized in a series of images despite relative movement of the Sun.

[0042] The process 10 further includes calculating a pixel intensity histogram of the selected region (see, e.g., FIG. 3B). At step 20, a computer-implemented algorithm (e.g. a stochastic funnel algorithm) is used to calculate initial fitting parameters for two skewed Gaussian curves to the pixel intensity histogram data. In the example of FIG. 3B, the skewed Gaussian curves are shown as the disordered and ordered curves, and the over-all fit is also shown in FIG. 3B. Referring again to FIG. 8, at step 22, the least squares method is used to calculate the Gaussian curve parameters using the initial fitting parameters from step 20.

[0043] At step 24, a root-finding algorithm (e.g., Newton's method) is used to find the) standard deviations that the intersection is away from the curve where the highest center point is calculated. The threshold value is set to the value of the highest center point value minus the floor of that number of standard deviations.

[0044] At step 28, a binary threshold is performed on the region of interest in the image using the threshold calculated in step 26. An example of the results of thresholding are shown in FIG. 3C. At step 30, the squared order parameter (S.sup.2) value of the region of interest is calculated by counting the bright pixels in the binary image and dividing this number by the total number of pixels contained with the region of interest. The method 10 then ends as shown at 32.

[0045] Although only a few embodiments of the present innovations have been described in detail in this disclosure, those skilled in the art who review this disclosure will readily appreciate that many modifications are possible (e.g., variations in the process used to extract order parameter (S) and squared order parameter (S.sup.2), the order of the steps, values of parameters, use of colors, orientations, etc.) without materially departing from the novel teachings and advantages of the subject matter recited. Also, as used herein, the term “order parameter” may generally refer to the order parameter (S), order parameter (S.sup.2), and/or other variables or values that incorporate or relate to the order parameter as described herein.

[0046] It will be understood that any described processes or steps within described processes may be combined with other disclosed processes or steps within the scope of the present device. The exemplary processes disclosed herein are for illustrative purposes and are not to be construed as limiting.

[0047] The above description is considered that of the illustrated embodiments only. Modifications of the process will occur to those skilled in the art and to those who use the process. Therefore, it is understood that the embodiments shown in the drawings and described above are merely for illustrative purposes and not intended to limit the scope of the invention, which is defined by the following claims as interpreted according to the principles of patent law, including the Doctrine of Equivalents.