© ESO, 2017

1. Introduction

Our ability to probe the solar chromosphere has seen significant advances in the recent past. The increased spatial, temporal, and spectral resolution of both ground and space based observatories has provided unique insights into the highly dynamic behaviour of fine-scale features. The observation of energy transfer by magnetohydrodynamic (MHD) waves is one area that has benefited from these improvements, with the identification of ubiquitous transverse waves – at least in chromospheric features associated with enhanced photospheric magnetic field concentrations that emanate from the magnetic network and plage, that is, spicules (e.g. De Pontieu et al. 2007) and fibrils (e.g. Pietarila et al. 2011; Morton et al. 2012a).

Estimates for the energy flux imply that the waves carry enough energy for plasma heating. However, current techniques for estimating the energy flux are generally crude, meaning that the actual energy content of the waves and their role in atmospheric plasma heating is subject to uncertainty (see e.g. Goossens et al. 2013; Van Doorsselaere et al. 2014). There is further uncertainty surrounding the mechanism(s) by which these waves can actually deposit their energy in the local plasma. By their very nature they are incompressible and difficult to dissipate. Currently, the most favoured mechanism for the damping of the kink mode is resonant absorption (Terradas et al. 2010; Verth et al. 2010), although this process is just a transfer of energy from the kink modes to the quasi-torsional m = 1 Alfvén modes, and requires a further mechanism to dissipate their energy, for example, phase mixing of the m = 1 Alfvén modes (Soler & Terradas 2015). An alternative mechanism may involve the generation of instabilities at the boundaries of the flux tube by the kink modes (Terradas et al. 2008; Antolin et al. 2015). In order to develop our understanding, detailed observations and analysis of MHD waves in the solar atmosphere are required, and one of the key objectives will be to determine the role the chromosphere plays in regulating the energy flow in the atmosphere.

In general, recent studies of chromospheric phenomena have focused on events that occur around network concentrations of magnetic flux and plage regions, with a wealth of individual features identified (e.g. Rapid blue-shifted events – RBEs, Type I and II spicules; Rouppe van der Voort et al. 2009; Pereira et al. 2012). This focus is primarily due to the magnetic connection of these features to the corona and the hypothesis that they contribute to the mass and heat flux required to sustain the corona (e.g. De Pontieu et al. 2011; Pereira et al. 2014; Henriques et al. 2016). In contrast, there has been much less attention given to the study of fibrils, which reach out from the network across the internetwork (IN), forming dense canopies that occupy a large volume of the visible chromosphere¹ (Rutten 2006). In addition to the network fibrils there exist shorter, absorbing Hα features that appear in the IN, which have gone somewhat unmentioned in literature (a discussion of these features is the subject of this paper).

Observations reveal that the fibrils display time-dependent behaviour, appearing, disappearing, and reappearing over tens of minutes, likely indicating a departure from a hydrostatic plasma. Such a thermodynamic cycle must then require an additional and time-dependent energy input into the system. It has been suggested that the fibrils are the signatures of heating events or at least the markers that heating events may have occurred at an earlier, recent instance in time (Rutten & Rouppe van der Voort 2017). However, there is no consensus about the processes that lead to the formation of fibrils. Therefore, it seems prudent to engage in detailed studies of the fibrils using the improved capabilities of modern instruments to gain insight into their nature, and, as such, shedding light on the heating of the chromosphere.

In anticipation of later discussion here, we summarise some salient features of IN magnetism here. On average the IN magnetic flux makes up 14% of the total Quiet Sun flux and appears to migrate towards the network boundaries (potentially due to supergranule flows), providing a substantial contribution to sustaining the network flux through mergers (e.g. Wang & Zirin 1988; Wang et al. 1995; Gošić et al. 2014). Typical estimates of photospheric IN magnetic field strengths are on the order of ~200 G but there is evidence for the existence of kGauss features (Domínguez Cerdeña et al. 2003; Orozco Suárez & Bellot Rubio 2012). The IN elements have lifetimes on the order of minutes but some survive longer in order to contribute to the network (Zhou et al. 2010). Extrapolations that include the IN field demonstrate that short, closed magnetic fields should be ubiquitous (Schrijver & Title 2003; Wiegelmann et al. 2010), with only a small fraction reaching the corona. Wiegelmann et al. (2010) suggest photospheric field strengths above 300 G contribute to around 90% of the chromospheric magnetic energy. Analogous low-lying magnetic fields have potentially been identified in recent high spectral resolution Ca II H data of active regions (Jafarzadeh et al. 2017a).

The magnetic fields in the photosphere are also observed to be highly dynamic. Study of the motions of magnetic bright points in G-band images (Berger & Title 1996; Nisenson et al. 2003; Chitta et al. 2012) and flux concentrations in magnetograms (Giannattasio et al. 2014b,a) has revealed that these small-scale magnetic elements migrate through the photosphere, whether due to some local flow or larger supergranule flow pattern. Superimposed on the longer timescale migration are random motions of shorter timescales due to granular buffeting (e.g. de Wijn et al. 2005). It is assumed that such motions can excite MHD waves, in particular transverse waves (e.g. Choudhuri et al. 1993). Stangalini et al. (2015) provided evidence that these motions excite kink waves in chromospheric network bright points in Ca II H (3969 Å). Morton et al. (2013, 2014) also demonstrate a correspondence between photospheric bright point motions and the transverse waves of long network fibrils in Hα (cf. Hillier et al. 2013, for prominences). There is no evident reason to believe that photospheric velocity fields cannot also excite waves along IN magnetic fields.

