Free Access
Issue
A&A
Volume 620, December 2018
Article Number A124
Number of page(s) 20
Section The Sun
DOI https://doi.org/10.1051/0004-6361/201833664
Published online 07 December 2018

© ESO 2018

1. Introduction

The solar chromosphere is made of a partially ionized magneto-fluid that is opaque in the core of Hα, conferring its characteristic fibrilar fine structure on disk and spicular appearance at the limb (e.g., De Pontieu et al. 2003, 2004, 2007; Rouppe van der Voort & de la Cruz Rodríguez 2013; Rutten & Rouppe van der Voort 2017; Rutten 2017). This layer is of special interest because it is not clear what are the most significant heating mechanisms among the numerous candidates such as acoustic and MHD waves, magnetic reconnection, microflares and spicules, that can explain the energy balance (Withbroe & Noyes 1977; Narain & Ulmschneider 1990).

In addition to Hα, one usually resorts to complicated diagnostics formed under non-local thermodynamic equilibrium (non-LTE) conditions such as the Mg II h and k and Ca II H and K resonance doublet lines, moreover affected by scattering and partial frequency redistribution (PRD). These have been used to observe the chromosphere and infer temperatures and velocities, while the Ca II 8542 Å infrared line has also been commonly used to probe chromospheric magnetic fields (see review by de la Cruz Rodríguez & van Noort 2017).

On the other side of the spectrum, sub-mm/mm continua have the advantage of being formed under LTE in the chromosphere, therefore the source function is Planckian and the intensity is a linear function of the electron temperature since the Rayleigh-Jeans limit applies (e.g. Rutten 2003). However, LTE does not hold to describe the ionization degree of hydrogen in the chromosphere where time-dependent non-equilibrium effects should be taken into account to compute opacities (Leenaarts & Wedemeyer-Böhm 2006). These are set by ther–mal bremsstrahlung processes in the quiet sun (Bastian 2002; Rutten 2017), namely the electron-proton (H0) and electron-neutral (H) free-free absorptions, having a quadratic dependence on wavelength. This in turn means that mm radiation becomes optically thick at increasing heights sampled by increasing wavelengths (Wedemeyer-Böhm et al. 2007).

In flaring active regions a contribution from gyrosynchrotron is expected (e.g., see review by Wedemeyer et al. 2016), whereas cyclotron radiation (gyro-resonance) would require unrealistically strong magnetic fields with strengths of 28 000 G and 3300–6000 G for this emission to be significant at 1.2 mm and 8 mm respectively (Brajša et al. 2009).

Radio observations are thus great alternative chromospheric probes, and can be used to test models of the solar atmosphere (Loukitcheva et al. 2008). The problem had been that such observations (e.g. White et al. 1992, 2006; Bastian et al. 1993; Loukitcheva et al. 2014; Iwai & Shimojo 2015; Iwai et al. 2016) have had low angular resolution (≫5″ ) for a long time, which is critical to observe the solar chromosphere that shows rapidly evolving fine structures at least down to the diffraction limit of today’s 1m-class optical telescopes.

Although that is promised to change with the advent of the Atacama Large Millimeter Array (ALMA; Wootten & Thompson 2009), providing for the first time a high spatial and temporal resolution view of the solar chromosphere in the radio, interferometric imaging is not a straightforward task. The existence of both extended and localized sources and their variability on very short time scales, and other practical issues, make it difficult to perform image reconstruction of the Sun (Shimojo et al. 2017). Even if that is achieved, ALMA observations will still benefit from coordinated observations with other ground facilities such as the Swedish Solar Telescope (SST; Scharmer et al. 2003) or GREGOR (Schmidt et al. 2012), and the upcoming Daniel K. Inouye Solar Telescope (DKIST), and space telescopes like Hinode (Kosugi et al. 2007) or the Interface Region Imaging Spectrograph (IRIS; De Pontieu et al. 2014).

Spectral lines can provide complementary information on the plasma parameters such as temperature, velocities and magnetic field. The latter is specially important because at least in the first years of ALMA solar observations (Cycles 4, 5, and 6), polarimetry is not commissioned yet. In the future, the study of the continuum circular polarization could also be an alternative method of determining the longitudinal component of the chromospheric magnetic field (e.g., Shimojo et al. 2017; Loukitcheva et al. 2017a).

Studies of synthetic mm images have been done in the past, but the combined diagnostic potential of ALMA mm-bands with visible and UV wavelengths remains largely unassessed. Wedemeyer-Böhm et al. (2005, 2007) were the first to study sub-mm/mm brightness synthesized from 3D radiation-MHD simulations of the quiet sun using an LTE equation of state (EoS) and radiative transfer treated in LTE and non-LTE. Their brightness maps show filamentary brightenings resulting from shock-induced thermal structure, and fainter regions in between, while the average brightness temperature (Tb) and relative contrast were found to increase with wavelength.

Loukitcheva et al. (2015) obtained Tb maps from a more realistic Bifrost simulation by Carlsson et al. (2016) which includes non-equilibrium hydrogen ionization. They reiterate that mm-brightness can be a measure of the chromospheric thermal structure at the height at which the radiation is formed as expected. However this height is not unique and unambiguous owing to the complex thermal structure. They find that Tb is a reasonable measure (albeit not perfect) of the gas temperature at the effective formation height, defined as the height corresponding to the centroid of the contribution function. Furthermore, their results indicate that although instrumental smearing reduces the correlation between brightness and temperature, the former can still be used to diagnose electron temperatures up to a resolution of 1″ with reasonable accuracy. Thus spatial resolution becomes a critical issue for observing the solar chromosphere with ALMA.

The first study combining actual solar ALMA data with IRIS data was carried out by Bastian et al. (2017) who find a correlation between Tb(1.2 mm) and the radiation temperature of Mg II h, with a slope smaller than one, highest in plage and lowest in umbra and quiet-Sun regions (see also Bastian et al. 2018). This is probably consequence of the local decoupling of Mg II due to the scattering component in its source function, although there could also be systematic differences in the formation heights of both diagnostics. In any case this study showed how complementary IRIS and ALMA observations are, and how useful their comparison can be, in line with predictions (Wedemeyer et al. 2016, and references therein).

Inferring physical properties of the plasma in the solar atmosphere from spectropolarimetric observations involves inverting the radiative transfer equation (RTE; del Toro Iniesta & Ruiz Cobo 2016, and references therein). In this paper we used the new STockholm Inversion Code (STiC; de la Cruz Rodríguez et al. 2016), which is capable of dealing with simultaneous, LTE/non-LTE, multiatom modeling of a large number of frequencies, to investigate whether mm bands could help constraining chromospheric temperatures while improving upon common inversion schemes, for example, using the Fe I 6301, 6302 Å Zeeman triplet, Ca II 8542 Å and H and K or Mg II h and k. To this end we performed multiple tests with a 3D radiation-MHD model and different combinations of spectral lines. We want to address the question: do radio observations help state-of-the-art inversions of the solar atmosphere?

2. Model

Temperature, velocities and magnetic fields have a clear imprint on the synthesized line strengths and shapes, therefore it is important to test a fairly large range of physical conditions not only to assess the response of the diagnostics to changes in the plasma parameters, but also to investigate where in the atmosphere the combination of multiwavelength inversions is most useful.

We used snapshot 385 of the radiation-MHD simulation of an enhanced network (“en024048_hion”) performed with the Bifrost code (Gudiksen et al. 2011). This snapshot has been extensively used in the literature to study the formation of chromospheric lines and mm-continua (e.g., Leenaarts et al. 2013a,b; Pereira et al. 2013; de la Cruz Rodríguez et al. 2013, 2016; Štěpán et al. 2015; Quintero Noda et al. 2016, 2017; Loukitcheva et al. 2015, 2017a).

The computational box has 504 × 504 × 496 grid points with a pixel scale of approximately 48 km in the horizontal direction, making a domain with 24 × 24 Mm2, and variable grid separation in the vertical direction from 19 km in the photosphere to 100 km in the corona.

The simulation includes optically thick radiative transfer in the photosphere and low chromosphere, and parameterized radiative losses in the upper atmosphere (Carlsson & Leenaarts 2012). It takes into account thermal conduction along magnetic field lines and an EoS that includes the effects of non-equilibrium ionization of hydrogen (Leenaarts et al. 2007). For additional references and a detailed description of this simulation we refer to Carlsson et al. (2016).

However, it lacks non-equilibrium ionization of helium and the effects of partial ionization on chromospheric heating by magnetic fields (Leenaarts et al. 2013b). Furthermore, comparisons between synthetic observables and real data by those authors showed that the simulated chromospheric lines tend to be relatively weaker and narrower. Even so, the overall agreement between the simulated features and the real data is good, implying that much of the relevant physics is included.

