© ESO 2020

1. Introduction

Strong absorption lines in the visible and infrared wavelengths are important spectral windows into the solar chromosphere. These lines are observable from the ground and contain useful information of chromospheric plasmas. The Hα line of hydrogen has been the most popular of these windows. This line is favored because it is strong, and broad enough for filtergraph observations. The Hα filter images of solar regions display a great variety of intensity structures (Rutten 2008; Leenaarts et al. 2012). Even though the image data of the intensity directly provide much useful (mostly morphological) information of the underlying plasma structures, they do not provide estimates of plasma parameters, which are crucial for understanding the physical conditions. The inference of plasma parameters requires the spectral data of the strong absorption lines and a successful spectral inversion.

Spectral inversion is the process of inferring the plasma parameters from the observed profile of a spectral line. Two types of spectral inversion have been popular in solar observations that assume the constancy of physical parameters. One is the Milne-Eddington inversion, and the other is the cloud model inversion (Beckers 1964). The Milne-Eddington inversion is based on the assumption that the spectral line is formed in a plasma layer of infinite optical thickness where the absorption profile is constant over optical depth and the source function varies with a constant gradient. This inversion has been used mainly to model spectral lines formed in the photosphere and to infer the magnetic fields from their Stokes profiles (Unno 1956; Skumanich & Lites 1987).

The cloud model inversion, on the other hand, assumes that the line is formed in a plasma layer of finite optical depth where the source function as well as the absorption profile is constant over optical depth. This model has been used mostly to infer the physical parameters of cloud-like features lying far above the solar surface (Tziotziou 2007), as was well illustrated in Fig. 1 of Heinzel et al. (1999). A number of variants have been proposed to generalize the original cloud model of Beckers (1964) by incorporating the varying source function (Mein et al. 1996; Heinzel et al. 1999; Tsiropoula et al. 1999), the presence of multiple clouds (Gu et al. 1996), the concept of the embedded cloud (Steinitz et al. 1977; Chae 2014), etc. Despite these variants, the usage of the cloud model inversion is still limited, and is often hampered by the difficulty of choosing the incident intensity profile. Because the incident intensity below the feature of interest cannot be determined from observations, it has to be assumed to be the same as that in its neighborhood, for instance.

Here we present a multilayer inversion for modeling the spectral profiles of strong absorption lines. This represents a combined generalization of the two types of spectral inversion. A strong line is formed over a wide height range of the atmosphere from the photosphere to the chromosphere. The formation of the line in the photosphere can be modeled by the Milne-Eddington model, and the formation in the chromosphere can be modeled by the cloud model inversion. When the source function is allowed to vary with optical depth, there is no fundamental difference between the two types of inversion. Thus we expect that the formation of a strong line can be modeled by the radiative transfer across a finite number of layers in each of which the absorption profile is constant and the source function varies with a constant gradient over optical depth. This is the multilayer inversion we aim to implement.

In this multilayer inversion, parameters other than the source function are kept constant in each layer, but they can vary from layer to layer throughout the solar atmosphere. In this regard, the multilayer inversion method is somewhat similar to the response-function-based inversion proposed by Ruiz Cobo & del Toro Iniesta (1992, 1994) to infer height-varying temperature, magnetic field, and line-of-sight velocity from Stokes profiles.

The multilayer spectral inversion is described in detail in the following section. Its specific version, the three-layer model, is applied to the spectral data of the Hα and Ca II 854.21 nm lines taken by the Fast Imaging Solar Spectrograph (FISS) of the Goode Solar Telescope (GST). We measure wavelengths in nm and pm, and lengths or distances on the Sun in km and Mm: 1 nm = 10⁻⁹ m = 10 Å, 1 pm = 10⁻³ nm = 10⁻² Å = 10 mÅ, 1 km = 10³ m, and 1 Mm = 10⁶ m.

2. Multilayer spectral inversion

2.1. Multilayer model of radiative transfer

We assume that the source function S is independent of wavelength over the spectral line, and is given as a function of height in the atmosphere, which is measured by the optical depth t_λ at a wavelength λ. If S is known as a function of t_λ, the intensity emergent out of the atmosphere I_λ is given by the solution of the radiative transfer equation,

$$\begin{array}{c}\hfill I_{\lambda }={\displaystyle {\int }_0^{\infty }}S(t_{\lambda })exp(-t_{\lambda })\mathrm{d}{t}_{\lambda }.\end{array}$$(1)

The atmosphere consists of the chromospheric layer of finite optical thickness τ₀ at the line center and the photospheric layer of infinite optical thickness, therefore we can rewrite Eq. (1) as

$$\begin{array}{c}\hfill I_{\lambda }=I_{\lambda ,\mathrm{c}}+exp(-{\tau }_{\lambda })I_{\lambda ,\mathrm{p}}\end{array}$$(2)

in terms of the chromospheric contribution I_λ, c and the photospheric contribution I_λ, p defined by