In this study, we begin to examine the nature and dynamics of chromospheric features located in the IN, focusing in particular on the apparent periodic transverse displacement of the IN fibrils (see Sect. 3 for definition and description). The periodic displacements of the IN fibrils axis are interpreted in terms of the kink wave. The data used here also has a very high temporal cadence (~1 s), which provides an ideal opportunity to measure the range of propagation speeds of kink waves in the chromosphere, something that has proved difficult to do in previous studies due to the large Alfvén speed and short-scales of the features (Jess et al. 2015). Furthermore, we exploit measured wave properties to investigate the nature of the IN fibrils.

Fig. 1 Hα line core and blue wing images are shown in the left and centre panels respectively. The two major magnetic flux concentrations that form part of the network and associated chromospheric rosettes are evident. Right hand panel: photospheric magnetogram. The levels have been clipped to B ~ ± 200 G to reveal some of the weaker fields in the IN. Contours highlight absolute magnetic field strengths greater than the cut-off, red showing positive flux and blue negative.

2. Observations and data reduction

The observations presented in the following are focused on the boundary of a coronal hole at the disk centre (Heliocentric Cartesian – X0′′, Y31.5′′). The dataset used here was taken by the Swedish Solar Telescope (SST – Scharmer et al. 2003) using the CRisp Imaging SpectroPolarimeter (CRISP – Scharmer 2006, Scharmer et al. 2008) between 09:06 and 09:35 UT on 2013 May 3 at the Roque de los Muchachos Observatory, La Palma in the Canary Islands. The main sequence is a spectral scan of Hα (6563 Å) including the following wavelength positions from line centre: (−0.91, −0.54,−0.36, 0.0, 0.36, 0.54, 0.91) Å, corresponding to a range of ±41 km s^-1 in velocity. The data was reconstructed with an extended multi-object multi-frame blind deconvolution (MOMFBD) process (van Noort et al. 2005; Henriques 2012; de la Cruz Rodríguez et al. 2015) and de-rotated, aligned, and de-stretched (Shine et al. 1994). Post-reconstruction, the cadence of the full spectral scan is approximately 1.3 s with 1150 frames, the spatial sampling is ~0.̋06 per pixel, and the resolution is 0.̋16. CRISP also undertook photospheric Fe I 6301 and 6302 Å spectral scans every five minutes apart over the same field of view (FOV) to obtain full Stokes profiles. The Stokes V component was used to construct line-of-sight (LOS) magnetograms (see Kuridze et al. 2015, for additional details). The scans of the iron lines then leave 30 s gaps in the Hα time-series. Other works that used this dataset include Samanta et al. (2016), who studied the impact of the lifetime of chromospheric structures on Fourier power spectra, and Henriques et al. (2016), who found a connection between RBEs and Rapid red-shifted events (RREs) and emission in SDO pass bands sensitive to transition region and coronal temperatures.

The highly dynamic nature of the chromosphere leads to Doppler shifts of line profiles and can cause problems when observing phenomena at a fixed wavelength position, mixing fluctuations in plasma conditions and velocities. To overcome this we also determine the intensity at the location of maximum absorption. Each spectral line profile is fit with a polynomial to determine the intensity at the minimum location and the wavelength position of the profile minimum².

3. Fibrils in the internetwork

The observed FOV is displayed in Fig. 1 and reveals a top-down view of the lower solar atmosphere. The Fe magnetograms and wing images from Hα show the state of the photosphere, revealing that there are two dominant regions of negative polarity, kilogauss magnetic fields that are likely part of the network. They also reveal the existence of widespread, small-scale patches of magnetic flux away from the network fields, among which exist a small number of compact magnetic concentrations with | B | > 200 G.

The network magnetic fields are visible as clusters of bright points in the photosphere and Hα line core images suggest that the network has a dominant influence on the visual appearance of the chromosphere. In this FOV, the chromospheric fine-scale structure predominantly originates from these network regions. Spicules are observed directly over the network, with an apparent orientation that is near perpendicular to the surface, and long fibrils are also evident around these kG fields, with clear structuring that forms two distinct rosettes. The fibrils extend out near-radially from the network to around 10−15 Mm into the IN. Although, this behaviour differs in regions where there is a strong positive polarity field (e.g. [18, 4] and [26, 15] Mm) and shorter fibrils exist.

Additionally, a number of relatively long-lived (>200 s), extended curvilinear features are found in the internetwork that show absorption in Hα line-core images. We will refer to these features as IN fibrils and examples are shown in Fig. 2. The IN fibrils are found to disappear and re-appear in approximately the same location over the length of the dataset and follow the same course. This behaviour would imply the persistence of an underlying magnetic field and its visibility may be just due to a variation in the Hα opacity. The existence of long-lived chromospheric magnetic fields with varying opacity, occupying similar spatial locations, is not at odds with the proper motion speed of photospheric IN magnetic elements (average speeds of ~0.2 km s^-1; Wang & Zirin 1988). Their visual appearance and behaviour are similar to the fibrils that protrude from the network, which would suggest that these features are dense plasma that outlines the IN magnetic fields. There is also some similarity in behaviour with slender Ca II H fibrils observed in an active region, reported recently from the balloon-borne Sunrise data (Gafeira et al. 2017). The elongated nature of the IN features in this FOV would imply, if they do outline the magnetic field, that the magnetic field is relatively low-lying in the atmosphere. This is opposed to being near-vertical features like spicules, which protrude significantly into the corona and predominantly have limited inclinations from the vertical (Tsiropoula et al. 2012; Pereira et al. 2012). This situation represents the view that the Hα chromosphere is essentially a corrugated fibrilar canopy, with fibrils generated in the 3D simulations found to exist up to heights considered to be coronal (Leenaarts et al. 2012). Furthermore, Hα limb observations (E. Scullion, priv. comm.) also appear to show fibrils rising above the “bulk chromosphere” (i.e. the non-fibrilar chromosphere described in Judge & Carlsson 2010). The IN fibrils may be a visible signature of low-lying loops that are suggested to dominate the IN by the magnetic field extrapolations of Wiegelmann et al. (2010). Alternatively, the IN fibrils’ curvilinear nature could be in the “plane” of the chromosphere, that is, the visible sections of the fibrils are horizontal and snake through the chromosphere. However, it is near impossible to determine this from a visual inspection of the current data set. In Sect. 4.4, we employ magneto-seismology to gain some insight into this aspect of the fibrils.

