Solar Wind Forecasting with Coronal Holes 1 2 2 S. Robbins , C. J. Henney and J. W. Harvey 1APS Department, University of Colorado, Boulder, CO 80309, USA 2National Solar Observatory, Tucson, Arizona, 85719, USA Abstract. Anempiricalmodelforforecastingsolarwindspeedrelatedgeomagnetic events is presented here. The model is based on the estimated location and size of solarcoronalholes.Thismethoddiffersfrommodelsthatarebasedonphotospheric 7 magnetograms (e.g., Wang-Sheeley model) to estimate the open field line configu- 0 ration. Rather than requiring the use of a full magnetic synoptic map, the method 0 presented here can be used to forecast solar wind velocities and magnetic polarity 2 from a single coronal hole image, along with a single magnetic full-disk image. The n coronal hole parameters used in this study are estimated with Kitt Peak Vacuum a Telescope HeI 1083 nm spectrograms and photospheric magnetograms. Solar wind J and coronal hole data for the period between May 1992 and September 2003 are 9 investigated.Thenewmodelisfoundtobeaccuratetowithin10%ofobservedsolar wind measurements for its best one-month periods, and it has a linear correlation 1 coefficient of ∼0.38 for the full 11 years studied. Using a single estimated coronal v hole map, the model can forecast the Earth directed solar wind velocity up to 8.5 5 7 days in advance. In addition, this method can be used with any source of coronal 2 hole area and location data. 1 0 7 0 1. Introduction / h p - Prediction of space weather near the Earth is a major goal of solar o research. An important aspect of attaining this goal is to accurately r st describe the solar drivers of space weather. The drivers are the solar a wind and the various phenomena that shape and modulate that wind. : v Among the earliest findings from space observations of the solar wind i X was that it consisted of recurrent low-speed, dense streams and high- speed tenuous streams, and that the latter were strongly associated r a with increased geomagnetic activity (Synder, Neugebauer, and Rao, 1963). Many suggestions were made that the high-speed solar wind streams might be associated with regions on the Sun having magnetic fields open to interplanetary space. When coronal holes were found to be regions likely to have open magnetic fields (Altschuler, Trotter, and Orrall, 1972), it was not long until Krieger, Timothy, and Roelof (1973), and several other investiga- tors, demonstrated a link between open-field coronal holes, high-speed solar wind streams, and enhanced geomagnetic activity. In an effort to strengthen this linkage, Sheeley, Harvey, and Feldman (1976) con- structed time-stacked diagrams of coronal holes, solar wind speed and (cid:13)c 2008 Kluwer Academic Publishers. Printed in the Netherlands. Paper_v10.tex; 4/02/2008; 0:35; p.1 2 Figure 1. AsampleNSO/KPVTcomputer-assistedhand-drawncoronalholeimage (left) and a EIT 19.5 nm Fe XII emission line image (left) for July 14, 2003 at approximately 17 UT. Note that the coronal hole regions appear dark in the EIT image. geomagnetic activity in one-rotation-long rows. The diagrams covering the years 1973-1975 strongly supported the linkage and the authors suggestedthatobservationsofcoronalholescouldbeusedtopredictthe arrival of high-speed streams and their associated geomagnetic activity a week in advance. Coronal holes are best seen against the solar disk as low-intensity regions in space observations of material at coronal temperatures. This can also be done from the ground using radio observations. Harvey et al. (1975) found that coronal holes could be seen faintly in ground- based images made with helium lines such as 587.6 and 1083.0 nm because the strength of these lines is partly controlled by the intensity of overlying coronal radiation (see, e.g., Andretta and Jones, 1997). A program of regular 1083 nm observations has been conducted by the National Solar Observatory (NSO) Kitt Peak Vacuum Telescope (KPVT) starting in 1974. Among the derived products are estimates of the locations and magnetic polarity of coronal holes. An example coronalholeestimateimagederivedfromaKPVTobservationisshown in Figure 1, along with a 19.5 nm Fe XII emission line image measured by the Extreme Ultraviolet Imaging Telescope (EIT) for comparison. Harvey and Recely (2002) describe how coronal holes are identified using KPVT He I 1083 nm observations. Predictions of solar wind speed at Earth are regularly made by several groups based on solar potential field extrapolations (e.g., http://solar.sec.noaa.gov/ws/, http://www.lmsal.com/forecast/, http://bdm.iszf.irk.ru/Vel.html, and http://gse.gi.alaska.edu) and in- Paper_v10.tex; 4/02/2008; 0:35; p.2 3 terplanetary scintillation (Hewish et al., 1964) observations (e.g., http://cassfos02.ucsd.edu/solar/forecast/index.html, and http://stesun5.stelab.nagoya-u.ac.jp/forecast/). Theformersetoffore- casts are based on extrapolation of photospheric longitudinal magnetic field measurements using a potential field assumption to locate open field lines. Using one of these models, a modified Wang and Sheeley (1990, 1992) flux-transport model, Arge and Pizzo (2000) studied a three-year period centered about the May 1996 solar minimum. They comparedpredictedsolarwindspeedandmagneticpolaritywithobser- vationsnearEarth.Theirthree-yearsampleperiodhadanoverallcorre- lation of ∼ 0.4 with observed solar wind velocities and an average frac- tionaldeviation,ξ,of0.15,whereξ = h(prediction−observed)/observedi. When excluding a 6-month period with large data gaps, they correctly forecast the solar wind to within 10 − 15%. Interplanetary magnetic field (IMF) polarity was correctly forecast ∼ 75% of the time. In this paper, we address the suggestion of Sheeley, Harvey, and Feldman (1976) that observations of coronal hole regions can be used to predict the solar wind speed at Earth as much as a week in advance. Inaddition,themodelpresentedhereisbasedonobservationsthatfind moderate and high-speed solar wind streams are associated with small andlargenear-equatorialcoronalholes,respectively(Nolteetal.,1976). Here we correlate the coronal hole percent area coverage of sectoral regions of the observed solar surface with solar wind measurements to derive a simple empirical model (discussed in Sections 2 through 5). As a measure of the merit of this model for solar wind forecasting, we comparepredictionswithobservationsandcontrastthistechniquewith the ones based on magnetic field extrapolations (Sections 5 and 6). 2. Model Input Data The coronal hole data used here are based on KPVT observations from May 28, 1992 through September 25, 2003 (i.e. the last half of cycle 22 and the first half of cycle 23). The coronal hole locations and area estimates are from computer-assisted, hand-drawn maps (see Figure1)basedupontheKPVTHeI1083nmimagesandphotospheric magnetograms (Harvey and Recely, 2002). For this investigation, the estimated coronal hole boundaries were mapped into sine-latitude and longitude to create heliographic images. The coronal hole region image pixels are set to a value of 1, whereas the background is defined as 0. For the time period analyzed here, the KPVT coronal hole maps have a 69% daily coverage. Paper_v10.tex; 4/02/2008; 0:35; p.3 4 The solar wind speed data utilized here was obtained from the OM- NIWebwebsite(http://nssdc.gsfc.nasa.gov/omniweb/) providedbythe National Space Science Data Center. Daily averages of the solar wind speed time series were created with the approximate cadence of the KPVT-based coronal hole maps. For the time period analyzed here, the solar wind speed time series has a 92% daily coverage. Data gaps in the time series are interpolated using a cubic spline. 3. Solar Wind Correlation Analysis Forcomparisonwiththesolarwindspeedtimeseries,eachheliographic coronal hole image was divided into 23 swaths (i.e. sectoral regions) 14 degree-wide in longitude overlapped by 7 degrees. The approximately 1-day-wide longitudinal window was selected to correspond with the temporal cadence of the KPVT observations. These sectoral samples are then summed, where each pixel corresponding to a coronal hole is equal to 1, to yield a percent coverage of that area by coronal holes. For each coronal hole image there may be no or only a few coronal hole regions observed for that time. For example, swath sectors with no coronal hole regions would yield a hole coverage of zero percent. This is repeated for each coronal hole image available in the 11-year period to form a coronal hole time series for each of the 23 sectoral samples. Each sectoral time series is then interpolated into the time frame of the solar wind velocity data. Thecorrelation andtime lag between the time series were estimated with weighted cross-correlations (e.g. Bevington and Robinson, 2003). The weighted cross-correlation simplifies the analysis by allowing the use of the continuous time series. The gap-filled data are given small weights to minimize their contribution while the measured or derived values are given equal and relatively large weight values. In addition, following Arge et al. (2004), periods of CME events were estimated using the plasma β value (obtained from the OMNIWeb data set) when β ≤ 0.1. For periods estimated to correspond to a coronal mass ejection (CME) event, solar wind speed values were given negligible weight values. 3.1. Longitudinal Forecasting Windows Twenty-threelongitudeswathswereexamined,rangingfrom77degrees east to 77 degrees west of the central meridian. Each of these was cross-correlated with the solar wind speed time series in the manner described above. Figure 2 exhibits the maximum cross-correlation co- efficient within a lag window of ±0.5 days for each longitude swath Paper_v10.tex; 4/02/2008; 0:35; p.4 5 0.5 nt e effici 0.4 o C n 0.3 o ati el orr 0.2 C m u +/-90 deg Latitude Range xim 0.1 +/-70 deg Latitude Range a +/-60 deg Latitude Range M 0.0 -70 -56 -42 -28 -14 0 14 28 42 56 70 Central Meridian Distance Figure2. LongitudeCross-Correlations:14degree-wideswathsofthesolardiskwere tiled across the solar surface with 7 degrees overlap between successive swaths for 3latituderanges.Themaximumcross-correlationswithinalagwindowof±0.5days (corresponding to a day offset described by Equation (1)) of the resulting data set with the solar wind velocity data are shown here. Negative and positive longitude values correspond to eastward (E) and westward (W) respectively. Note the two relatively good forecast windows centered around 63 degrees east and the central meridian. The correlations continue to improve with model time series using de- creasinglatituderanges,peakingwith thelatituderangeof±60degrees.Thetrend reverses for latitude ranges narrower than ±60 degrees. investigated with the latitude ranges ±90 degrees, ±70 degrees, and ±60 degrees with equal latitude weights. Note that besides the central meridian peak, the correlation also peaks towards the east and west limbs. This is a result of better coronal hole detection from the KPVT He I 1083 nm spectroheliograms towards the image limb. In addition to the three latitude bands shown in Figure 2, the correlation values were also determined for ±50 and ±40 deg cases. These two latitude bands are similar to the ±60 degree case but have lower correlation values eastward of -56 degrees. For the ±60 deg case, two preferred forecasting windows centered at 63 degrees east and on the central meridian are clearly visible. The estimated swath time series lags corresponding to maximum cross-correlation coefficients are found to be linear, and can be ex- pressed as: d = (3.69±0.02)−(0.07386±0.0005)θ, (1) where d is the time delay in days and θ is the center of the longitude swath in degrees as measured from the central meridian (east is neg- Paper_v10.tex; 4/02/2008; 0:35; p.5 6 0.5 a m xi 0.4 a M n 0.3 o ati el orr 0.2 C e of 0.1 d plitu 0.0 CMD C Sewntaetrhe dP oats i6ti3o°n E m A Centered at 0° -0.1 N N+Eq Eq S+Eq S Latitude Bin Figure3. Thetwolongitudeforecasts’ cross-correlationswiththesolardiskdivided intofivelatitude regions described in §3.2. Seetext for discussion. ative, west is positive) and the uncertainty values are 1-σ. Note that Equation (1) is valid within the central meridian distance range of -80 to 80 degrees, but it only has physical meaning (days forecast) for the centralmeridiandistancerangeof-80to40degrees.FromEquation(1), the two preferred forecasting windows centered at 63 degrees east and on the central meridian (see Figure 2), the time lags of these forecasts correspond to 8.3 days and 3.7 days respectively. In other words, the delay between detected solar wind variations at Earth and a coronal hole region at the central meridian (θ = 0 degrees) is approximately 3.7 days. A coronal hole region observed at 63 degrees east central meridian distance (θ = −63 degrees) is associated to solar wind speed variations at Earth approximately 8.3 days later. Theobserved delay is as expected, and is best explained as the result of a corotating stream of plasma moving nearly radially outward from the sun (e.g. Gosling, 1996). 3.2. Latitude Weighting Analysis In addition to the longitudinal correlation analysis above, the heli- ographic coronal hole images were divided into three latitude bins: Northern (90 degrees north to 30 degrees north), Equatorial (30 de- grees north to 30 degrees south), and Southern (30 degrees north to 90 degrees south) regions. Combinations of Northern with Equatorial and Southern with Equatorial were used for a total of five bins. The cross-correlation coefficient values of each latitude bin were calculated for the two longitudinal windows discussed in §3.1. Paper_v10.tex; 4/02/2008; 0:35; p.6 7 The amplitude of the correlation coefficient maxima shown in Fig- ure 3 reveals a bias towards the southern hemisphere. This bias in the cross-correlation between the hemispheres is most likely due to the significantly greater number of coronal holes detected in the southern hemispherefor theperiodinvestigated (e.g. Henney andHarvey, 2005). Assumingthatthehemisphericasymmetryisaresultofthelimiteddis- tribution sample of coronal holes, and following the analysis discussed in the previous section, the latitude range of ±60 degrees was used for the forecasting analysis done below. The maximum weighted cross- correlations for the East window was 0.376 and the central meridian had a correlation of 0.378 (see Figure 2). Amodelsolar windtimeseries, Vmod, is created by firstdetermining thearea percentage of the 14-degree widesectors that is acoronal hole, Is.Thesecoronalholepercentagevaluesforeach longitudinalswathare rescaled to agree with the observed solar wind speed, Vobs, using the linear scaling coefficients α and γ, where Vmod = α+γIs. The linear scaling coefficients were determined by minimizing the mean of the absoluteaverage fractionaldeviation,ξ,ofthemodelfromtheobserved values, where ξ = h(Vobs −Vmod)/Vobsi. (2) Using the above criteria and a latitude range of ±60 degrees, the av- erage linear scaling values for each longitude window are found to be: α = 330 km/s and γ = 930 km/s. With these scaling factors and the best weighting as discussed above, the resulting ξ for the entire 11- year data set for the central meridian swath is ∼16% with a standard deviation of ±20%. 4. Magnetic Activity Cycle Dependence Sheeley and Harvey (1981) reported a dependencebetween the coronal hole and solar wind correlation and the sunspot cycle. This was quan- titatively explored for the time period spanned by the coronal hole data set - the last half of cycle 22 and the first half of cycle 23. The full time series was divided into six time series of approximately 690 days, illustrated in the upper graph in Figure 4 with sunspot number timeseries(sunspotcountdatausedherewasobtainedfromtheNGDC website(http://www.ngdc.noaa.gov/stp/), maintainedbytheNational Geophysical Data Center). The date range and percent coverage of the six periods are outlined in Table I. Figure4highlightsastrongdependenceofthetimeseriescorrelation values on the phase of the sunspot cycle, similar to the qualitative Paper_v10.tex; 4/02/2008; 0:35; p.7 8 250 200 er mb 150 Bin 1 Bin 2 Bin 3 u N ot p 100 s n u S 50 Bin 4 Bin 5 Bin 6 0 1988 1990 1992 1994 1996 1998 2000 2002 2004 Time (Year) 0.5 a m xi 0.4 a M n o ati 0.3 el orr C of 0.2 e d u plit 0.1 Centered at 63° E m Centered at 0° A 0.0 1 2 3 4 5 6 Time Series Bin Figure4. Monthly(soliddots)andyearly(line)sunspotnumberaveragesareshown in the upper panel along with six time bins into which the 11-year data set was divided.Thesixbinsareindicatedbytheroundedrectangles andeach representsa durationofapproximately 690days.Daterangesforeach period arelisted in Table I.Showninthelowerpanelarethecross-correlation valuesofthecoronalholedata withthesolarwindforthe6intervalsillustratedintheupperfigure.Thecorrelation is found to be best during the declining phase of the sunspot cycle and worst just after solar minimum. results of Sheeley and Harvey (1981). The correlation between the coronal hole and solar wind time series is best during the declining phase of the cycle, and it is worst during solar minimum and the beginning ascending phase of the cycle. Some of the lack of correlation can be attributed to the observed latitudinal distribution of coronal holes with respect to the solar cycle. During the minimum phase of the Paper_v10.tex; 4/02/2008; 0:35; p.8 9 Table I. Time series subdivisions used in the cross-correlation analysis relative to thesunspot numbertime series shown in Figure 4. Period Bin DateRange Interval(days) Completeness 1 May 28, 1992 - Apr17, 1994 690 66.1% 2 Apr18, 1994 - Mar 7, 1996 690 68.4% 3 Mar 8, 1996 - Jan 26, 1998 690 72.2% 4 Jan 27, 1998 - Dec17, 1999 690 73.8% 5 Dec18, 1999 - Nov 6, 2001 690 64.6% 6 Nov 7, 2001 - Sep 25, 2003 688 66.3% 1 - 6 May 27, 1992 - Sep 25, 2003 4138 68.6% solar cycle, fewer low-latitude coronal holes are observed which means few fast streams are observed by spacecraft in the ecliptic plane (e.g. Woch et al., 1997; Kojima et al., 2001; Kojima et al., 2004). However, the source of the poor correlation during solar minimum is also likely a result of a noted difficulty in the determination of coronal hole regions using He I 1083 nm spectroheliograms during periods of low magnetic activity. Both longitude windows vary over approximately the same range throughout the cycle, and even though during the worst forecasting period the correlation is far below the full time series’ average, it is still well above statistical significance. Correlation significance, the probability that two uncorrelated random sets of variables withagiven numberof observations would give similar correlation values, was calculated following Appendix C in Bevington and Robinson (2003). The correlation is considered highly significant and nominally significant if the probability of chance occurrence is less than 1% and 5%, respectively (Taylor, 1997). The time delays between coronal hole observation and the effects seeninthesolarwindvelocity areslightly longerduringsolarminimum and slightly shorter during solar maximum than the averages quoted in §3.1. The lags range over 8.38±0.34 and 3.82±0.58 days. 5. IMF Polarity Forecast Besides area and position information, the Harvey and Recely (2002) coronalholeboundaryimagesincludemagneticpolarity.So,inaddition to forecasting the solar wind velocity, the KPVT estimated coronal hole maps can be used to predict the interplanetary magnetic field (IMF)polarityatEarth.Thelongitudinalswathcenteredatthecentral Paper_v10.tex; 4/02/2008; 0:35; p.9 10 meridian, the 3.7-day forecast, was used in the following analysis. The heliographic coronal hole maps are scaled so that each pixel with a positive polarity hole has a value of +1 and each pixel with a negative polarity hole has a value of −1, whereas non-coronal hole regions are set to 0. The average value of all the pixels in the 14 degree wide swath wastaken,witharangebetween±60degreesinlatitude.Theaveraging of the coronal hole polarity, albeit simple, is treated here a as baseline for the polarity forecasting when using only coronal hole regions. Excludingonly thedays that didnothave bothvelocity andcoronal hole data, leaving 63.0% of the comparison period, the IMF polarity was correctly forecast 57.9% of thetime. When excludinganadditional 17.9% of the days where the average model magnetic polarity was within 0.1% of 0, the IMF polarity is correctly forecast for 63.5% of the time. This is approximately 10% lower than the values reported by Arge and Pizzo (2000); however, we expect improvement with future models.Thoughthecurrentmodelhas inherentinaccuracies as aresult of the oversimplification of the magnetic field structureassociated with coronal holes and the resulting solar wind, the use of higher signal-to- noise ratio magnetograms is expected to improve the polarity forecast. Inaddition,thismodelispartlybasedontheassumptionthatthesolar wind velocity is related linearly with the size of the coronal hole. How- ever, it has been shown that there is a critical scale size for which the wind velocity is independent of coronal hole size (Kojima et al., 2004). In future models, we plan to include coronal hole size and additional topology-related parameters to potentially improve forecasts. 6. Forecast Comparison and Discussion Figure5illustrates two sampleforecastperiodsthatcover twoCarring- ton rotations (CR) each: CR 1862 and 1863 (top), and CR 1955 and 1954 (bottom). To objectively find periods of good and bad forecasts, both weighted cross-correlations for 90-day periods as well as absolute average fractional deviations, see Equation (2), were performed. For each time series, only valid data (non-gap-filled data values) are shown in the figures. Ranges quoted in this section are from the two different forecast windows (discussed in §3) centered at 63 degrees east and 0 degrees relative to the central meridian. The top forecast comparison shown in Figure 5 is a sample period when the modeltime series is well-correlated with the solar wind speed data. For the two forecast windows, the weighted correlation coeffi- cients range from 0.626 to 0.698, and the absolute average fractional deviation, ξ, is found to be 0.098. In comparison, the best one-month Paper_v10.tex; 4/02/2008; 0:35; p.10