$$\begin{array}{cc}& I_{\lambda ,\mathrm{c}}\equiv {\displaystyle {\int }_0^{{\tau }_{\lambda }}}S(t_{\lambda })exp(-t_{\lambda })\mathrm{d}{t}_{\lambda }=L(0,{\tau }_{\lambda })\hfill \\ \hfill & I_{\lambda ,\mathrm{p}}\equiv {\displaystyle {\int }_{{\tau }_{\lambda }}^{\infty }}S(t_{\lambda })exp(-(t_{\lambda }-{\tau }_{\lambda }))\mathrm{d}{t}_{\lambda }=L({\tau }_{\lambda },\infty ).\hfill \end{array}$$

Here the operator L(t₁, t₂) is defined as

$$\begin{array}{c}\hfill L(t_1,t_2)\equiv {\displaystyle {\int }_{t_1}^{t_2}}S(t)exp(-(t-t_1))\mathrm{d}t.\end{array}$$(3)

We can decompose the chromosphere into N layers of equal optical thickness δ₀ = τ₀/N at the line center (as illustrated in Fig. 1 in the case N = 2). Then we obtain the expression

$$\begin{array}{c}\hfill I_{\lambda ,\mathrm{c}}={\displaystyle \sum _{l=1}^N}exp(-t_{\lambda ,l-1})L(t_{\lambda ,l-1},t_{\lambda ,l}),\end{array}$$(4)

where

$$\begin{array}{c}\hfill t_{\lambda ,l}-t_{\lambda ,l-1}=r_{\lambda ,l}{\delta }_0,\end{array}$$(5)

which assumes constancy of r_λ ≡ χ_λ/χ₀ over optical depth in each layer, where χ_λ is the absorption coefficient at wavelength λ and χ₀ is the absorption coefficient at the central wavelength λ₀ of the line. The optical thickness of the chromospheric layer is then given by

$$\begin{array}{c}\hfill {\tau }_{\lambda }=t_{\lambda ,N}={\displaystyle \sum _{l=1}^N}{r}_{\lambda ,l}{\delta }_0.\end{array}$$(6)

Fig. 1. Three-layer model of radiative transfer that consists of two chromospheric layers and one photospheric layer.

We assume that the source function has a constant gradient dS/dt in each layer. Then the integration yields the expression

$$\begin{array}{c}\hfill L(t_1,t_2)=S_1(1-e^{-\mathrm{\Delta }})+\frac{\mathrm{d}S}{\mathrm{d}t}(1-e^{-\mathrm{\Delta }}-\mathrm{\Delta }{e}^{-\mathrm{\Delta }}),\end{array}$$(7)

written in terms of the source function at top S​₁ and Δ = t₂ − t₁.

In the chromospheric layers, the optical thickness is finite, therefore the gradient can be written in terms of the source function at the top S​_l − 1, and the one S​_l at the bottom of each layer indexed by l, and it follows

$$\begin{array}{cc}\hfill L(t_{\lambda ,l-1},t_{\lambda ,l})=& S_{l-1}(1-e^{-\mathrm{\Delta }})\hfill \\ \hfill & +(S_l-S_{l-1})(1-e^{-\mathrm{\Delta }}-\mathrm{\Delta }{e}^{-\mathrm{\Delta }})/\mathrm{\Delta },\hfill \end{array}$$(8)

with Δ = r_λ, lδ₀.

The absorption profile is assumed to be constant in each layer. In a chromospheric layer, the values of the line-of-sight velocity v_l (defined to be positive in downward motion) and Doppler width w_l are assumed to be constant over optical depth. The damping parameter and continuum absorption are set to zero. The absorption coefficient in the lth layer is described by the Gaussian function

$$\begin{array}{c}\hfill {\chi }_{\lambda ,l}={\chi }_{\mathrm{p},l}exp(-u_{\lambda ,l}^2),\end{array}$$(9)

with the peak absorption coefficient χ_p, l, and

$$\begin{array}{c}\hfill u_{\lambda ,l}\equiv \frac{\lambda -(1+v_{\mathrm{l}}/c){\lambda }_0}{w_l}.\end{array}$$(10)

The absorption profile at wavelength λ₀ is given by

$$\begin{array}{c}\hfill {\chi }_{0,l}={\chi }_{\mathrm{p},l}exp(-u_{0,l}^2),\end{array}$$(11)

with

$$\begin{array}{c}\hfill u_{0,l}\equiv -\frac{v_{\mathrm{l}}}{c}\frac{{\lambda }_0}{w_l},\end{array}$$(12)

and we obtain

$$\begin{array}{c}\hfill r_{\lambda ,l}=exp(-u_{\lambda ,l}^2+u_{0,l}^2).\end{array}$$(13)

In the photosphere, Δ is infinite in Eq. (7), and the gradient can be expressed as