Fig. 2 Examples of internetwork fibrils. Left hand column: intensity images determined from the profile minimum, revealing the chromospheric absorption features. Middle panels: corresponding wing images; and last column: photospheric magnetic field strength clipped between ± 50 G.

Fig. 3 Panel a: one of IN fibrils in the data set, which can be seen at (8, 20) Mm in the line core image in Fig. 1. Panel b: a parcel of plasma is seen at the upper end of the fibril, significantly broadening the apparent width of the feature. The dark clump can be seen to move along the apparent axis of the feature in panels c and d (from upper right to lower left).

From examining the data, it proves difficult to determine the photospheric footpoints of the IN fibrils. Both in Hα wing images, which show a photospheric scene, and magnetograms, there are no obvious magnetic field concentrations that can be directly associated with the visible endpoints of the IN fibril, if we compare, for example, absorption features seen in the profile minimum to line wing and magnetogram images in Fig. 2. Taking the example in the second row, the left hand-end of the IN fibril lies close to a negative small patch of negative polarity (2, 2.5) Mm that is potentially its footpoint. However, there is no clear positive polarity patch at the other visible end of the feature, ~(6, 8) Mm. It has been noted previously by Reardon et al. (2011) that network fibrils typically only have one clearly identifiable footpoint, that being the network patch it originates in, while the other footpoints are in the IN with their exact location significantly harder to establish.

A number of on-disk chromospheric features located around the network are observed to have a significant blue wing enhancement (Rouppe van der Voort et al. 2009) that is thought to be related to an initial heating event (Rutten & Rouppe van der Voort 2017). A cursory investigation found no obvious heated precursors in the IN fibrils that we examine.

However, we observe an example of absorbing material flowing along the length of one of these IN features after formation (Fig. 3). The observation sequence begins at ~09:06 and this particular IN fibril is present then, hence has been in existence for an unknown length of time. At ~09:07 a broad absorption feature is seen at the upper right of IN fibril and propagates along the axis of the fibril over the next two minutes. The material does not have a large enough velocity component along the line of sight to show a distinct signature in either the line core Doppler shift or the line wings. It is also difficult to track through individual images, so no precise estimate of its velocity is possible. A crude estimate puts the apparent velocity somewhere between 10−20 km s^-1. We suggest this observed case is evidence of a flow aligned with the elongated axis of the absorption feature, that provides some support for the presence of a guiding magnetic field over an extended period of time.

The IN fibrils are also observed to support quasi-periodic transverse displacements and the remainder of this article focuses on them. Such motion can (and has previously) been interpreted as MHD wave behaviour, in particular the kink mode, and has been observed in network Hα fibrils (e.g. Morton et al. 2012a) and calcium active region fibrils (Pietarila et al. 2011; Jafarzadeh et al. 2017b).

4. Data analysis

4.1. Estimating the uncertainty

In order to make reliable measurements of the transverse displacement of fibrils and the associated errors, an estimate of the measurement errors associated with the data is required. The SST and CRISP provide high signal-to-noise observations with low photon noise. However, there are a number of other sources of noise associated with the data.

Distortion due to atmospheric seeing is a major source of noise. While the extended MOMFBD procedure and de-stretching of the data make a significant impact on reducing this source of data noise, uncertainty will evidently remain on both the intensity values and the physical location of features of interest. Ideally, some form of statistical re-sampling method would be employed, modelling the noise associated with the pre-processed data and re-running the data reduction pipeline numerous times in order to see how those uncertainties are propagated through the pipeline and impact upon the relative variations in spatial location and intensity. However, this is currently impractical due to the amount of processing time required for each data set.

CRISP is a double Fabry Pérot Instrument (FPI) set up to finely tune its spectral observations, allowing for high resolution spectral imaging. However, imperfections on the surface of the etalons of the FPI result in “cavity” errors (de la Cruz Rodríguez et al. 2015), leading to field-dependent shifts in the central wavelength of the transmission profile. This leads to uncertainties in measured intensity at high-gradient spectral regions such as the wings of the line, which are corrected intensity-wise by the pipeline and wavelength-wise on a per pixel basis when computing Doppler maps (or other situations where high accuracy in the wavelength is necessary) by applying the pipeline calculated cavity shifts. Such uncertainties should be extremely small at the core of the line where the gradient is low, but will nonetheless contribute.

In order to gain an estimate for the uncertainty associated with the data, we turn to multi-scale image processing techniques. We begin by assuming that we can estimate the uncertainty on the measured intensity values from the noise within the processed data and we would like to establish a relationship between the intensity and the noise. To do this we first create a “noise image” from the whole FOV of the data. The data is cropped to remove artefacts at the edges left over from the reduction pipeline.

The data is then filtered in space by applying a multi-scale filter utilising the à-trous algorithm with a 2D B3 spline filter (e.g. Stenborg & Cobelli 2003; Starck & Murtagh 2006), and then we apply unsharp masking to the highest-frequency scale. This process ensures the data is reduced to variations on the smallest spatial scale and the residuals can be taken to be indicative of the data noise (e.g. Olsen 1993). The root mean square (rms) of forty successive noise frames in time is then taken to estimate the standard deviation per pixel. Figure 4a shows the estimated rms noise for this data set. There is clear evidence of fixed pattern noise in the rms noise image. It is revealed as a grid pattern and is likely due to small differences between the different MOMFBD sub-fields. The presence of the grid pattern in this rms image is likely exacerbated due to one or two frames taken during poor seeing conditions, which have relatively worse resolution and a somewhat smoother intensity profile. Hence, there is less intensity variation across these regions and a reduction of the average noise in the regions. The grid pattern only occupies a minority of the pixels across the image and will likely only impact minimally on the following analysis of the data noise.