Rutten (2017) noted that Hα fibrils synthesized from this simulation are less numerous and opaque than real solar fibrils, which in turn should be as opaque or even more opaque at mm-wavelengths. This is also not how the synthetic ALMA images obtained by Loukitcheva et al. (2015) from the same Bifrost simulation snapshot look like. Therefore one has to keep in mind that the simulation has certain limitations.

The computational cost of inversions scales with the number of pixels and it is currently prohibitive to perform several cycles of inversions using different combinations of lines for the entire simulation box treating each of the 254 016 pixels individually since each one can take between ∼ 2–6 h to converge depending on the atmosphere, the exact numerical settings and the number of frequency points and depth nodes used. For this reason we restrain the analysis to a more modest box of 74 × 74 pixels (3.5 × 3.5 Mm2) that encloses one of the magnetic patches and it is large enough for some structure to be recognizable while providing sufficient statistics. This particular patch was chosen to contain both quiet-Sun like atmospheres and more strongly magnetic regions with a magnitude reaching up to 2870 G at the bottom of the domain. However, |B| is overall weak with a median value and median absolute deviation of |B ∼ 50(±30) G at z = 0 km, and a maximum of |B| ∼ 90 G at z ∼ 2000 km. In fact, strong magnetic fields (|B| > 1000 G) only occur in about 3% of the pixels. We do not find any significant difference in the inversion results depending on the pixel location on the selected patch, so we do not expect that to occur elsewhere in the whole box. Since neither the Ca II nor the Mg II lines are sensitive to temperatures much higher than 16 000 K (Carlsson & Leenaarts 2012), the coronal-like part was cropped out as well as the domain below z < −0.2 Mm, keeping only formation heights of the studied lines.

The two leftmost panels in Fig. 1 show the temperature stratification in the snapshot as function of optical depth. Throughout this paper we abbreviate the optical depth at 500 nm (log τ500) simply as log τ. The black square illustrates the selected box for the calculations. At log τ = 0 we clearly see the granulation in the photosphere with hot granules with temperatures between T ∼ 6200−7000 K, whereas the darker intergranular lanes are cooler (T ∼ 5800 K), and overall the average temperature is 6300 K. Going higher up to log τ = − 2, which corresponds to z ∼ 250−350 km, the average temperature drops to approximately 4950 K and the reverse granulation becomes visible featuring dark granules with temperatures below 4500 K and brighter intergranular patches with temperatures between (T ∼ 5000−6000 K). In the chromosphere a given optical depth can correspond to very different geometrical heights, depending on the underlying density stratification and, to some extent, the magnetic field from pixel-to-pixel. For example, log τ = −6 corresponds to a broad height range usually within z ∼ 1300−1800 km. The gas temperature at this depth shows a filamentary structure connecting the magnetic patches of opposite polarity with cold pockets of only T ∼ 2500 K plasma juxtaposed with hotter structures up to T ∼ 11 700 K.

For the purpose of comparing response functions only we also extracted a diagonal slice running through the two magnetic patches defined by [x0 = 50, y0 = 300, x1 = 450, y1 = 150] in pixel coordinates that is approximately 21 Mm in physical length in the simulation snapshot (Sect. 4.5).

thumbnail Fig. 1.

Temperature structure and emerging intensities from a Bifrost simulation. Left panel: gas temperature as function of optical depth; central panel: subregion selected for the inversions with temperature scale in units of 103 K; the magnetic field strength at log τ = 0 is represented by thick (2 kG) and thin (0.1 kG) green isocontours. Right panels: emerging intensities at selected wavelengths given in ångströms except where indicated; for each spectral line, the top and bottom images correspond to core and blue wing wavelengths respectively.

Open with DEXTER

3. Calculations

We performed the calculations using the STiC code (de la Cruz Rodríguez et al. 2016), which is based on the RH code (Uitenbroek 2001), to synthesize and invert the spectral data. STiC is capable of simultaneously inverting the Stokes profiles of multiple LTE/non-LTE lines (and the mm-continuum) including PRD effects using node parametrization. It uses third-order DELO-Bezier splines to integrate the RTE (de la Cruz Rodríguez & Piskunov 2013) and it assumes plane-parallel geometry and hydrostatic equilibrium. Further details on more recent updates on the code can be found in de la Cruz Rodríguez et al. (2018).

3.1. Forward calculations

The spectral lines and the radio continuum were synthesized column-by-column (1.5D approximation) using an LTE EoS while the statistical equilibrium is solved for non-LTE populations of a five-level plus continuum hydrogen atom model, a five-level plus continuum model for Ca II with PRD in the H and K lines and complete frequency redistribution (CRD) in the 8542 Å line, and ten-level plus continuum model for Mg II with PRD in the h and k lines and CRD in the subordinate UV triplet. For Fe I we used a fifteen-level plus continuum model assuming LTE and CRD. The radiation is assumed to come from the solar disk center.

In the solar chromosphere full 3D radiative transfer is important for strongly scattering lines (Leenaarts et al. 2009, 2012, 2013a). Moreover, PRD effects have to be taken into account. This has been done in full 3D for the first time in Sukhorukov & Leenaarts (2017) for Mg II h and k, and in Bjørgen et al. (2018) for Ca II H and K based on the same Bifrost simulation we use in this paper. This sophisticated treatment comes with an enormous increase in computational cost, and remains prohibitive in inversion codes. In rigor such computations would mostly be necessary in the cores of Mg II h and k between the emission peaks (Leenaarts et al. 2013b), which we neglect here. In practice this results in Mg II intensity maps of higher contrast than reality at core wavelengths. This is also true for the analogous Ca II resonance lines. For the same reason we do not include Hα, the classic chromospheric line, because a 3D treatment of the RTE is crucial due to its highly scattering nature (Leenaarts et al. 2012). This is not yet possible with current inversion codes. Finally, these effects are not as strong in the Ca II 8542 Å (de la Cruz Rodríguez et al. 2012) because it has a deeper formation height, so it can be modeled in 1.5D to good approximation.

We also included in the synthesis a weak photospheric line of Ni I 2814.4 Å in the IRIS spectral window which, according to de la Cruz Rodríguez et al. (2016), helps to infer vertical velocities in the photosphere.

This study is not exhaustive in finding the best lines to invert together with ALMA data. Rather, it concerns the most common diagnostics, and in particular those that can be observed by IRIS and SST, and how ALMA helps inverting those. Other combinations of lines that could have been interesting to investigate include Mg I b2 5173 Å and Na I D1 5896 Å, forming in the upper photosphere (Briand & Solanki 1995; Leenaarts et al. 2010), or even Ca II 8498 Å which is also temperature-sensitive in the chromosphere (e.g., Quintero Noda et al. 2017).

We did not take into account gyrosynchrotron emission at mm-wavelengths since it would only originate from very energetic positrons and electrons primarily at flaring sites (Wedemeyer et al. 2016), which is not the scope of the model atmospheres we were interested in.

The polarization is assumed to come solely from the Zeeman effect, which means that we do not include the circular polarization of the mm-continuum in the calculations, so we cannot investigate in this paper whether ALMA helps constraining magnetic fields.

In real observations it is customary to introduce an additional parameter to explain the observed line broadening that does not correspond to thermal or instrumental origin, but could be due to a non-thermal contribution or simply unresolved motions in the same pixel. This micro-turbulence broadening is assumed to have a Gaussian line-of-sight velocity distribution such that the total Doppler broadening is a sum in quadrature of the thermal () and turbulent components. We investigated two cases regarding the micro-turbulence velocity (νturb): (1) a constant νturb(log τ)= 2 km s−1 which has a line broadening effect, and (2) νturb(log τ) = 0.

Table 1

List of inversion schemes.

The near-UV lines were convolved with a Gaussian with full width at half maximum (FWHM) of 0.053 Å (to simulate IRIS observations), and the optical lines between 0.052–0.1 Å (SST/CHROMIS and CRISP). In the mm-continuum we used only eight points centered on the four sub-bands of each currently-available ALMA band: [1.20, 1.21, 1.29, 1.30] mm and [2.80, 2.86, 3.16, 3.22] mm. In principle real ALMA observations could achieve higher spectral resolution1 but we do not expect significant opacity or source function changes within shorter wavelength intervals. We did not add random noise to the data to ensure that the results do not depend on an arbitrary noise level, and we did not take into account the issue of different spatial resolutions, or in other words we assume each telescope perfectly resolves horizontal structures.

3.2. Inverse calculations

We tested several inversion schemes using different combinations of lines plus mm-continuum. For clarity we shall call them i1 (IRIS), i2+i7 (IRIS+ALMA), i3 (SST/CRISP), i4 (SST/CRISP+ALMA), i5 (SST+IRIS) and i6 (SST+IRIS+ALMA) as summarized in Table 1. The reason behind this separation is twofold: firstly we want to be able to discern the constraining power of ALMA apart from the lines, so we need to test inversion schemes with and without chromospheric lines, and secondly, such co-temporal and co-spatial data spanning the entire electromagnetic spectrum is not often available for a given solar feature of interest, but we still want to know to what extent does ALMA help constraining whatever common line observations one disposes, in particular when observations from other ground telescopes are not possible, whereas IRIS from space might still be able to record data.