$$\begin{array}{c}\hfill \frac{\mathrm{d}S}{\mathrm{d}{t}_{\lambda }}=\frac{\mathrm{d}S}{\mathrm{d}{t}_{\mathrm{c}}}\frac{\mathrm{d}{t}_{\mathrm{c}}}{\mathrm{d}{t}_{\lambda }}=(S_{\mathrm{p}}-S_N)\frac{\mathrm{d}{t}_{\mathrm{c}}}{\mathrm{d}{t}_{\lambda }},\end{array}$$(14)

where t_c is the continuum optical depth and S​_p is the source function at the level of t_c = 1. Thus Eq. (7) leads to the expression

$$\begin{array}{c}\hfill L({\tau }_{\lambda },\infty )=S_N+(S_{\mathrm{p}}-S_N)\frac{\mathrm{d}{t}_{\mathrm{c}}}{\mathrm{d}{t}_{\lambda }}\end{array}$$(15)

in the photosphere, where S​_N is the source function at the bottom of the lowest layer indexed by N.

The absorption profile is constant over depth in the photosphere as well. Here the continuum absorption is not negligible and the collisional broadening also has to be taken into account. The ratio of line center-to-continuum absorption η as well as the dimensionless damping parameter a is taken to be constant over optical depth. The absorption coefficient in the photosphere is then given by

$$\begin{array}{c}\hfill {\chi }_{\lambda ,\mathrm{p}}={\chi }_{\mathrm{c}}+{\chi }_{\mathrm{p},\mathrm{p}}H(u_{\lambda ,\mathrm{p}},a),\end{array}$$(16)

where χ_p, p is the peak absorption coefficient in the photosphere, and H is the Voigt function normalized to satisfy H(0, a)=1 with

$$\begin{array}{c}\hfill u_{\lambda ,\mathrm{p}}=\frac{\lambda -(1+v_{\mathrm{p}}/c){\lambda }_0}{w_{\mathrm{p}}}·\end{array}$$(17)

Thus we have

$$\begin{array}{c}\hfill \frac{\mathrm{d}{t}_{\mathrm{c}}}{\mathrm{d}{t}_{\lambda }}=\frac{{\chi }_{\mathrm{c}}}{{\chi }_{\lambda ,\mathrm{p}}}=\frac{1}{1+\eta H(u_{\lambda ,\mathrm{p}},a)},\end{array}$$(18)

with the line center-to-continuum opacity ratio η ≡ χ_p, p/χ_c.

In this work, we specifically consider the three-layer model of radiative transfer where the chromosphere is assumed to consist of two layers (N = 2), as illustrated in Fig. 1. We are mainly interested in the chromosphere. Because it covers a wide range of heights over which the physical conditions vary, it is necessary to describe it with a model of at least two layers. The three-layer model is the simplest model that includes the photospheric layer and characterizes both the height variation of the absorption profile and the nonlinear variation of the source function inside the chromosphere.

This model is fully specified by a total of 13 parameters. The optical thickness of the chromosphere is specified by τ₀. The variation in source function is specified by the four parameters S​₀, S​₁, S​₂, and S​_p, and the absorption profile is described by v_p, w_p, a, and η in the photosphere, by v₂ and w₂ in the lower chromosphere, and by v₁ and w₁ in the upper chromosphere.

2.2. Model fitting

We fit each profile of a strong absorption line with the model specified by the independent parameter vector, which is defined as

$$\begin{array}{cc}\hfill \mathit{p}\equiv & (v_{\mathrm{p}},log\eta ,log{w}_{\mathrm{p}},loga,log{S}_{\mathrm{p}},log{S}_2,log{\tau }_0,\hfill \\ \hfill & v_2,v_1,log{w}_2,log{w}_1,log{S}_1,log{S}_0),\hfill \end{array}$$(19)

where each element has a real value and log refers to the common logarithm to the base of 10. The adoption of the logarithmic values at the model parameters automatically guarantees the positivity requirement for η, w_p, a, S​_p, S​₂, τ₀, w₂, w₁, S​₁, and S​₀.

For the regularized model fitting, we employ the technique of constrained least-squares fitting that minimizes the functional

$$\begin{array}{c}\hfill H\equiv {\displaystyle \sum _j}\frac{{(y_j-x_j(\mathit{p}))}^2}{{\sigma }_y^2}+{\displaystyle \sum _i}\frac{{(p_i-p_{i,e})}^2}{{ϵ}_i^2}+{\displaystyle \sum _{\mathit{lm}}}\frac{{(p_l-p_m)}^2}{{ϵ}_l^2+{ϵ}_m^2}·\end{array}$$(20)

Here the first sum is the classical χ² term where y_j is the data, σ_y is the standard noise of data, and x_j is the model with parameters p_i. The second sum is the constraint on the individual components p_i, where p_i, e and ϵ_i are the expectation value of p_i and its standard deviation, respectively, that are to be known from the a priori information. The third sum is the similarity constraint forcing absorption profiles in the two chromospheric layers to be similar to each other as far as the data allow. Specifically, we select (p_l, p_m)=(v₁, v₂), (logw₁, logw₂), and (logS​, logS​₀).