Figure 4d displays the joint probability distribution function (JPDF) of the average intensity (averaged over the same 40 frames) and the rms noise. The JPDF shows a discernible trend of increasing noise with intensity. In order to characterise this trend, the rms noise values are binned as a function of average intensity, and the 1D PDFs of rms noise show an approximate log-normal behaviour (e.g. Fig. 4b). The values of the mean and standard deviation of the log-normal rms noise are determined for each average intensity bin. Figure 4c shows the average intensity versus the mean values of rms noise and a quadratic function is fit to the mean values to establish a relationship between the two quantities. The found parameters give: (1)This equation will be used to provide our estimates of uncertainty associated with measured intensity. The resulting uncertainty is approximately 0.1% of the intensity, which is in line with expected noise levels estimated from theoretical arguments by de la Cruz Rodríguez et al. (2012) for CRISP-like FPIs operating on similar optical systems to the SST.

Fig. 4 Noise estimation process. The rms residual intensity or noise (natural log units) is obtained across the FOV panel a, which has the joint probability distribution shown in panel d when compared to average intensity. The bin widths used in panel d are 0.01 units in the rms noise and five units in the intensity. In panel b, an example of the 1D histograms of the rms noise estimate for the intensity bin 400−405 is shown. The over-plotted red curve in panel b is a fitted Gaussian and the vertical red dotted line signifies the location of the mean. Panel c demonstrates the relationship between average intensity and the mean value of the log rms noise estimate. The polynomial fit to the mean values is over-plotted in red, the equation of this curve (Eq. (1)) will serve as the estimate of noise. The same data points are over-plotted on the JPDF panel d as a red line.

4.2. Transverse motion

The IN fibrils are observed to undergo a quasi-periodic displacement transverse to the elongated direction. This motion is revealed in time-distance diagrams created from slits placed perpendicular to the observed axis of the elongated direction (Fig. 5).

The fibrils’ cross-sectional flux profile is approximately an inverse Gaussian, and the central location of the fibrils axis is defined as the location of greatest absorption, that is, the minimum of the flux profile. Hence, from the time-distance diagrams, the central locations of the fibril axis, x(t_i) = x_i, were found using the NUWT (Northumbria University Wave Tracking) code (Morton et al. 2016; further details are given Morton et al. 2012a; Morton 2014). NUWT requires the uncertainties in intensity (Eq. (1)) to establish the central location with sub-pixel accuracy and also provide meaningful uncertainties on this position (σ_i). Examples of the IN fibrils and the resulting tracked axis in the time-distance diagrams are displayed in Fig. 5, with the x_i’s forming a time-series of the transverse displacement. Time-series are obtained in this manner for multiple cross-cuts placed at equally spaced intervals along the axis of the fibril, enabling a picture to be built up of how the wave evolves as it propagates.

In order to characterise the apparent periodic nature of the displacements, the individual displacement time-series ((t_i,x_i), i = 0,...,n) are fit by least-squares (Markwardt 2009), using a model of the form, (2)where a₃ is the amplitude, a₄ is the period, and a₅ is the phase. The linear terms represent any potential long-term drift of the fibril over the course of the quasi-periodic motion.

Fig. 5 Examples of IN fibrils are shown in the left-hand column. White lines are over-plotted to show locations of cross-cuts taken along the structure. The longitudinal line serves as a guide line along the fibril axis. Right-hand panels: examples of time-distance diagrams from the IN fibrils, revealing the transverse displacements of the IN fibril. The white stars highlight where the measured central locations of the IN fibrils axis are from the fitting routine. Every 10th point is plotted for clarity.

4.3. Propagation speed

To measure the propagation speed of the transverse waves supported by the fibril, the time-series obtained at multiple locations along the fibril are cross-correlated. The time-series from the central cross-cut is used as a reference series and the preceding and proceeding series are cross-correlated with the reference, enabling a measurement of the time-delay between time-series. We note this procedure assumes that there is no significant curvature of the fibril perpendicular to our view of it. Any such curvature would lead to an underestimate of the distance that the wave propagates, leading to slower measured propagation speeds.

Shorter sections of the time-series are chosen for cross-correlation due to the possibility that neighbouring peaks in displacement are related to counter-propagating waves (e.g. Morton et al. 2012a). Additionally, upon examining the displacement time-series, there are data points that show relatively large variations in position from the local “average”. This feature is observed across all time-series for a particular feature and occurs when the data quality drops in that particular region, which we associate with seeing variations across the image. The prominent nature of these excursions can impact upon the cross-correlation results. Hence, a cautious approach is applied to selecting segments of the time-series to cross-correlate, avoiding these outliers. This is currently a rather subjective process and future studies will require more advanced objective methodology to counter such issues, for example, wavelet or empirical mode decomposition may provide this ability, or an initial cleaning of the data before feature measurement occurs. Once a suitable portion of the time-series has been identified, the series are cross-correlated and a polynomial function is fitted around the local maximum of the correlation function to achieve sub-cadence accuracy on the lag value. Figure 6 shows an example of the different time-series taken along a fibril used for cross-correlation.

In an attempt to determine the uncertainties on the lag values, a parametric re-sampling technique is employed. In brief, the idea behind the re-sampling methodology is to try and determine the distribution of an estimator in order to establish indicators of the estimators reliability, for example standard errors or confidence intervals. The parametric re-sampling makes use of an estimated or hypothesised model for the uncertainties in the original data, and propagates these errors via Monte Carlo simulation. This methodology eliminates the need for simplifying assumptions, which are required to obtain approximate analytic formulae describing how uncertainties propagate through complicated operations on the data. Here the estimator of interest is the lag value, and we would like to know how the uncertainties associated with the measured locations of the IN fibrils’ central axis influence the lag value. This methodology enables us to incorporate the heteroscedastic nature of the measurement uncertainties, σ_i, which we assume are normally distributed about the measured location, . The methodology is then to generate repetitions of the displacement time-series with different values of measurement noise and to cross-correlate these re-sampled series.