Although we synthesized the full Stokes vector in each grid point, we fitted only Stokes I because in this paper we are primarily interested in the temperature stratification and not in magnetic fields. Therefore we run the inversions in the “NO_STOKES” mode. This makes the computations faster and less prone to failed convergence. Because most pixels are weakly magnetic and |B| quickly decays with height (see Sect. 2), the response of the spectral lines to the magnetic field is mostly small, with temperature being a far more important parameter, specially at the line cores which are formed higher, so this approximation does not critically affect the conclusions. Also we note that IRIS and SST/CHROMIS cannot observe in full-polarization mode, therefore in real data we would not fit the full Stokes vector of the Mg II h and k and the Ca II H and K lines, but only of the Ca II 8542 Å and Fe I 6301, 6302Å lines. In Sect. 4.4 we show that this approximation does not allow us to correctly account for the significant Zeeman splitting in the iron lines, especially above strong magnetic concentrations, but that has no impact on the inferred chromospheric temperatures.

The fitting weights were set in a arbitrary absolute scale which means the value of the merit function (χ2) is meaningless, but we assigned relatively more weight to the line cores compared to the optical and UV continuum. The mm-wavelength points were assigned even stronger weights (eight times higher than the line cores) to make sure the code always fits the mm-continuum.

We run the inversions in two initial cycles with an increasing number of nodes to improve convergence, and with regularization switched off. Bad pixels are eliminated and replaced with interpolated values at all optical depths using a 2D Gaussian convolution kernel. Then, we performed a third cycle with regularization. At all cycles the best-fit in each pixel was chosen by selecting the model with the lowest χ2 among the results of a few randomizations of the initial guess.

We extensively tested the dependency of the results on the placement and number of nodes, but it was impossible to find an optimal setup that works for any pixel in the atmosphere (see Sect. 5.2). On the other hand, it is necessary to minimize the differences between the numerical settings of each set of inversions to ensure that the results are dominated by physics, not by the inversion method.

This being said, we started with [6, 3, 1] nodes in [T, νlos, νturb] respectively, and we increased to [10, 6, 1] in cycle 2 and up to [14, 8, 1] in cycle 3 adding a modified Tikhonov regularization (de la Cruz Rodríguez et al. 2018) in T and νlos. The first cycle was run with equidistant node placement, but the two following cycles were run with a more sparse placement at higher optical depths and a more dense sampling at lower optical depths to ensure that the accuracy of the inverted chromosphere is not limited by the number of nodes.

We also performed one additional test running with 0 nodes in νturb on data synthesized with νturb = 0.

3.3. Response functions

To better understand the results of the inversions it is instructive to compute response functions (Rλ) as function of wavelength and optical depth since they reflect the different sensitivities of the diagnostics to changes in the plasma parameters (Beckers & Milkey 1975; Landi Degl’Innocenti & Landi Degl’Innocenti 1977), in particular to temperature as we are interested here. We numerically computed Rλ to temperature by introducing a small first-order perturbation in the atmosphere and comparing the resulting spectra with the previously obtained. This can be expressed by the following formula:

(1)

where Iλ is the intensity at a given wavelength and ΔT = 10 K is the magnitude of the temperature perturbation.

The exact shape of Rλ as function of wavelength and optical depth is obviously not only model dependent, but also dependent on the magnitude of the perturbation to some extent. The perturbation should be kept small, but large enough to beat the numerical noise (Mili´c & van Noort 2017). Numerical tests have shown that the value adopted is robust in the sense that it ensures the linear regime.

4. Results

4.1. Emergent spectra

The emerging intensities from the model resemble the temperature maps at a given height (see Fig. 1). For example, the granulation in the photosphere is very clearly seen at wing-wavelengths of Ni I 2814.4 Å and Fe I 6301.5 Å. The layers above can be probed for example with Ca II 8542 Å where we see the reversed granulation at λ0 − 1Å, but also if one scans sufficiently far way from the core of Ca II H and K and Mg II h and k. The images in the cores of the K and k lines look very similar, showing a somewhat higher layer of the atmosphere, and correlate with the temperature maps between log τ ∼ [−6.0, −5.0]. It has been suggested that the range of formation heights of the continuum at 1.2 mm and 3 mm is known to be similar to that of the cores of Mg II h and k (e.g., Bastian et al. 2017, and references therein), and in fact one can see that our synthetic maps show analogous structures.

Interestingly we also found a correlation between Tb(1.2 mm) and the radiation temperature of the core of Mg II h with a slope smaller than one as in Bastian et al. (2017, 2018) from actual observations. Here we find Pearson’s correlation coefficients between 1.2 mm and Mg II temperatures approximately equal to 0.78 for the k line and 0.79 for the h line.

The average spectra from this simulation patch shows Ca II H, K and 8542 Å in and UV triplet in absorption, whereas the h and k lines of Mg II are in emission with double peaks. The mean brightness temperatures (with standard deviation) in the radio are Tb(1.2 mm) = 5800 (± 880) K and Tb(3 mm) = 7250 (± 1140) K. These values are consistent with recent ALMA single-dish mapping observations of the quiet Sun. White et al. (2017) reported Tb(1.3 mm) = 5900 (± 100) K and Tb(3 mm) = 7300 (± 100) K, while Alissandrakis et al. (2017) obtained Tb(1.3 mm) = 6180 (± 100) K and Tb(3 mm) = 7250 (± 170) K.

4.2. Inverted temperatures

Figure 2 shows the results of six inversion schemes at five selected optical depths at which we compared model temperature versus inverted temperature. These maps correspond to the final inversion results and were not smoothed to remove bad pixels. We show the Pearson’s correlation coefficient as a rough measure of association between Tmodel and Tinv. This figure is complemented by Fig. 3 which illustrates the kernel density distributions2 of temperature residuals as function of optical depth.

thumbnail Fig. 2.

Inverted temperatures compared to the model. Each column corresponds to a given optical depth and each row to a different inversion scheme (not smoothed) which is to be compared to the model in the top row. The colorbars are in units of 103 K. The three markers in the top row correspond to the spectra in Fig. 6. r is the Pearson’s correlation coefficient.

Open with DEXTER

thumbnail Fig. 3.

Temperature residuals at selected optical depth nodes. Each panel compares the probability density distributions of the temperature residuals of two inversion schemes (with and without mm-continuum), with the median indicated by the white dotted lines while the 16th and 84th percentiles are shown by dashed lines.

Open with DEXTER

Overall, we find major performance differences most notably in the chromosphere. In qualitative terms we find: (1) a poor temperature reconstruction in the low photosphere in i1 and slightly better in i2, but very good in the other schemes; (2) a better inversion of the chromosphere when mm-bands are included in i2, i4 and i6; (3) a generally better inversion at all heights when more spectral data is used together, but (4) no substantial improvement in the low chromosphere when mm-wavelengths are added to schemes already solving for Mg II h and k and Ca II H, K and 8542 Å simultaneously. The first point was expected since i1 does not contain any information of the deeper layers of the atmospheres, whereas i3–6 contain the photospheric Fe I lines. Note, however, how the temperature map at log τ = −5.6 of i2 shows much more structure compared to i1 which tends to overestimate temperatures in this layer. The choice of this particular optical depth is not fortuitous. Rather, it corresponds to the average optical depth of the maximum response to temperature of the 3 mm band as we will show in Sect. 4.5.

The difference between i3 and i4 in this layer is even more striking. The temperature in the higher chromosphere above log τ ∼ −5.0 is mostly unconstrained when we use only Ca II 8542 Å, but when mm bands are added we suddenly gain access to information that is critical to reconstruct the stratification at these heights. However, at log τ = −4.0 we see that the i4 map is a lot more noisy than i3 which is a consequence of difficulties in simultaneously fitting both mm-bands and the wings of Ca II 8542Å.

We obtained generally better inversions at all optical depths when we combined visible and UV lines as depicted in the i5 row. The region spanning from the middle of the photosphere to the low chromosphere is more accurately reconstructed than in any of the other schemes, which is shown by the much smoother and noise-free temperature maps at, for example, log τ = −2.0, −4.0, −5.2. Higher up in the chromosphere at log τ = −5.6 we also see that the map shows more structure compared to i1 and i3, and it has cleaned up most of the “inversion noise” seen in i1. This means that the Ca II H and K lines offer important constraints in these layers. However, there are still some patches of systematically overestimated temperatures (by more than 1500 K). These are greatly reduced in i6 with the addition of mm bands, but not always.

4.3. Error analysis