Depending on the value of ϵ_i, the parameter p_i can be categorized as either fixed or free. If ϵ_i is set to be very small, p_i is practically fixed to p_i, e. The value of v_p is inferred from the center of the Ti II 655.958 nm line in the Hα band and the Si I 853.6165 nm line in the Ca II band. It is thus given as an input, and treated as a fixed parameter in the model fitting. In addition, we found after several experiments that most line profiles can be fairly well fit with logη, logw_p, a, and τ₀ being fixed to the values listed in Table 1.

We fix the values of logη, logw_p, and a, because this reduces the degree of freedom and facilitates fitting the portions of the line profile that is formed in the photosphere. We also fix the value of τ₀ because it is not uniquely determined from the fitting. We have examined the performance of the fitting by varying τ₀. As a result, we found that in both lines, the fitting is fairly good, regardless of τ₀, because it varies from 1 to 10. The fit becomes worse when τ₀ becomes larger. Thus we fix τ₀ to 10, the highest value that can yield a sufficiently good fitting.

There remain eight free parameters that have to be determined from the fitting itself. To determine these values, we first apply the fit to a large number of data sets without constraints (by setting ϵ_i = ∞). Then all p obtained with the best fit form an ensemble of p. The mean of p_i over this ensemble is then identified with p_i, e. The value of ϵ_i is set to the standard deviation of p_i in the ensemble or, if necessary, to a higher value. The values of p_i, e and ϵ_i are listed in Table 2.

All the fits are made with respect to the clean wavelengths, that is, the wavelengths where the blending by other lines is negligible. We apply the fitting in two stages. In the first stage, the far wings of the line profile (with wavelength offsets larger than 1.2 Å in the Hα line and 0.8 Å in the Ca II line) are fit by the model

$$\begin{array}{c}\hfill I_{\lambda }=I_{\lambda ,\mathrm{p}},\end{array}$$(21)

assuming τ_λ is negligibly small. This first fitting produces the estimates of two parameters logS​_p and logS​₂, and the construction of I_λ, p over the whole line. With I_λ, p being determined, we can construct the contrast profile,

$$\begin{array}{c}\hfill C_{\lambda }\equiv \frac{I_{\lambda }-I_{\lambda ,\mathrm{p}}}{I_{\lambda ,\mathrm{p}}}=\frac{I_{\lambda ,\mathrm{c}}}{I_{\lambda ,\mathrm{p}}}-(1-exp(-{\tau }_{\lambda })).\end{array}$$(22)

This contrast profile is fit by the corresponding model at the clean wavelengths. This second fitting yields the estimates of the other six parameters: v₂, v₁, logw₂, logw₁, logS​₁, and logS​₀. The special condition of v₂ = v₁, w₁, and S​₂ = S​₁ = S​₀ ≡ S reduces the above equation to that of the classical cloud model.

The goodness of fit is measured by the standard error ϵ defined by the root mean square of the difference between the observed C_λ and the model C_λ, where the average is taken over the clean wavelengths where the fitting is applied.

3. Data and reduction

3.1. Data

We used the spectral data taken with the FISS of the GST at Big Bear Solar Observatory (Chae et al. 2013). The FISS is a dual-band echelle spectrograph that usually records the Hα band and the Ca II 854.2 nm band simultaneously using two cameras. The cameras record the spectral ranges of 0.97 nm (Hα) and 1.29 nm (Ca II) with a pixel sampling of 1.9 pm and 2.5 pm, respectively. With a slit width of 0.16′′, the spectral resolving power (λ/δλ) is estimated at 140 000 and 130 000. We use the PCA-compressed spectral data where noise is suppressed very well (Chae et al. 2013). The signal-to-noise ratio of the compressed intensity data is 700 (Hα) and 250 (Ca II) at the line centers. The height of 40′′ is covered with the sampling of 0.16′′. The two-dimensional imaging of the FISS is achieved by scanning the slit across the field of view with a step size of 0.16′′.

Figure 2 presents the monochromatic images of an observed region constructed from the raster-scan observation of a quiet-Sun region outside active regions. The time of observation was 17:14:06 UT on 2017 June 14. The field of view is 24′′ × 40′′ or 17.4 Mm × 29 Mm. Even though small, the field of view includes both network features and intranetwork features.

Fig. 2. Monochromatic images of a quiet region constructed at ten selected wavelengths. The two symbols mark the position of a intranetwork feature (x) and a network feature (+) selected for the illustration of the observed line profiles and the model fitting.

3.2. Reduction

The average spectral profile of an observed solar region was taken as the reference profile. All the profiles were normalized by the maximum intensity (a proxy of the continuum intensity) of this reference profile. Each profile was then corrected for stray light in two steps following Chae et al. (2013). The observed line profile was corrected for spatial stray light by subtracting 0.027 times the reference profile from it and by then dividing it by 0.973. The corrected line profile was then corrected for spectral stray light by subtracting 0.065 times the maximum value of the line profile from it and by dividing it by 0.935.