The following steps describe this implementation.

• 1.

  Generate random noise for each displacement time-series data point (t_i,x_i) from a normal distribution, .

• 2.

  Add this random noise to the original time-series, creating a re-sampled time-series with a different realisation of the data noise.

• 3.

  Generate five hundred replicates in this manner for each time-series.

• 4.

  Cross-correlate the five hundred reference signal replicates with its counterpart from different spatial locations.

• 5.

  Establish a distribution of time-lags from the five hundred correlations, calculating the mean lag value and its standard deviation (e.g. Fig. 7).

A linear function is then fit to the lag values as a function of distance along the fibril, with the gradient providing a measure of the propagation speed (Fig. 7). It is found that the lag values show little evidence of deviations from a straight line, with the linear fits having a reduced χ² ~ 1. In certain cases, the χ² is substantially less than one, potentially implying that the associated uncertainties with the lag values are over-estimated. However, at present we remain on the side of caution with our interpretation, and the constant time-lag with height implies there is little variation of propagation speed as a function of distance along the IN fibrils.

4.4. Magneto-seismology

In our visible inspection of the data in Sect. 3, it was noted that the relative inclination of the IN fibrils within the atmosphere was near impossible to establish from visual examination alone. We now outline the magneto-seismological methodology required to estimate variations of plasma parameters along the fibrils, which may provide information on whether the fibrils are lying horizontally in the chromosphere or not.

There is a substantial volume of literature dedicated to the use of kink waves as an inversion tool to diagnose the local plasma properties. In general, the theoretical development of the atmospheric magneto-seismology has contained very strict assumptions about the nature of the plasma (e.g. that it is hydrostatic), and the wave evolution (e.g. that there is no damping during propagation). There are attempts to relax some of these assumptions, for example Morton & Erdélyi (2009b), Ruderman (2011), and Soler et al. (2011)³, but this still does not capture the complex physical effects observed during the evolution of the kink waves along chromospheric wave-guides (e.g. Morton 2014). Despite these constraints, it is useful to compare observed wave behaviour to the current theory in order to aid interpretation.

Fig. 6 Set of time-series taken from a single IN fibril. Each time-series corresponds to a measurement of the central location of the fibril from a time-distance diagram. In this example, the time-series are measured at separations of 85 km taken along an IN fibril feature. The red dotted lines shows the time interval used for cross correlation.

Fig. 7 Propagation speed measurement. Top panel: lag values as a function of the time-series position. The zero time-series number corresponds to the middle time-series in the range. A linear fit is shown by the red line, the gradient of which determines the propagation speed. Bottom panel: typical distribution of time lags between two time-series from the re-sampling technique. The red line is a fitted Gaussian with the mean shown by the dashed red line. The mean and sigma values of the Gaussian are used as the lag and uncertainty on the lag respectively.

Verth et al. (2011) demonstrated that inversion of the properties (amplitude, phase speed) of a propagating kink wave along a magnetised wave-guide can inform us about the variation of the local plasma density and magnetic field strength. We recall that the WKB solution to the wave equation governing a kink wave propagating along a stratified flux tube under the thin flux tube approximation (Morton et al. 2012b) is given as⁴, (3)Here ξ(z) is the displacement amplitude as a function of distance z, c_k(z) is the kink phase speed, ω is the angular frequency, R(z) is the radius as of the flux tube as a function of height, and C is a constant. Starting from Eq. (3), we derive relations for the normalised density, radius, and magnetic field following Morton (2014), which can give an indication of the quantities’ relative evolution of the stratified flux tube, Here, ⟨ ρ(z) ⟩ is the density averaged over the internal and ambient plasma and ⟨ B(z) ⟩ is the local average magnetic field. These relations show that the above quantities rely solely on the amplitude (ξ) and kink phase speed (c_k). We emphasise that these relations only hold under an assumption of no wave damping and a hydrostatic plasma. As mentioned in Sect. 4.3, a linear function provides a reasonable fit to the slope of the time-lags, with no evidence that a higher order polynomial is required (within the current uncertainty estimates). As such, the propagation speed of the transverse wave can be considered constant along the wave-guide, which means the propagation speed terms in Eqs. (4)–(6) drop out, leaving relations entirely dependent on amplitude.

Finally, depending on the nature of the IN fibrils, the variation in magnetic field and density will be dependent upon the inclination of the fibril to the vertical (Sect. 3). In a hydrostatic atmosphere, a straight flux tube will have an exponential decrease in density, modified by any inclination. However, should the IN fibril be loop-like, the expected variation in density in a hydrostatic atmosphere is given by (7)assuming the loop is semi-circular (e.g. Dymova & Ruderman 2006), where z is the distance along the loop, H is the scale height, and L is the loop length. Further modifications to this relationship are present should the loop geometry deviate from this (e.g. elliptical; Morton & Erdélyi 2009a). These gravitationally stratified profiles for the variation in density will lead to upward propagating waves showing amplification and downward propagating waves displaying attenuation (Eq. (4)).

One must also recognise that any estimated density profile will be an average of internal fibril and ambient plasma densities, hence, variations of the ambient plasma with height can also contribute to the measured gradients (Morton 2014). If, for example, we were to assume that the ambient atmosphere around the IN fibrils takes on a profile of an approximately gravitationally stratified atmosphere, and the IN fibrils are features composed of non-hydrostatic chromospheric plasma that protrude into the corona (although not as far as spicules), then the external density may drop off rapidly compared to the internal density, leading to significant changes in wave amplitude.