We now define the relative temperature error (δT) and the absolute error (Θ) of the inversion as:

(2)

(3)

where TBIF is the model temperature at a given optical depth.

Figure 4 shows more quantitative estimates of the error in the inversions. In the top panel we plotted the Pearson’s correlation coefficient between TBIF and Tinv as a function of optical depth. We note that r is aspatial in the sense that it disregards the spatial distribution of the pair of values and assumes that they are independent and non-autocorrelated, which is not true. As a consequence, normal tests of significance cannot be performed. In spite of this, r is still a useful and simple measure of association. Since r is very affected by strong outliers we made a 4σ-clipping on Tinv a priori. The exact value of r is not so important, but rather how it evolves with height. Qualitatively this behavior is what one would expect after looking at the results in Fig. 2. We find that r tends to decrease with increasing height. These diagrams are specially useful in revealing the loss of correlation for optical depths below log τ ∼ −4.5 which is more severe in the model i1 and i3, but more moderate in the other schemes. We find that schemes with mm-wavelength points maintain a higher correlation with increasing heights compared to the ones without them. However, in the upper photosphere this behavior reverses since the inverted maps at those heights are typically more noisy (see i2 and i4 in Fig. 2).

thumbnail Fig. 4.

Temperature error of the inversions and correlation coefficient as a function of optical depth.

Open with DEXTER

In the bottom panel of Fig. 4 we show δT as a function of optical depth. The solid and dashed lines are the median of δT for schemes with and without ALMA respectively, whereas the filled contour is a robust measure of its spread given by where IQR is the interquartile range3 and erf is the Gauss error function. The dotted vertical lines delimit a the range of optical depths that contains 95% of the response functions at 1.2 mm and 3 mm bands (see Fig. 10, discussed further on).

The plots confirm that the error of the inversion is usually smaller in the photosphere and centered around zero, and larger in the chromosphere where it begins to show a systematic deviation with Tinv being more often underestimated. We find that the median absolute error between log τ = [−4, −2] is approximately 1–2% in all schemes. In the photosphere at log τ > −1 i1 and i2 have a much larger error (>50%) as expected, but the remaining schemes have δT ≲ 3%. Between log τ ∼ [−6, −4.5] is where we find the most interesting differences in the error distributions. For example, comparing i1 to i2, we see that the median δT oscillates around zero but with a lower dispersion in i2 than in i1 with an IQR that is almost twice as narrow. In this range, i2 has a consistently smaller median absolute error than i1 by 1–5%. Below log τ = −6, δT greatly increases in both schemes and the Tinv are not trustworthy. The difference in δT between i3 and i4 in the chromosphere is remarkable. While in i3 the median absolute error peaks at ∼ 250%, in i4 δT it is typically ≲15%. As we had seen in Fig. 2 this is a major improvement. This is not as clear if we compare i5 to i6 because the dispersion in δT is similar (≲10%). Between log τ ∼ [−6, −4] the difference of median absolute errors of i5 and i6 is less than 0.5%.

These numbers do not tell the complete picture of the error characterization and inversion performance of the schemes since they are robust measures of scale and location which are blind to the tails of the distributions. In fact, despite the difference in the median errors in the chromosphere being arguably small (except for i3/i4, see Fig. 3), the addition of mm-wavelengths is very useful to remove the outliers of the error distributions and substantially reduce the inversion noise.

Figure 5 shows the joint probability densities of temperature for the six inversion schemes at two selected optical depths in the chromosphere. On the top side of each panel we show the respective histograms. This figure emphasizes that the differences accentuate with decreasing optical depths, while the dispersion in Tinv is reduced with the inclusion of mm- data in the inversions, specially in i1 and i3. This is not only shown by the more compact contours of i2, i4 and i6 around the locus Tinv = TBIF, but also by their corresponding temperature histograms which, with respect to the histograms of TBIF, have smaller relative entropy than i1, i3 and i5.

thumbnail Fig. 5.

Joint probability densities of the temperature in the Bifrost simulation and inverted temperatures at two different optical depths. Each panel compares two inversion schemes (color code as in Fig. 4), and the dashed and solid lines enclose 90% of the distributions. The side plots show the corresponding histograms: top row: i1 vs i2, middle: i3 vs i4, and bottom: i5 vs. i6. The black dotted lines are the TBIF = Tinv locus.

Open with DEXTER

4.4. Examples of fitted spectra

Figure 6 shows three examples of fitted spectra at the different locations marked in the top row of the mosaic in Fig. 2. Overall the spectra are well reproduced with the fit residuals being usually slightly larger in the line cores than in the wings.

thumbnail Fig. 6.

Example of fits to the spectra from the Bifrost simulation snapshot. These correspond to the points marked in Fig. 2. Panels ac: overview of Mg II h and k, UV triplet lines; panels di: zoom into the cores of Mg II h and k, UV triplet; panels j, l, n: Ca II K; panels k, m, o: Ca II 8542 Å; panels p, r, t: Fe I 6301, 6302 Å; and panels q, s, u: ∼ 1–3 mm continua. Bottom panels: model temperature (dotted line) versus inverted temperatures.

Open with DEXTER

For example, on top of the kilogauss magnetic element (red diamond marker), we obtained very good fits except in the h and k cores. There is also a slight misfit in the cores of the Fe I lines because we did not include magnetic fields in the inversions, therefore the code cannot correctly reproduce these profiles simply with temperature and velocity changes because they also have a strong response to the magnetic field strength. Although we obtained an excellent fit to the Ca II 8542 Å line in this column, the absolute temperature residual in i3 is approximately Θ ∼ 114 K on average between log τ ∼ [−4.8, −3.0] which is larger than the residuals in i4 (Θ ∼ 80 K), or even i1 and i5 which are accurate to (Θ ∼ 30−40 K). Adding the mm points in i4 greatly improved the inversion compared to i3 from log τ = −4.8 and it is in fact the scheme that performed the best at the nodes log τ = −5.2 and log τ = −5.6. At these locations the schemes that performed the worst in the chromosphere were i2 and i6 with an average error of Θ ∼ 1750 K between log τ ∼ [−6.0, −5.0] despite the mm-continuum being well fitted.

At the weakly magnetic location (|Bmax| ∼ 128 G) given by the blue square, a 150 K bump above the average value in Tmodel at log τ ∼ −4.8 originates two emission peaks on both sides of the absorption core of the K line. The Mg II spectra are best fitted using the i1 and i5 schemes but they are not the ones that best recover Tmodel across the atmosphere. In fact, i6 is the scheme that performs the best overall (at least up to log τ ∼ −6.4) with an average error of Θ ∼ 80 K between log τ = [−5.0, −4.0] and ∼ 250 K at log τ = −5.4, but it fails to reproduce the temperature peak at log τ = −5.6 which is underestimated by ∼ 1630 K. The scheme i1 was the one that best recovered this feature with an error of Θ ∼ 650 K but it shows a much larger error Θ ∼ 430 K between log τ = [−5.0, −4.0]. In this example, i2 clearly provides a more accurate temperature at log τ = −5.2 than i1 having reduced the discrepancy from ∼ 500 K (in i1) to ∼ 300 K, which is also approximately the same residual level of i6. In this case the Ca II 8542 Å line is also very well fitted both with i3–6 but the addition of mm-points and/or extra lines is necessary to improve the accuracy in the layers above log τ ∼ −4.8. We see that in i4 the temperature between log τ = [−6.0, −5.0] is much better constrained, but the addition of the mm-continuum made it overestimate the temperatures in essentially all the upper photosphere between log τ = [−4.2, −3.0].

Finally, at the dark inter-granule location marked by the green cross we obtained excellent fits to all the spectral points with all the schemes except i2 which was not able to properly fit the continuum between the h and k peaks. Despite that, the accuracy of the inversions in the chromosphere is very different across the schemes with i1, i3 and i5 performing worse in general than their counterparts with added mm-points. At log τ = −5.2, for example, i1 has a residual of Θ ∼ 500 K, whereas i2 has a residual of Θ ∼ 300 K, and i5 has a residual of Θ ∼ 900 K whereas i6 has a residual of Θ ∼ 400 K. This emphasizes that a good fit to a single or even several spectral lines does not guarantee a high-fidelity reconstruction of the corresponding chromosphere due to the high degeneracy of the problem. Increasing the wavelength coverage by adding the mm-continuum did help retrieving a more accurate representation of the real temperature profile, but there is still a systematic underestimation of a few hundred kelvin. We should note that, not only at this particular location but also in general, the scheme i4 using only Ca II 8542 Å and ALMA to infer the properties of the chromosphere provides a temperature estimate that is of the same order of accuracy of i6 which makes use of a lot more lines. However, the larger combination of lines in i5/i6 provides additional constraints on νlos (see Fig. A.3), so overall we consider the latter schemes more powerful.

4.5. Response functions