The reference profile was also used to calibrate the wavelength precisely. In the Hα band, the centers of the Hα line and the Ti II 655.9580 nm line in the reference profile were determined in pixel units and were used to convert the wavelength pixels into physical wavelength values. In the Ca II band, the Ca II line itself could not be used as a reference because the center of this line averaged over a solar region is known to be offset from its laboratory wavelength. Instead, we used the pair of the Fe I 853.80152 nm line and the Si I 853.6165 nm line for the wavelength calibration. Thus the Doppler shift in each band was measured with respect to the average photosphere of the observed region.

Before applying the model fitting, we corrected each observed line profile for the slightly nonuniform pattern that may exist and may cause the asymmetry between the far blue wing and the far red wing. We fit the intensities at the far red and blue wings > 0.4 nm by a first-order polynomial. This fit, after normalizing by its mean value, was used as the estimate of the nonuniform pattern. The profile corrected for this pattern then had the required symmetry.

The observed spectral profiles of the Hα and the Ca II lines are contaminated by weak spectral lines in the same band, some of which are solar lines and others are terrestrial lines (mostly H₂O lines). Fortunately, the terrestrial lines in our data are very weak, because the observatory is located in a dry land. The model fitting described below uses only the spectral data that are less contaminated by these satellite lines.

4. Results

4.1. Model fitting of line profiles

We applied the multilayer fitting to the profiles of the Hα and Ca II 854.2 nm lines. It took 0.018 s on average to fit each Hα line profile, and 0.046 s to fit each Ca II line profile when we used the Python 3.7 software installed on a laptop computer with a 1.80 GHz CPU and Windows 10 operating system.

We found that the multilayer model fitting is reasonably good in both lines. The fitting of the line profiles of an intranetwork (IN) feature is illustrated in Fig. 3, and that of a network (NW) feature in Fig. 4. The standard error of the fitting is ϵ = 0.007 (IN) and 0.007 (NW) in the Hα line, and 0.013 (IN) and 0.006 (NW) in the Ca II line. The two figures show that the NW feature is different from the IN feature in the spectral characteristics. The NW feature has a broader profile of the Hα line and a higher core intensity of the Ca II line than the IN feature. For this reason, the model fitting produces a larger w₁ of the Hα line and a higher S​₀ of the Ca II line in the NW feature than in the IN feature.

Fig. 3. Three-layer model fitting of the line profiles taken from an intranetwork feature marked by the cross in Fig. 2. The reference profile is the average of all the profiles over the observed region.

Fig. 4. Three-layer model fitting of the line profiles taken from a network feature marked by a plus in Fig. 2.

4.2. Spatial distribution of the model parameters

We find from Fig. 5 that the parameter maps of the Hα line are very regular. There are few noticeable irregularities. The map of ϵ indicates that the fitting is better than ϵ = 0.01 in most of the spatial pixels in the field of view. This shows the effectivity of our spectral inversion. The model parameters can be used for a reasonable derivation of chromospheric plasma parameters.

Fig. 5. Maps of the Hα line parameters. The spatial averages have been subtracted in the Doppler velocity maps.

The maps of the source function provide information on temperature or radiation field at the different atmospheric levels. The map of S​_p is very similar to the near-continuum intensity image at 0.4 nm off the line (see Fig. 2) and corresponds to the spatial distribution of the temperature at the continuum-forming level of the photosphere. The map of S​₀ is practically a reproduction of the core intensity image of the line (Fig. 2). This is expected because it is well known that the core of a very strong line is formed at the outer part of the atmosphere, and its intensity is expected to be very close to the source function at the formation level. The maps of S​₂ and S​₁ show the structure of the source function at the top level of the photosphere and at the lower level of the chromosphere, respectively. The map of S​₂ displays the inverse convective pattern consisting of dark cells and bright lanes, and the map of S​₁ is similar to that of S​₀.

Our spectral inversion produces the Doppler velocity maps at the three atmospheric levels (Fig. 5). These are very useful probes of the atmospheric dynamics. The map of v_p shows the velocity pattern associated with granulation and photospheric oscillations, while the maps of v₂ and v₁ mostly show the velocity pattern linked to jet-like features and chromospheric oscillations. The similarity of v₂ and v₁ (with the Pearson correlation of 0.56) suggests that the difference in the oscillation phase between the lower chromosphere and the upper chromosphere is not large, implying that the wavelength of the associated waves may be longer than the height difference between the two layers.

The maps of either w₂ or w₁ of the Hα line (Fig. 5) provide a convenient way of distinguishing between the network regions and the internetwork regions. The network regions have high values of Doppler width and the internetwork regions have low values. The map of w₁ practically corresponds to the temperature map in the upper atmosphere because the hydrogen atom is light and its thermal speed dominates the Doppler broadening.