Fig. 8 Histograms displaying the measured amplitudes, periods, phase speeds, and velocity amplitudes from the 28 IN fibrils measured.

5. Results

In this section we present the results from the measurements of the MHD kink waves observed in 28 IN fibrils. An overview of the measured wave properties is presented in the histograms in Fig. 8, showing the distribution of displacement, period, velocity amplitude, and propagation speed. The values shown in the histograms are the weighted means of the measured properties along fibrils, providing in essence a summary of the wave properties. The velocity amplitude is calculated utilising the standard equation, (8)where v is the velocity amplitude, and ξ_wm and P_wm are the weighted mean displacement amplitude and period of each feature, respectively. The mean, median, and standard deviation of these results are given in Table 1. The results are consistent with previous Quiet Sun measurements of network fibrils (e.g. Table 3 in Jess et al. 2015). We note that the gaps in the data may have an impact on our ability to measure waves with periods on the order of, and longer than, 300 s, although Morton et al. (2014) did not find many examples of such periods in both quiet and active region fibrils.

The results for the propagation speeds of kink waves along the IN fibrils suggest that, in general, they are comparable to the limited measurements of propagation speeds along network fibrils and mottles (Jess et al. 2015). The results are notably greater than those measured for Ca II slender fibrils in active regions by Jafarzadeh et al. (2017b). Moreover, there are a number of IN fibrils with propagation speeds that exceed 500 km s^-1, which are perhaps unexpectedly large. This is discussed further in the final section.

Fig. 9 Magneto-seismological inversions for the relative change in density along the IN fibrils. The profiles are separated whether the wave was measured to be propagating away from the apparent endpoint of the fibril (left) or towards it (right). The black dashed line is the expected relative change in density for a hydrostatic atmosphere with scale height 250 km.

The magneto-seismology inversions for the density profiles are shown in Fig. 9. To interpret the results, we work under the assumption that the fibrils are closed magnetic field lines, that is low lying loops, in line with our expectations from the visual inspection of the data in Sect. 3. The results are separated into two plots dependent on whether the wave is observed to propagate away from the observable endpoint of IN fibril (“upward” – left panel) or towards it (“downward” – right). To aid comparison of density changes along the fibrils, an exponential density profile for a hydrostatic atmosphere is over plotted, with a scale height of 250 km (Uitenbroek 2006). This represents the maximum change in density through the atmosphere for this value of scale height, with inclined and loop-like features having less variation as inclination to the vertical increases (Eq. (7)). The magnetic field variation can be inferred from these plots as it is proportional to the square root of the density variation (Eqs. (4) and (6)).

In general, the measurements reveal evidence for trends in the wave amplitude as they propagate along the IN fibrils, although the measurements have a high variability (with correspondingly sizeable error bars). The variation of amplitude suggests that the fibrils and ambient plasma possess some stratification of density and magnetic field, indicating that the IN fibrils are not merely horizontal features in the chromosphere. In both the “upward” and “downward” propagating waves, there are examples whose variation in amplitude is roughly in line with that expected from a hydrostatic density profile.

The blue and pink “upward” profiles would suggest the density increases along those fibrils. However, this is potentially a signature of wave damping, which effectively works against the density to attenuate any amplitude variation for waves propagating upwards along the wave-guide (the interplay of longitudinal inhomogeneities and wave damping was discussed in Morton 2014). There is no evident reason why all observed waves should not also be damped to some degree during propagation, with the inhomogeneous transverse structuring of the local plasma enabling resonant absorption (e.g. Terradas et al. 2010) or the driving of the Kelvin-Helmholtz instability (e.g. Antolin et al. 2015) to act upon the observed waves. For the upwardly propagating waves, this would imply the average density of the fibrils and external atmosphere are more stratified than suggested by the majority of the magneto-seismology profiles in Fig. 9. The potential signature of damping is also present in the “downward” propagating waves, with the damping potentially enhancing the effect of increasing density to attenuate the amplitude, as seen with the yellow and purple profiles.

Finally, there is potential evidence for additional wave amplification, the “upward” cyan profile (with an alternative explanation as a highly stratified fibril) or the gold “downward” profile. It is not known to us what may cause an additional amplification of the wave as it propagates, although a scenario can be imagined involving a suitable combination of internal and ambient density profiles, for example Morton (2014) reported rapid amplification of a propagating kink wave along a spicule and demonstrated the required density profiles that could recreate this behaviour.

The vast range of possible scenarios makes it difficult, at present, to extract actual gradients of the density along the observed waveguides. It is likely that an advanced inversion scheme is required to exploit the observed propagating kink waves, potentially utilising the power of Bayesian analysis and requiring additional information about the plasma. In the hope that such a scheme is developed, the trends have been fitted with an exponential profile in order to provide a measure of the variation (results given in Table 2). While we refer to this measure as a scale-height, we emphasise that the value will incorporate the influence of internal and external variations in density, as well as wave damping (or even amplification).

6. Discussion and conclusion

Here we have begun to examine the properties and dynamic behaviour of internetwork fibrils observed in Hα, with the long term goal of trying to understand their formation and the antecedent energy deposition. The relatively isolated nature of these features, as opposed to the numerous network fibrils that occupy the rosettes, may mean that identifying the mechanism(s) responsible for depositing energy is a simpler task. Before concluding, we highlight that the number of IN fibrils in this study is small and urge caution in assuming that these properties are representative of all IN fibrils.