Figure 7 shows the response functions to temperature perturbations of selected frequencies at the locations marked in Fig. 2. The plots show that Rλ is stronger in the line wings but maximal at higher optical depths, and weaker in the cores where it reaches lower optical depths. We also see that only a few wavelength points close to the line core are sensitive to the temperature variations in the chromosphere. In addition, none of the wavelength points is specially sensitive in the narrow temperature minimum/steep rise layers usually around log τ ∼ [−4.8, −4.0]. This is why we only obtained an accurate Tinv at these depths when we fitted multiple lines (i5 and i6) simultaneously, which together have more constraining power in this region.

thumbnail Fig. 7.

Normalized response functions to temperature perturbations at selected frequencies and pixels in the simulation box. Top row: response functions to temperature of Mg II h and k and Ca II 8542 Å at the same locations given in Fig. 6.; the solid line traces the maximum at each column. Bottom row: normalized responses at 1.2 mm (solid line) and 3 mm (dashed line), and model temperature (dotted line) at those locations.

Open with DEXTER

For example, at the location marked by the blue square, the h and k lines of Mg II are the ones whose Rλ is sensitive to the widest range of optical depths from log τ ∼ −2 to approximately log τ ∼ −6.0, whereas the core of Ca II 8542 Å is sensitive to lower layers with a maximum response up to log τ ∼ −5.2 (similar to Ca II H and K, not displayed). The response of the 1.2 mm-continuum has an absolute maximum at log τ ∼ −5.3 and a secondary peak at log τ ∼ −4.9, but the response of the 3 mm band is single-peaked at log τ ∼ −5.6. The response of the mm bands is practically contained within log τ ∼ [−6, −4], which encloses the range where we find accuracy improvements with the addition of those wavelengths.

In Fig. 6 we showed that none of the schemes was able to recover the steady temperature rise in the higher chromosphere above log τ ∼ −5.4 at the red diamond location, even though the spectra are very well fitted. This can be explained by the response functions at this particular location which occur deeper in the atmosphere than the average. The response of the Mg II and Ca II lines for log τ≲ −5.0 is negligible, but the 3 mm continuum has its peak at log τ = −5.4. Therefore, Tinv is poorly constrained beyond log τ ∼ −4.8 in the i3 scheme since the response to perturbations in higher layers is practically zero, and only a few points close to the very line core convey information of the layers close to log τ ∼ −5.0. The inclusion of the mm-continuum improved the inverted temperature just until log τ = −5.4 above which the constraints are very weak.

The diverse density and temperature structure in the simulation generates corrugated surfaces with overlapping distributions and significant responses at distinct heights. To further illustrate this point, we computed Rλ along a diagonal cut running through the two patches of opposite polarity.

Figure 8 shows the temperature stratification along the 2D cut. On top of the temperature map we overlay the maximum and the half-maximum level of R1.2 mm and R3 mm. In the bottom panel we show the magnetic field strength for context. We see that the optical depth at which R3 mm is maximal is always lower than the maximal R1.2 mm, but we also see that these surfaces are highly corrugated depending very much on the local atmosphere from column to column. For example, the maximum R1.2 mm spans two units of optical depth and has two peaks at log τ ∼ −4.6 and log τ ∼ −5.2. By fitting a Gaussian mixture model to the corresponding distributions of geometrical heights we find a main peak at z1 ∼ 700 ± 80 km and a secondary one at z2 ∼ 1200 ± 270 km, with 64% of the maximum R1.2 mm lying below 1000 km. The maximum R3 mm does not so evidently show a double peak in optical depth scale, but it does so in z-scale with peaks at 780 ± 84 km and 1380 ± 186 km, although in this case for 79% of the columns R3 mm peaks at z > 1000 km.

thumbnail Fig. 8.

Temperature, response functions to temperature at mm-wavelengths and magnetic field strength in a 2D-cut of a Bifrost simulation. Panel a: temperature (clipped at 15 000 K) as a function of optical depth; the solid and dashed lines are the maximum response in each column at 1.2 mm and 3 mm respectively (smoothed for visibility), and the filled contours enclose the Rλ at half-peak; the peak response the core of Mg II h is shown for comparison (dash-dotted). Panel b: magnetic field strength for context; panels IIV: normalized response functions of the 1.2 mm- (solid) and 3 mm- (dashed) continuum and temperature profile (dotted line) at the selected columns in panels a,b).

Open with DEXTER

The bimodality at 1.2 mm seems to be correlated with the underlying magnetic field since R1.2 mm tends to peak at lower optical depths in columns with strong kilogauss magnetic fields which can be found, for example, around x ∼ 4.5 Mm and x ∼ 13 Mm, but not always. We also see that over the magnetic field concentrations R1.2 mm and R3 mm are more localized in height, whereas in more quiet areas (e.g., x ∼ 14−18 Mm) the distributions are broader. These aspects are illustrated on the right panels in Fig. 8 that show how the response varies depending on the location along the slice.

Figure 9 compares mass density and electron-proton number density at two locations marked in Fig. 8, along with the respective responses to temperature perturbations of the 1.2 mm and 3 mm continuum and the core of Mg II h. These UV and radio frequencies respond differently according to changes in the opacity. The chromosphere above strong kG-magnetic elements is evacuated, that is, it has lower mass density than more quiet locations which makes the Mg II h and k lines respond to perturbations in deeper layers than the average. On the contrary, the mm-continuum, responds to changes in relatively higher layers because the dominant source of opacity χepnenp, where ne and np are the electron and proton densities, respectively, is greatly increased compared to more quiet locations. In weakly magnetic locations the difference between the sensitivity of the h core and the mm-continuum is much larger, with the former responding to changes in higher layers at log τ ≳ −5.6, whereas the 1.2 mm, for example, has maximum response at approximately log τ ∼ −4.6 at the selected location.

thumbnail Fig. 9.

Mass density, electron-proton number density and response functions to temperature of the mm-continuum and Mg II h. Top: mass-density and electron-proton number density at the magnetic (I) and non-magnetic locations (IV) marked in Fig. 8. Bottom: normalized response functions of the 1.2 mm (solid) and 3 mm-continuum (dashed) and Mg II h core (dashed-dot) at location I (left) and IV (right).

Open with DEXTER

There are significant differences between the layers contributing to the response at the core of the h line and the mm-continuum which accentuate in areas of lower magnetic field strength. Above kG-magnetic-field concentrations the response functions of the 1.2 mm and 3 mm continua generally peak at hotter (T ∼ 6500−8000 K) layers with higher electron densities than in more quiet areas (see Fig. A.4).

Figure 10 shows the distributions of maximum response functions to temperature perturbations for the selected cube from the Bifrost simulation. In the mm-continuum this corresponds to selected mm-wavelengths whereas for the lines it refers to some wavelength at the line center where Rλ occurs at the lowest optical depth, that is, Doppler effects have been taken into account. Apart from the ALMA bands used in the inversions (1.2 and 3 mm) we also computed for completeness Rλ at 4.5 mm and 8.5 mm that would correspond to future ALMA bands 1 and 2 under development at the moment of writing. These would allow the study of the upper chromosphere and transition region according to, for example, Loukitcheva et al. (2015), Wedemeyer et al. (2016). We find that the distributions of maximum Rλ move to higher heights with increasing wavelength as expected, although there is considerable overlap between wavelengths.

thumbnail Fig. 10.

Distributions of the maximum response functions to temperature perturbations at selected wavelengths. For the spectral lines we show the distributions at the wavelength where Rλ attains the greatest height in each pixel. The values shown on the right correspond to the median, the 2.5th and the 97.5th percentiles.

Open with DEXTER

In the selected portion of the Bifrost simulation, R1.2 mm and R3 mm peak at approximately the same range of optical depths as the cores of the Ca II and Mg II lines with 95% of the maximum R1.2 mm lying between log τ = [−5.8, −4.5], while the maximum R3 mm lies between log τ = [−6.0, −5.1].

The 1.2 mm continuum responds to temperature perturbations higher up than the core of Ca II 8542 Å (and the H and K lines) in about 46% of the pixels, whereas for the 3 mm continuum that amounts to 82% of the pixels. The core of Mg II lin es forms slightly higher than that of the Ca II lines and generally above the 1.2 mm continuum since the latter responds to temperature perturbations above the h (and k) line in only 30% (and 20%) of the pixels. In 58% and 49% of the pixels the 3 mm continuum responds higher than the h and k line respectively. In the majority of the pixels (94%) the 1.2 mm continuum responds to temperature changes higher up than the UV triplet between the h and k lines.