Figure 6 shows the maps of the Ca II model parameters. Because the maps of S​_p, S​₂, and S​₀ are very regular, they contain information on the spatial distribution of the physical parameters at the corresponding atmospheric level. The map of S​₀ is very similar to the core intensity image (Fig. 2). It is similar to the map of the Hαw₁ in Fig. 5. On the other hand, the map of S​₁ is not regular in that it contains a number of discontinuities at different spatial scales.

Fig. 6. Maps of the Ca II line parameters. The spatial averages have been subtracted in the Doppler velocity maps.

The maps of v_p, v₂, and v₁ of the Ca II band look very similar to those of the Hα line. It is expected that the Ca II band v_p and the Hα band v_p have very similar patterns (with the Pearson correlation of 0.88) because they were derived from two weak lines that formed in the photosphere. The similarities of the Ca IIv₂ and the Hαv₂ (correlation = 0.72) and that of the Ca IIv₁ and the Hαv₁ (correlation = 0.76) are more interesting. This supports the notion that the formation layer of the Ca II 854.2 nm line significantly overlaps that of the Hα line.

The maps of the Ca IIw₂ and w₁ are more complicated than those of the Hαw₂ and w₁. There is a tendency for w₂ and w₁ to be larger in the NW regions than in the IN regions. This tendency, however, is not as strong as in the Hα line. We find that numerous small patches of enhanced value of w₁ are found to be scattered in the IN regions as well as in the NW regions.

We can estimate the random fitting errors by assuming that the physical conditions and systematic errors vary very smoothly with position so that their second-order spatial derivatives are approximately zero. Then we assume that the nonzero second-order finite difference of three neighboring points is contributed by the random errors. We applied this second-order derivative method to the variation in parameters along the slit direction at each slit position, which provides the estimates of the random error in each parameter. By repeating this process at a number of slit positions and by taking the average, we obtained the representative random errors of the mode parameters as listed in Table 3.

4.3. Inferring temperature and nonthermal speed

It is possible to separately determine the hydrogen temperature T_H and the nonthermal speed (or microturbulence speed) ξ of the upper chromosphere by combining the values of w₁ of the two lines. By denoting w₁ of the Hα line by w_H and that of the Ca II line by w_Ca, we can derive the expressions

$$\begin{array}{cc}& T_{\mathrm{H}}=8100\mathrm{K}{\left[\frac{w_{\mathrm{H}}}{0.025\mathrm{nm}}\right]}^2(1-0.59{\left[\frac{w_{\mathrm{Ca}}}{w_{\mathrm{H}}}\right]}^2)\hfill \end{array}$$(23)

$$\begin{array}{cc}& \xi =5.40\mathrm{km}{\mathrm{s}}^{-1}\frac{w_{\mathrm{Ca}}}{0.015\mathrm{nm}}{(1-0.042{\left[\frac{w_{\mathrm{H}}}{w_{\mathrm{Ca}}}\right]}^2)}^{1/2},\hfill \end{array}$$(24)

which indicate that T_H is mostly determined by w_H because of the light mass of the hydrogen atom and ξ, mostly by w_Ca because of the high mass of a Ca II ion.

In the IN feature of Fig. 3, we have w_H = 0.025 nm and w_Ca = 0.015 nm, which leads to the estimates T = 6500 K and ξ = 5.1 km s⁻¹. In the NW feature of Fig. 4, we have w_H = 0.036 nm and w_Ca = 0.023 nm in the Ca II line, which leads to the estimates T = 130 00 K and ξ = 7.8 km s⁻¹.

Figure 7 shows the spatial distribution of T_H and ξ determined from the Hα line w₁ in Fig. 5 and the Ca II line w₁ in Fig. 6. The map of T_H is clearly very similar to the Hα line w₁, and that of ξ, to the Ca II line w₁. Making use of the second derivative method described above, we estimated the errors of T_H and ξ at 140 K and 0.1 km s⁻¹.

Fig. 7. Top: spatial distributions of TH and ξ. Bottom: number distributions of TH and ξ.

The number distribution of T_H is approximately a double Gaussian. The peak around 7000 K represents the typical temperature of the upper chromosphere in the IN regions and the other peak at 11 000 K, that in the NW regions. In contrast, the number distribution of ξ is singly peaked at 6.1 km s⁻¹, which represents the typical ξ of the upper chromosphere, regardless of the specific region.

We confirm from Fig. 8 that S​₀ of the Ca II line is strongly correlated with T_H. The Pearson correlation coefficient is 0.75. This means that the value of S​₀, or equivalently, the core intensity, is sensitive to temperature. This result supports the study of Cauzzi et al. (2009). The temperature sensitivity of the core intensity of Ca II line originates from the property of the Ca II line. In this line, the collisional excitation by electrons significantly contributes to the Ca II line emission. In contrast, the Hα line S​₀ or the core intensity is not correlated with T_H at all, in agreement with Leenaarts et al. (2012), which means that the collisional excitation by electrons is less important in the Hα line emission. The Hα core intensity is correlated with the average formation height, with the lower intensity corresponding to the higher average formation height (Leenaarts et al. 2012).