In general, we find that the IN fibrils display similar behaviour to network fibrils in Hα. The fibrils are found to appear, disappear, and re-appear over periods of five to ten minutes suggesting regular but not continuous depositions of energy. The reason for their pattern of visibility in Hα is unclear at present, although it may be the result of the fibrils being in a near constant state of thermal non-equilibrium, undergoing cycles of heating and cooling (such cycles have been extensively studied to explain coronal loop behaviour, e.g. Klimchuk et al. 2010). The following is a crude description of how such a cycle could impact on visibility. An initial energy deposition occurs along the magnetic fields, heating the local plasma, potentially – and eventually – leading Hα visibility through delayed post Saha-Boltzmann extinction (Rutten 2016; Rutten & Rouppe van der Voort 2017). The heating could drive up-flows of denser material from lower atmospheric heights, which would also lead to enhanced Hα extinction (Leenaarts et al. 2012). The fibril material then cools and condenses, evacuating the magnetic flux tube; both cooling and density decrease leading to less Hα opacity. Another heating event occurs and the cycle repeats. For coronal loop modelling at least, this behaviour leads to over-dense loops, hence, larger than expected scale-heights.

Such a heating scenario may require the existence of a semi-stable chromospheric magnetic flux tube that supports denser plasma than its ambient environment, which may not be guaranteed on theoretical grounds in the partially ionised chromosphere (e.g. Martínez-Sykora et al. 2016). However, we have observed flows of dense (i.e. absorbing) material along the IN fibril sometime after the initial formation, which would indicate that a guide for the plasma exists. Moreover, we observe signatures of quasi-periodic, propagating transverse displacements of IN fibrils, which we interpret as MHD kink waves. The presence of kink waves requires an over-dense magnetic flux tube to be present, with the field orientation perpendicular to the direction of displacement. The mere fact fibrils are visible as Hα absorption features could indicate they are denser than their surroundings, with their appearance due to increased column mass along the line-of-sight that leads to a greater than average formation height (Leenaarts et al. 2012). The presence of the propagating kink modes along the elongated axis supports the idea that the observed absorption features outline magnetic fields.

At present there is some discussion of whether the observed chromospheric features are tracers of the magnetic field. Early work by de la Cruz Rodríguez & Socas-Navarro (2011) using Ca II 8542 Å spectropolarimteric observations suggested a certain degree of misalignment between fibrils and magnetic fields. Recently, more robust analysis implies the misalignment degree is much less than previously suggested, although the dispersion of misalignment angle increases in less magnetised regions (Asensio Ramos et al. 2017). Moreover, observations with helium 10830 Å find little evidence for misalignment (Schad et al. 2013).

On the other hand, Martínez-Sykora et al. (2016) demonstrated that 2D simulations of partially ionised chromospheric plasmas including the effect of ambipolar diffusion also show misalignment between the magnetic fields and temperature (density) structures. However, if these features undergo significant heating during their formation, as has been suggested for some spicules and RBEs (Pereira et al. 2014; Skogsrud et al. 2015; Rouppe van der Voort et al. 2015) and also long fibrils (Rutten & Rouppe van der Voort 2017), then it appears crucial for models to require the inclusion of a sluggish, non-equilibrium hydrogen ionisation and to take into account the history of heating events in the chromosphere (R. Rutten, priv. comm.)⁵. This would then naturally lead to a significantly greater fraction of ions and electron densities in the post-heating plasma (e.g. see post (inter)shock regions in Leenaarts et al. 2007), and a potential reduction in the influence of ion-neutral related effects. Additionally, 2D simulations are unlikely to capture the necessary physics for the heating of chromospheric phenomena, missing torsional motions that can input additional energy and momentum into the chromospheric plasma (Matsumoto & Shibata 2010; Iijima 2016). A recent investigation of the Bifrost simulation by Leenaarts et al. (2015) suggests that their model chromosphere shows a mixed picture of alignment.

Returning to the measurements made here, we have been able to measure the propagation speed of a number of kink waves along the fibrils. Somewhat surprisingly, we find no evidence for variation in the propagation speed as a function of distance. The propagation speed measurements allow us to place rather coarse constraints on the local Alfvén speeds, hence, magnetic field values associated with the fibrils in the upper chromosphere. The following equation describes the relationship between the propagation speed (kink speed), plasma quantities, and the Alfvén speed, (9)where, μ₀ is the magnetic permeability of free space, ζ = ρ_e/ρ_i, and the subscripts i, e refer to internal and external plasma quantities. Here we have assumed B_i ≈ B_e. For estimated chromospheric densities of 10^-9−10^-10 kg m^-3 along absorbing Hα features (taken from simulations of a model quiescent chromosphere; Leenaarts et al. 2012), this would imply chromospheric magnetic field strengths somewhere in the range of 3−200 G. The larger values of propagation speeds (>500 km s^-1) are significantly greater than previous measurements, and would correspond to magnetic fields of B ~ 25−200 G. While the lower end of this range seems feasible, it is unclear whether the larger values of field strength would exist in the IN chromosphere. These larger values are certainly at odds with the LOS field strengths observed in and around the IN fibrils’ apparent footpoints. However, we note that no correction for stray light has been attempted on the magnetograms, with stray light potentially leading to an underestimation of the field strength up to as much as half for isolated magnetic elements, due to the contribution to the Stokes V profiles of surrounding non-vertically-magnetised regions (e.g. the Milne-Eddington inversions of similarly processed data by Narayan 2011). Furthermore, the larger values given here are in excess of estimated magnetic field strengths for spicules at the limb (e.g. ~50 G Centeno et al. 2010), which are generally associated with the stronger network magnetic field. Moreover, recent results from Asensio Ramos et al. (2017) find median values of ~60 G from Ca II 8542 Å observations of plage regions, although their estimated distribution of values has a long tail to larger values.

The larger propagation speed measurements here typically have larger errors, which may account for some of the apparent spread in the distribution of speeds. An alternative may be that there are standing modes supported by the IN fibrils, which would lead to spurious measurements of fast propagation speeds. Similar results were obtained in spicule measurements (Okamoto & De Pontieu 2011) and are still unexplained.