These results show that these ALMA bands are sensitive to temperature changes from just above the temperature minimum up to the middle-high chromosphere. The range of optical depths log τ ∼ [−6.0, −4.5] is where we find the most significant improvements with the addition of the mm-continuum to the inversions not surprisingly . It is clear now why Tinv substantially deviates from Tmodel from optical depths below log τ = −6.0 since temperature perturbations in the top of the chromosphere have no impact on the mm-spectra, although that is not the case for the Mg II h and k spectra whose Rλ has a longer tail to even lower optical depths. However, Mg II h and k (and Ca II H and K) already show a strong scattering decline at these depths meaning that the core of those lines carries little information on the local plasma temperature.

4.6. The degeneracy between temperature and micro-turbulence

In multiatom inversions it is, in principle, possible to separately determine both the thermal and microturbulent velocity components because all elements share the same turbulent broadening but undergo different (mass-dependent) thermal broadening. In our inversion schemes T and νturb are fitted simultaneously.

Figure 11 shows the empirical cumulative distributions functions (ECDFs) of micro-turbulence velocity residuals for the different inversion schemes. We find that νturb is more often underestimated than overestimated in all schemes, but the constraints significantly depend on the number of the fitted lines of different atoms. For example, in i1 νturb is underestimated in 56% of the pixels, and about 83% of pixels have a misfit smaller than 0.5 km s−1 . In contrast, in i5 νturb is underestimated in 43 % and overestimated in 38%, meaning that 19% of the pixels have a residual identically zero. In i5, 89% of the pixels have a misfit smaller than 0.5 km s−1.

thumbnail Fig. 11.

Empirical cumulative distribution functions of the micro-turbulence velocity residuals for different inversion schemes.

Open with DEXTER

This is interesting because the temperatures are also more often underestimated (see Fig. 3) in the chromosphere between log τ ∼ [−6.0, −4.5] where the cores of the Mg II and Ca II lines respond to temperature perturbations, so it seems that the underestimated temperatures in these layers are not due to overestimated νturb. Furthermore, since the absolute residual in micro-turbulence is mostly smaller than ∼0.5 km s−1, this parameter is not expected to have a large impact on the temperature accuracy.

However, we see that schemes fitting the mm-continua (i2, i4 and i6) typically result in larger underestimation of νturb than their counterparts without them, despite improving the inverted temperatures between log τ ∼ [−6.0, −4.5].

To further investigate whether our temperature accuracies were limited by the micro-turbulence degeneracy, we synthesized the same spectra setting νturb = 0, repeated the inversions for the schemes i5 and i6 with exactly the same setup, and compared to the results presented in Sect. 4.2.

Figure A.1 shows the temperature residuals of the inversions for the schemes i5 and i6 with νturb = 0 at selected optical depth nodes. We find that the temperature residuals are essentially the same in the two inversion modes, that is, having a constant νturb as an extra free parameter does not greatly impact the determination of temperature. This is because νturb itself is well constrained (better than 0.5 km s−1) given the number of lines of different atoms that we simultaneously fitted in i5/i6, although this could not be tested for the other schemes. In what regards i5 vs i6, we still found that the mm-continuum provides more accurate chromospheric temperatures when there is no micro-turbulence in the forward calculations.

5. Discussion

5.1. Temperature constraints from multiwavelength inversions

We have used a patch of a Bifrost simulation snapshot to investigate whether the addition of mm-spectra to current inversion schemes improves the reconstruction of the temperature stratification upon solely fitting the profiles of some visible and UV spectral lines.

We find that mm-wavelength bands can potentially be used to constrain chromospheric temperatures in combination with other diagnostics. This is especially noticeable when comparing the schemes i3 vs i4 which work as “control schemes” to test whether the mm-continuum has any impact on the inversions since, has shown in Fig. 2, the higher layers of the atmosphere are completely unconstrained if one uses only Ca II 8542 Å. Although we were able to accurately reproduce such line profiles with i3, the addition of mm-wavelength points in i4 enable us to invert the temperature structure in the chromosphere above log τ = −5.0 where Ca II 8542 Å loses sensitivity.

The mm-continuum also allows a more accurate inversion of the Mg II lines. The 1.2 mm band is especially helpful in setting the gradient between the temperature minimum and the chromospheric plateau, while the 3 mm continuum, which is formed slightly above, provides stronger constraints at higher layers of the chromosphere where the h and k line-cores are only weakly coupled to the local conditions. This happens because both 1.2 mm and 3 mm emission have stronger responses to temperature changes at approximately the same heights as the cores of the Ca II and Mg II doublet lines in a relatively narrow layer of the chromosphere, but the source function of the mm-continuum can be described in LTE by the Planck function (Bλ) which means we are able to retrieve more information of the local conditions in the higher layers of the chromosphere.

Our results show that the combination of IRIS and ALMA observations is a powerful tool to diagnose plasma temperatures practically across the entire atmosphere up to the upper chromosphere, specially if one includes photospheric lines in the IRIS spectral window. Inversion codes should make use of mm-data to improve the accuracy of the inverted chromospheric temperatures from Mg II lines, while revealing more structure otherwise hidden among the inversion noise. This is of great interest to observers working with ALMA since IRIS usually guarantees observing support.

Interestingly, we also found that a combination of a greater number of UV+optical spectral lines allows inferring reasonably accurate (average error ∼30−300 K) temperatures up to the mid-chromosphere (log τ ∼ −5.2) where the the Ca II line cores attain maximum response. This emphasizes the importance of the Ca II H and K observations with the new SST/CHROMIS instrument (Scharmer 2017) whose potential was first shown in Rouppe van der Voort et al. (2017), since these lines not only greatly improved the inversion of temperatures, but also vertical velocities in the low-mid- chromosphere which were not discussed here (see supplementary Fig. A.3). However, this implies that schemes simultaneously fitting a large number of chromospheric lines are not greatly improved by the addition of the 1.2 mm continuum which is sensitive to temperature perturbations at approximately the same heights. The fact that the i6 maps show more structure with higher spatial correlation with respect to the model at log τ ∼ −5.6 is mostly due to the 3 mm continuum whose response peaks around such optical depth.

We find that the response functions to temperature at mm-wavelengths change spatially depending on the underlying temperature/density structure and on the magnetic field to some extent. R1.2 mm spans a larger range of optical depths with two clear peaks at log τ = −4.8 and log τ = −5.2, which is lower down in the atmosphere than the typical maximum response of the core of the Mg II lines (log τ ∼ −5.5). We also found an offset in the optical depth of the maximum response of the h (and k) core and the mm-continuum depending on the underlying magnetic field strength.

This behavior is consistent with the effective formation heights defined as the centroid of the contribution functions of the mm-continuum shown in Loukitcheva et al. (2015) using the same enhanced network model. The authors found that on average the contribution maximum at λ = 1 mm and λ = 3 mm occurs at 900 km and 1500 km respectively, which are slightly higher than the values reported by Wedemeyer-Böhm et al. (2007) who found 730 km and 960 km, respectively, from a quiet-Sun (without magnetic field) simulation. The resemblance between contribution functions and response functions to temperature disturbances happens because temperature is the dominant parameter in mm-continuum formation.

We note that collapsing the entire contribution function to a single value hides the fact that there are in fact multiple heights contributing the the mm-radiation with some overlap between wavelengths. This means that in a quiet-atmosphere ALMA bands provide constraining power not at a single well-defined height in the chromosphere but at a wide range of heights anywhere between just above the temperature minimum and the upper chromosphere depending on the opacity.

In the Bifrost snapshot there are systematic differences between the responses of the Mg II line-cores and the 1.2 mm radiation with the latter being more influenced by perturbations in lower layers of the chromosphere than the former in the majority of the pixels which correspond to quiet atmospheres. On top of strong magnetic elements that trend can reverse, with the line cores responding to perturbations deeper in the atmosphere as consequence of opacity changes. At such locations the mm-continuum responds to temperature perturbations at lower optical depths in tendentially hotter and electron-denser layers of the chromosphere.

In principle, including the 4.5 mm and 8.5 mm bands (not available yet for ALMA observations) could help constraining even higher layers of the atmosphere since Rλ peaks at higher heights at those wavelengths, although we have not performed inversion experiments with those. Similarly, including bands 4 and 5 sampling between 1.4–2.4 mm (not available for solar observations) could also offer additional temperature constraints, while allowing the so-called tomographic or “volume imaging” of the chromosphere in the words of Wedemeyer-Böhm et al. (2007).

5.2. Caveats

We note that the estimated uncertainties are model dependent. The actual accuracy gain in real observations will depend on a number of factors and it is hard to quantify. These include slight temporal/spatial mismatches between different observatories, the signal-to-noise ratio, the absolute flux calibration uncertainties and the target on the solar surface.

In fact, strictly speaking these results refer to a “enhanced-network-like” region of the Sun so it is unclear whether the same results would be found in less quiet regions since we have shown that the response functions to temperature perturbations very much depend on the parameters of the atmosphere. It could be that in active regions and plage the hotter and denser plasma ensures a better coupling of the source functions to temperature, thus permitting better inversions. On the other hand, it might not be possible to ignore synchrotron contribution at mm-wavelengths at certain locations in active regions.