Fig. 8. Scatter plot of TH vs. the Ca II line S​0.

4.4. Temporal variations

It is interesting to determine the variation in model parameters with time at a fixed point. Figure 9 presents the line profiles of each line that were taken from an IN point and stacked over time (approximately one hour). The time epoch was set to the instant of the highest downflow speed. We find from the figure that the observed data (λ − t intensity maps) are fairly well reproduced by the corresponding models. The mean value of the fitting error ϵ was estimated at 0.0062 in the Hα line and at 0.014 in the Ca II line. The residual intensity maps show some systematic patterns that are correlated with the intensity distribution. Systematic errors seem to contribute significantly to the fitting errors. Some vertical patterns represent telluric lines and weak photospheric lines. These patterns are outside the clean wavelengths and were excluded in the inversion.

Fig. 9. Wavelength-time maps of the Hα and Ca II intensity stacked over time – observation (left) and model (middle), and maps of residual intensity normalized by the continuum intensity (right). The grey scale of the normalized residual intensity spans from −0.03 to 0.03.

Figure 10 clearly indicates that all the model parameters as well as T_H and ξ fluctuate with time. The standard deviations of the fluctuations are 2.4 km s⁻¹ in the Hα line v₁, 2.0 km s⁻¹ in the Ca II line v₁, 1100 K in T_H, and 0.9 km s⁻¹ in ξ. These are far larger than the estimates of the random errors 0.06 km s⁻¹, 0.04 km s⁻¹, 140 K, and 0.1 km s⁻¹, respectively. Thus the fluctuations shown in Fig. 10 may represent the real variations in the physical conditions. The variations of v₁ of the two lines represent chromospheric three-minute oscillations. The temporal variations of T_H and ξ may be mostly attributed to the three-minute oscillations as well. The proper interpretation of the variations of T_H and ξ in terms of chromospheric oscillations, however, is not straightforward, and requires further investigation to understand how the three-minute oscillations affect the formation of the two lines and their model parameters of the multiline spectral inversion. The three-minute oscillations in the upper chromosphere are now understood as highly nonlinear waves of very long wavelength that have sharp discontinuities, for example, shock fronts (Chae & Litvinenko 2017; Chae et al. 2018).

Fig. 10. Temporal variations of some of the model parameters as well as TH and ξ.

5. Discussion

We proposed a multilayer spectral inversion in order to infer the physical condition of chromospheric plasma from strong absorption lines, and we specifically investigated the three-layer model in detail. This model is fully specified by 13 parameters. By fixing some parameters and determining one parameter from a weak line that is formed in the photosphere, we reduced the number of free parameters to 8 and determined these free parameters by applying the constrained nonlinear least-squares fitting. We found that the three-layer model reproduces most of the observed profiles of the Hα line and the Ca II 854.2 nm line taken from a quiet-Sun region fairly well. The random fitting errors are much smaller than the intrinsic spatial variations, therefore the maps of most model parameters look very regular. Thus, we conclude that our implementation of the spectral inversion is successful in determining the model parameters, and we have reached the first goal of successfully implementing the three-layer spectral inversion.

The ultimate goal of our investigation is to infer the physical conditions from the determined model parameters. This goal is distinct from the first goal and requires further investigation. The physical condition at each time is to be specified by the height variations in temperature, electron density, velocity, and nonthermal speed, but the determined model parameters are only a limited number of height-averaged source function, Doppler width, and line-of-sight velocity in each layer. Moreover, the spectral inversion does not provide any information on the height range of each layer. Hence, adopting some assumptions is indispensable to derive the physical conditions from the model parameters. The physical plausibility of the derived physical condition depends on the validity of the adopted assumptions.

For instance, hydrogen temperature T_H and nonthermal speed ξ of the upper chromosphere were derived from the model parameters: w_H and w_Ca, the Doppler widths of the Hα and the Ca II 854.2 nm lines based on the assumption that the upper chromosphere seen through the Ca II line is the same as the upper chromosphere seen through the Hα line. Under this assumption, the same values of T_H and ξ contribute to both w_H and w_Ca. In fact, T_H is mostly determined by w_H because a hydrogen atom is lighter and has a higher thermal speed than any other atoms or ions, as was previously noted by Cauzzi et al. (2009) and Leenaarts et al. (2012).

We determined T_H and ξ in all the pixels of the observed quiet-Sun region and found that the distribution of the hydrogen temperature peaks at around 7000 K in the intranetwork regions, and around 11 000 K in network regions. The mean value of nonthermal speed is found to be 6 km s⁻¹ regardless of intranetwork regions and network regions. Our measurements may be compared with the estimates of Cauzzi et al. (2009). They determined T_H and ξ in a quiet region from the core widths of the two lines. Each core width was defined very like an FWHM and was directly measured from the line profile without taking the effect of radiative transfer into account. After subtracting the “intrinsic contribution” from the line widths that probably represents the opacity contribution to the line width, they obtained T_H ranging from 5000 K to 60 000 K, and ξ ranging from 1 km s⁻¹ to 11 km s⁻¹. These ranges include and are broader than the corresponding ranges we obtained. Cauzzi et al. (2009) also noted that the core width of the Hα line is larger in the network patches, supporting the notion that network regions are heated more strongly than intranetwork ones.