Finally, combining the measured wave amplitudes with propagation speeds, we are able to apply current solar magneto-seismology theory to the observations. We do so bearing in mind that the theory is quite conservative in terms of what is assumed about the behaviour of the local plasma (Sect. 4.4). The data reveals that the amplitude of the waves is non-constant along the waveguides, showing evidence of amplification and damping as a function of distance. The amplification is readily explained in terms of longitudinal stratification of the local plasma, combining together changes in internal and ambient plasmas. Observations and simulations suggest fibrils rise out of the surrounding “bulk chromosphere” and protrude some way into the corona (leading in part to a larger column mass), hence, we expect that there will be large changes in the external density as the ambient atmosphere rapidly transitions from chromosphere to corona. This situation would lead to the ambient density variations dominating the observed variation in amplitude, and masking the influence of internal density stratification. However, at present, it is not clear how to disentangle the two density profiles from the observed variation. In addition, it is expected that the kink waves are subject to damping of some form, which will also be entangled in the observed amplitude profile. Some of the profiles shown in Fig. 9 support the idea of the changes in external plasma dominating their behaviour, with measured variations in amplitude suggesting density scale-heights similar to that of a vertical gravitationally stratified atmosphere.

In conclusion, we provide here the first study of internetwork fibrils. The IN fibrils display a pattern of repeated Hα visibility suggesting repetitive heating events along a quasi-stable magnetic flux tube. We also find evidence of flows along the longitudinal axis of the fibril and observe quasi-periodic transverse propagating displacements that we interpret as the MHD kink wave. We suggest that both the flows and waves indicate the correspondence of an underlying magnetic field with the absorption feature visible in Hα.

Acknowledgments

K.M. acknowledges support from Northumbria University Research Development Fund and the Royal Astronomical Society. R.M. is grateful to the Leverhulme Trust for the award of an Early Career Fellowship. All authors acknowledge IDL support provided by the Science & Technologies Facilities Council. The Swedish 1-m Solar Telescope is operated on the island of La Palma by the Institute for Solar Physics (ISP) of Stockholm University in the Spanish Observatorio del Roque de los Muchachos of the Instituto de Astrofísica de Canarias. The observations were taken as part of work supported by the SOLARNET project (www.solarnet-east.eu), funded by the European Commissions FP7 Capacities Program under the Grant Agreement 312495.

References

All Tables

All Figures

Fig. 1 Hα line core and blue wing images are shown in the left and centre panels respectively. The two major magnetic flux concentrations that form part of the network and associated chromospheric rosettes are evident. Right hand panel: photospheric magnetogram. The levels have been clipped to B ~ ± 200 G to reveal some of the weaker fields in the IN. Contours highlight absolute magnetic field strengths greater than the cut-off, red showing positive flux and blue negative.

In the text

	Fig. 2 Examples of internetwork fibrils. Left hand column: intensity images determined from the profile minimum, revealing the chromospheric absorption features. Middle panels: corresponding wing images; and last column: photospheric magnetic field strength clipped between ± 50 G.
In the text

	Fig. 3 Panel a: one of IN fibrils in the data set, which can be seen at (8, 20) Mm in the line core image in Fig. 1. Panel b: a parcel of plasma is seen at the upper end of the fibril, significantly broadening the apparent width of the feature. The dark clump can be seen to move along the apparent axis of the feature in panels c and d (from upper right to lower left).
In the text

Fig. 4 Noise estimation process. The rms residual intensity or noise (natural log units) is obtained across the FOV panel a, which has the joint probability distribution shown in panel d when compared to average intensity. The bin widths used in panel d are 0.01 units in the rms noise and five units in the intensity. In panel b, an example of the 1D histograms of the rms noise estimate for the intensity bin 400−405 is shown. The over-plotted red curve in panel b is a fitted Gaussian and the vertical red dotted line signifies the location of the mean. Panel c demonstrates the relationship between average intensity and the mean value of the log rms noise estimate. The polynomial fit to the mean values is over-plotted in red, the equation of this curve (Eq. (1)) will serve as the estimate of noise. The same data points are over-plotted on the JPDF panel d as a red line.

In the text

Fig. 5 Examples of IN fibrils are shown in the left-hand column. White lines are over-plotted to show locations of cross-cuts taken along the structure. The longitudinal line serves as a guide line along the fibril axis. Right-hand panels: examples of time-distance diagrams from the IN fibrils, revealing the transverse displacements of the IN fibril. The white stars highlight where the measured central locations of the IN fibrils axis are from the fitting routine. Every 10th point is plotted for clarity.

In the text

	Fig. 6 Set of time-series taken from a single IN fibril. Each time-series corresponds to a measurement of the central location of the fibril from a time-distance diagram. In this example, the time-series are measured at separations of 85 km taken along an IN fibril feature. The red dotted lines shows the time interval used for cross correlation.
In the text

Fig. 7 Propagation speed measurement. Top panel: lag values as a function of the time-series position. The zero time-series number corresponds to the middle time-series in the range. A linear fit is shown by the red line, the gradient of which determines the propagation speed. Bottom panel: typical distribution of time lags between two time-series from the re-sampling technique. The red line is a fitted Gaussian with the mean shown by the dashed red line. The mean and sigma values of the Gaussian are used as the lag and uncertainty on the lag respectively.

In the text

	Fig. 8 Histograms displaying the measured amplitudes, periods, phase speeds, and velocity amplitudes from the 28 IN fibrils measured.
In the text

	Fig. 9 Magneto-seismological inversions for the relative change in density along the IN fibrils. The profiles are separated whether the wave was measured to be propagating away from the apparent endpoint of the fibril (left) or towards it (right). The black dashed line is the expected relative change in density for a hydrostatic atmosphere with scale height 250 km.
In the text