On the technical side we should mention that it is not currently possible to perform simultaneous ALMA observations at two distinct bands as we have assumed, since there are calibration overheads that typical amount to at least 30 min between bands at the moment of writting, although it is not clear whether this could change in the future. Therefore, it might not be possible to benefit from the information of both bands on the same target at the same time. This is definitely a problem for the most dynamic, short-lived phenomena in the chromosphere such as UV-bursts and flares occurring on much smaller time scales of a few minutes. However, in some cases it may be useful to perform some sort of spatial averaging to simply measure the average brightness temperature over larger scales.

In fact, some solar features have enough structural stability over longer periods of time which may justify the usage of such time sparse data. For example, as in the sunspot analysis of Loukitcheva et al. (2017b) who combined ALMA observations with both Band 3 and Band 6 taken in two different days to constrain the average properties of a sunspot’s umbra and penumbra while comparing to model predictions.

Nonetheless, we performed an additional inversion test (i7) with IRIS+ALMA synthetic data, this time using only the continuum around 1.2 mm (Band 6). The inverted maps are shown in the supplementary Fig. A.2. We find that even if we include only this band, we still obtain clear improvements on the inverted temperatures in the low-mid- chromosphere. This test could have been done using only Band 3 with a similar outcome, although they actually respond to temperature perturbations at slightly different layers (see Fig. 10). Band 6 was preferred since it has a higher spatial resolution, currently at 0.63″ for a maximum baseline of 500 m (at the moment of writing), which is closer to the IRIS NUV resolution of approximately 0.33″ , therefore it can more realistically provide reliable results on real data at the current ALMA capabilities. These are expected to improve in the future.

The node placement was extensively investigated and ultimately ruled out with the inclusion of higher number of nodes densely sampling the chromosphere. We investigated this issue in greater detail with the schemes i1 and i2 as shown in Fig. 12. The ECDFs of the absolute temperature residuals (Θ) for different inversion cycles show that the accuracy of the inversions generally increases with the increase of number of nodes, but the bulk of the pixels is not very affected by increasing from ten to fourteen nodes in temperature, so we would not expect to obtain substantially different results with even more nodes since the improvement in the residuals already seems to be damped. For example at log τ = −4.8 the median of the absolute residuals in i1 decreased by 1300 K from six to ten nodes, but it only further decreased another 100 K from ten to fourteen nodes. In i2, from ten to fourteen nodes the median accuracy gain was only 15 K, while the changes in the tails of the distributions are insignificant.

thumbnail Fig. 12.

Temperature residuals versus number of nodes for different inversion schemes. Top row: empirical cumulative distribution functions of the absolute residuals Eq. (3) of the schemes i1 (solid green), i2 (dashed blue) and i7 (dashed-dot red) for different number of temperature nodes at three selected optical depths in the chromosphere; bottom row: corresponding median (circle markers) and the interval (vertical lines) containing 90% of the distribution of absolute residuals.

Open with DEXTER

Figure 12 also shows that inversions with mm-wavelength points are always more accurate than inversions of only Mg II lines at the selected optical depths in the chromosphere. In particular, for a low number of nodes the accuracy gain by simply adding ALMA bands is quite striking (1000 K in the median). Progressive inversion cycles improve both i1 and i2, and eventually with a sufficiently high number of nodes the accuracy of these inversion schemes is of the same order in the layers log τ ∼ [−5.0, −4.0]. However, higher up, for example, at log τ = −5.6, i2 is still more accurate by some 130 K on average, and the temperature maps are a lot smoother because the tail of the distribution of residuals for the overestimated temperatures is significantly reduced (by ∼ 1700 K) as shown in the corresponding bottom panel of Fig. 12.

The ECDFs for i7 show that even the addition of just the continuum around 1.2 mm helps the inversions to reduce the errors across the chromosphere. As expected, in the low chromosphere one does not lose a great degree of accuracy by not simultaneously using the 3 mm continuum since it has much weaker response at log τ > −5.0. For example, at log τ = −4.8 i2 and i7 provide essentially the same result for 6 nodes in temperature which is significantly better than i1, but for 14 nodes there is no significant difference in the accuracy of the three schemes with i2 having only a minor (∼50 K) improvement on the median error compared to i7. However, higher up in the chromosphere, for example, at log τ = −5.6 we do see quite a noticeable difference between i7 and i2 because the 3 mm continuum provides stronger constraints at those optical depths. In this case, increasing the number of nodes from 10 to 14 actually caused a larger error in a portion of the pixels in i7 since, for example, in cycle 2 75% of the pixels had a residual smaller 1100 K whereas in i7 that fraction is 69%, although this could be attributed to the inversion settings.

It is clear that the 3 mm continuum is necessary to recover temperatures at even lower optical depths where the response of the 1.2 mm continuum and UV/optical lines is lower. Therefore we speculate that had we chosen to recompute i4 or i6 with only the 1.2 mm-band, we would have arrived to the same conclusion. This test confirms the intuition provided by the responses functions: the 1 mm continuum is useful to set the gradient between the temperature minimum and the low-chromosphere while improving the accuracy of the inverted temperatures at lower layers around log τ ∼ −5.0, whereas the 3 mm continuum acts in the mid-chromosphere since it responds to temperature perturbations slightly above.

5.3. Conclusions

We have used a snapshot of a 3D radiation-MHD simulation with Bifrost to investigate whether ALMA can help the inversion of commonly used non-LTE spectral lines, while essentially allowing a better representation of the real temperature structure of the solar atmosphere. Inversions of a combination of lines and mm-continuum can retrieve with reasonable accuracy, higher in the photosphere and lower in the chromosphere, the temperature structure of a solar-like atmosphere.

The mm-continua helps the inversions of commonly used non-LTE lines namely Mg II h and k and Ca II H, K and 8542 Å to retrieve a more accurate representation of the temperature stratification in the mid-upper chromosphere where those lines are only weakly coupled to local conditions. The addition of mm-wavelengths is very useful to remove outliers from the error distributions and to substantially reduce the inversion noise, revealing more spatial structure.

The 1.2 mm continuum responds to temperature perturbations typically around log τ ∼ −5.0, and thus is useful to set the gradient between the temperature minimum and the low-chromosphere, whereas the 3 mm continuum has its maximal response slightly above at log τ ∼ −5.6 on average, and extending up to log τ ∼ −6.0, constraining the mid-upper chromosphere.

The response of the mm-continua changes spatially depending on the underlying temperature, density and structure and magnetic field strength to some extent, since we find that on top of kG-magnetic field elements the peak response of the ALMA bands moves to lower optical depths at hotter layers of the chromosphere with higher electron densities than in more quiet areas. This is usually the opposite behavior of the Mg II line cores which have deeper responses at such locations where the mass density in the chromosphere is lower. There are systematic differences in the heights contributing to the response to temperature perturbations between the Mg II lines and the mm-continua. In quiet (weakly magnetic) areas there is a larger gap between the optical depth of the peak response of the h and k cores and the 1.2 mm radiation with the former responding to temperature changes well above the 1.2 continuum, whereas above strong magnetic elements this behavior reverses.

The 1.2 mm band is expected to give the most reliable results at the current ALMA spatial resolution capabilities (0.6″ ), and we have shown that it helps inversions of the Mg II lines and even Ca II 8542 Å alone. However, simultaneous inversions of Mg II and Ca II lines are not as greatly improved by the addition of the 1.2 mm continuum since the latter is sensitive to temperature changes in approximately the same layers as the cores of the H, K and 8542 Å lines in this Bifrost simulation. Provided that enough inversion cycles are performed with sufficiently high number of nodes, together with regularization to improve convergence, inversions of these lines are able to attain the same order of accuracy up to log τ ∼ −5.2, at least in our idealized scenario. This emphasizes the importance of the Ca II H and K observations with the new SST/CHROMIS instrument since these lines not only greatly improved the inversion of temperatures, but also vertical velocities in the low-mid- chromosphere. The 3 mm band is needed in any case to constrain even higher layers of the chromosphere where the lines show weaker responses and their source functions are dominated by scattering.

Depending on specific science goals it may not be possible to run multiple inversions cycles on large patches for different time steps due to the prohibitive computational costs. In that case, inversions with the mm-continuum provide relatively more accurate temperatures for a smaller computational effort since they achieve good a degree of accuracy with less cycles/nodes. The actual accuracy gain in real observations will depend on a number of factors such as slight temporal/spatial mismatches between different observatories, the signal-to-noise ratio, the absolute flux calibration uncertainties and the target itself.

The next decade is likely to see improvements in our understanding of the solar chromosphere now that the community has access to the full electromagnetic spectrum of the Sun from the radio to the ultraviolet in high resolution. This will necessarily require coordinated, multiwavelength campaigns along with inversion codes, such as STiC, to infer the thermodynamical state of the plasma from such observations.