The assumption used for the derivation of T_H and ξ is reasonable, but it also has certain limitations. When we determine the volume-averaged values of T_H and ξ in the upper chromosphere, the assumption seems to be good enough, as described above. If we were to determine the values in a fine Hα structure such as a fibril, the assumption is not satisfactory. The Hα fibrils are often invisible in the Ca II 854.2 nm line, which means that these plasma structures are transparent in the Ca II line. This may partly explain the lack of regional dependence on ξ we found above. ξ is mostly determined by the Doppler width of the Ca II line, but fibrils are transparent in the Ca II line, whereas they are clearly visible in the Hα line. If fibrils have higher temperature and higher nonthermal speed than the low-lying layers where the Ca II line is formed, it is likely that ξ in network regions was underestimated, and T_H was overestimated.

Making better use of the determined model parameters requires a good understanding of the line formation. In this regard, one can learn much from the forward modelling of the non-LTE radiative transfer. The forward modelling is different from the spectral inversion in that it initially adopts the height variation of temperature, density, and velocity, and calculates the source function by solving the rate equations for level populations and can determine the formation height. The core intensity of a line is close to the source function in the outermost part of the line formation region. It is well known that the Hα source function in the outer layers is mostly determined by the radiation field and is not sensitive to local temperature, but sensitive to height (Leenaarts et al. 2012). In contrast, the Ca II line source function in the outer layers is still affected by the collisional excitation and is sensitive to the local electron temperature. This explains the strong correlation shown in Fig. 8 and the similarly strong correlation between the Hα core width and the Ca II 854.2 nm core intensity reported by Cauzzi et al. (2009).

Finally, we would like to mention that the effect of isotopic splitting may have to be investigated in the multilayer spectral inversion of the Ca II 854.2 line in the future. Leenaarts et al. (2014) investigated this effect on the bisector and inversions. They showed that the line core asymmetry and inverse C-shape of the bisector of the Ca II 854.2 nm line can be explained by the isotopic splitting: the larger difference in line-of-sight velocity difference of more than 2 km s⁻¹ can result from this.

We conclude that the multilayer spectral inversion successfully infers the model parameters from the observed profiles of strong absorption lines. The model parameters can be used to derive the physical parameters of chromospheric plasma such as the temperature, when physically plausible assumptions are made and if the line formation is well understood. Determining subtle variations of physical parameters in space or in time requires further careful investigations. Combining the model parameters of several lines would be of much help in determining the height variation of physical parameters. We expect that the multiline multilayer spectral inversion will serve as a powerful tool to infer the physical parameters of chromospheric plasma from observations.

Acknowledgments

JC greatly appreciates Juhyung Kang’s assistance in the implementation of the method using the Python software. This research was supported by the National Research Foundation of the Korea (NRF-2019H1D3A2A01099143, NRF-2020R1A2C2004616), and by the Korea Astronomy and Space Science Institute under the R&D program(Project No. 2020-1-850-07) supervised by the Ministry of Science and ICT.

References

All Tables

All Figures

	Fig. 1. Three-layer model of radiative transfer that consists of two chromospheric layers and one photospheric layer.
In the text

	Fig. 2. Monochromatic images of a quiet region constructed at ten selected wavelengths. The two symbols mark the position of a intranetwork feature (x) and a network feature (+) selected for the illustration of the observed line profiles and the model fitting.
In the text

	Fig. 3. Three-layer model fitting of the line profiles taken from an intranetwork feature marked by the cross in Fig. 2. The reference profile is the average of all the profiles over the observed region.
In the text

	Fig. 4. Three-layer model fitting of the line profiles taken from a network feature marked by a plus in Fig. 2.
In the text

	Fig. 5. Maps of the Hα line parameters. The spatial averages have been subtracted in the Doppler velocity maps.
In the text

	Fig. 6. Maps of the Ca II line parameters. The spatial averages have been subtracted in the Doppler velocity maps.
In the text

	Fig. 7. Top: spatial distributions of TH and ξ. Bottom: number distributions of TH and ξ.
In the text

	Fig. 8. Scatter plot of TH vs. the Ca II line S​0.
In the text

	Fig. 9. Wavelength-time maps of the Hα and Ca II intensity stacked over time – observation (left) and model (middle), and maps of residual intensity normalized by the continuum intensity (right). The grey scale of the normalized residual intensity spans from −0.03 to 0.03.
In the text

	Fig. 10. Temporal variations of some of the model parameters as well as TH and ξ.
In the text