2

Computed using the Scikit-learn library (Pedregosa et al. 2011).

3

The difference between 75th and 25th percentiles.

Acknowledgments

We thank the referee for the useful remarks and suggestions to improve the manuscript. The Institute for Solar Physics is supported by a grant for research infrastructures of national importance from the Swedish Research Council (registration number 2017-00625). The computations were performed on resources provided by the Swedish National Infrastructure for Computing (SNIC) at the High Performance Computing Center North (HPC2N) at UmeÅ University. This research made use of Astropy4, a community-developed core Python package for Astronomy (Astropy Collaboration 2013, 2018). This research has made use of NASA’s Astrophysics Data System Bibliographic Services JdlCR is supported by grants from the Swedish Research Council (2015-03994), the Swedish National Space Board (128/15) and the Swedish Civil Contingencies Agency (MSB). This project has received funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation program (SUNMAG, grant agreement 759548).

References

Appendix A: Supplement figures

thumbnail Fig. A.1.

Temperature residuals of two inversion schemes at selected optical depths for synthetic spectra with null micro-turbulence. In each group of panels, the top row correspond to the fit to synthetic data with constant νturb = 2 km ,s−1 and the bottom row is fit to synthetic data with νturb = 0. The colorbars are in units of 103 K. Rightmost panels: relative error for the two inversion modes as in Fig. 4.

Open with DEXTER

thumbnail Fig. A.2.

Model and inverted temperatures as function of optical depth in a patch of a Bifrost simulation snapshot. Top row: temperature in the model at selected optical depths in the chromosphere; we compare the schemes i1 and i7 for cycle 2 (10 temperature nodes) and cycle 3 (14 temperature nodes). The colorbars are in units of 103 K. r is the correlation coefficient.

Open with DEXTER

thumbnail Fig. A.3.

Model and inverted line-of-sight velocity as function of optical depth in a patch of a Bifrost simulation snapshot. The colorbars are in units of km s−1.

Open with DEXTER

thumbnail Fig. A.4.

Distributions of temperature and electron density at the optical depth of the peak response at two mm-wavelengths. Left panels: compare the temperature histograms along the lines drawn in panels a/b of Fig. 8 for all magnetic field strengths (gray-filled), weak-magnetic elements (|Bmax| < 50 G, dashed orange) and strong magnetic concentrations (|Bmax| > 1000 G, orange-filled). Right panels: analogous comparison in terms of electron densities.

Open with DEXTER

All Tables

Table 1

List of inversion schemes.

All Figures

thumbnail Fig. 1.

Temperature structure and emerging intensities from a Bifrost simulation. Left panel: gas temperature as function of optical depth; central panel: subregion selected for the inversions with temperature scale in units of 103 K; the magnetic field strength at log τ = 0 is represented by thick (2 kG) and thin (0.1 kG) green isocontours. Right panels: emerging intensities at selected wavelengths given in ångströms except where indicated; for each spectral line, the top and bottom images correspond to core and blue wing wavelengths respectively.

Open with DEXTER
In the text
thumbnail Fig. 2.

Inverted temperatures compared to the model. Each column corresponds to a given optical depth and each row to a different inversion scheme (not smoothed) which is to be compared to the model in the top row. The colorbars are in units of 103 K. The three markers in the top row correspond to the spectra in Fig. 6. r is the Pearson’s correlation coefficient.

Open with DEXTER
In the text
thumbnail Fig. 3.

Temperature residuals at selected optical depth nodes. Each panel compares the probability density distributions of the temperature residuals of two inversion schemes (with and without mm-continuum), with the median indicated by the white dotted lines while the 16th and 84th percentiles are shown by dashed lines.

Open with DEXTER
In the text
thumbnail Fig. 4.

Temperature error of the inversions and correlation coefficient as a function of optical depth.

Open with DEXTER
In the text
thumbnail Fig. 5.

Joint probability densities of the temperature in the Bifrost simulation and inverted temperatures at two different optical depths. Each panel compares two inversion schemes (color code as in Fig. 4), and the dashed and solid lines enclose 90% of the distributions. The side plots show the corresponding histograms: top row: i1 vs i2, middle: i3 vs i4, and bottom: i5 vs. i6. The black dotted lines are the TBIF = Tinv locus.

Open with DEXTER
In the text
thumbnail Fig. 6.

Example of fits to the spectra from the Bifrost simulation snapshot. These correspond to the points marked in Fig. 2. Panels ac: overview of Mg II h and k, UV triplet lines; panels di: zoom into the cores of Mg II h and k, UV triplet; panels j, l, n: Ca II K; panels k, m, o: Ca II 8542 Å; panels p, r, t: Fe I 6301, 6302 Å; and panels q, s, u: ∼ 1–3 mm continua. Bottom panels: model temperature (dotted line) versus inverted temperatures.

Open with DEXTER
In the text
thumbnail Fig. 7.

Normalized response functions to temperature perturbations at selected frequencies and pixels in the simulation box. Top row: response functions to temperature of Mg II h and k and Ca II 8542 Å at the same locations given in Fig. 6.; the solid line traces the maximum at each column. Bottom row: normalized responses at 1.2 mm (solid line) and 3 mm (dashed line), and model temperature (dotted line) at those locations.

Open with DEXTER
In the text
thumbnail Fig. 8.

Temperature, response functions to temperature at mm-wavelengths and magnetic field strength in a 2D-cut of a Bifrost simulation. Panel a: temperature (clipped at 15 000 K) as a function of optical depth; the solid and dashed lines are the maximum response in each column at 1.2 mm and 3 mm respectively (smoothed for visibility), and the filled contours enclose the Rλ at half-peak; the peak response the core of Mg II h is shown for comparison (dash-dotted). Panel b: magnetic field strength for context; panels IIV: normalized response functions of the 1.2 mm- (solid) and 3 mm- (dashed) continuum and temperature profile (dotted line) at the selected columns in panels a,b).

Open with DEXTER
In the text
thumbnail Fig. 9.

Mass density, electron-proton number density and response functions to temperature of the mm-continuum and Mg II h. Top: mass-density and electron-proton number density at the magnetic (I) and non-magnetic locations (IV) marked in Fig. 8. Bottom: normalized response functions of the 1.2 mm (solid) and 3 mm-continuum (dashed) and Mg II h core (dashed-dot) at location I (left) and IV (right).

Open with DEXTER
In the text
thumbnail Fig. 10.

Distributions of the maximum response functions to temperature perturbations at selected wavelengths. For the spectral lines we show the distributions at the wavelength where Rλ attains the greatest height in each pixel. The values shown on the right correspond to the median, the 2.5th and the 97.5th percentiles.

Open with DEXTER
In the text
thumbnail Fig. 11.

Empirical cumulative distribution functions of the micro-turbulence velocity residuals for different inversion schemes.

Open with DEXTER
In the text
thumbnail Fig. 12.

Temperature residuals versus number of nodes for different inversion schemes. Top row: empirical cumulative distribution functions of the absolute residuals Eq. (3) of the schemes i1 (solid green), i2 (dashed blue) and i7 (dashed-dot red) for different number of temperature nodes at three selected optical depths in the chromosphere; bottom row: corresponding median (circle markers) and the interval (vertical lines) containing 90% of the distribution of absolute residuals.

Open with DEXTER
In the text
thumbnail Fig. A.1.

Temperature residuals of two inversion schemes at selected optical depths for synthetic spectra with null micro-turbulence. In each group of panels, the top row correspond to the fit to synthetic data with constant νturb = 2 km ,s−1 and the bottom row is fit to synthetic data with νturb = 0. The colorbars are in units of 103 K. Rightmost panels: relative error for the two inversion modes as in Fig. 4.

Open with DEXTER
In the text
thumbnail Fig. A.2.

Model and inverted temperatures as function of optical depth in a patch of a Bifrost simulation snapshot. Top row: temperature in the model at selected optical depths in the chromosphere; we compare the schemes i1 and i7 for cycle 2 (10 temperature nodes) and cycle 3 (14 temperature nodes). The colorbars are in units of 103 K. r is the correlation coefficient.

Open with DEXTER
In the text
thumbnail Fig. A.3.

Model and inverted line-of-sight velocity as function of optical depth in a patch of a Bifrost simulation snapshot. The colorbars are in units of km s−1.

Open with DEXTER
In the text
thumbnail Fig. A.4.

Distributions of temperature and electron density at the optical depth of the peak response at two mm-wavelengths. Left panels: compare the temperature histograms along the lines drawn in panels a/b of Fig. 8 for all magnetic field strengths (gray-filled), weak-magnetic elements (|Bmax| < 50 G, dashed orange) and strong magnetic concentrations (|Bmax| > 1000 G, orange-filled). Right panels: analogous comparison in terms of electron densities.

Open with DEXTER
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.