Evaluating Non-LTE Spectral Inversions with ALMA and IBIS

We present observations of a solar plage in the millimeter-continuum with the ALMA and in the Ca 8542 and Na 5896 spectral lines with the Interferometric BIdimensional Spectrometer (IBIS). Our goal is to compare the measurement of local gas temperatures provided by ALMA with the temperature diagnostics provided by non-LTE inversions using STIC. In performing these inversions, we find that using column mass as the reference height scale, rather than optical depth, provides more reliable atmospheric profiles above the temperature minimum and that the treatment of non- LTE hydrogen ionization brings the inferred chromospheric temperatures into better agreement with the ALMA measurements. The Band 3 brightness temperatures are higher but well correlated with the inversion-derived temperatures at the height of formation of the Ca 8542 line core. The Band 6 temperatures instead do not show good correlations with the temperatures at any specific layer in the inverted atmospheres. We then performed inversions that included the millimeter continuum intensities as an additional constraint. Incorporating Band 3 generally resulted in atmospheres showing a strong temperature rise in the upper atmosphere, while including Band 6 led to significant regions of anomalously low temperatures at chromospheric heights. This is consistent with the idea that the Band 6 emission can come from a range of heights. The poor constraints on the chromospheric electron density with existing inversion codes introduces difficulties in determining the height(s) of formation of the millimeter continuum as well as uncertainties in the temperatures derived from the spectral lines.


INTRODUCTION
In recent years, solar physics has gained a new tool for the study of the chromosphere-high-resolution observations in the millimeter wavelengths from the Atacama Large Millimeter/submillimeter Array (ALMA, Wootten & Thompson 2009). The solar radiation in the wavelength regime of 1 to 10 mm is believed to mostly arise from the continuum emission in the chromosphere; see Wedemeyer et al. (2016) for discussion. Previous observations in this regime were performed with smaller arrays, and thus had limited resolution; one such example is the Berkeley-Illinois-Maryland Array (BIMA) (Loukitcheva et al. 2009), which has a resolution of 10 arcseconds, sufficient to reveal the shape of magnetic networks, but incapable of resolving the finer chromospheric structures that are apparent at shorter wavelengths (Wedemeyer et al. 2016). With ALMA, spatial resolutions of 1 arcsecond or better are now achievable, fine enough to resolve distinct features in the solar atmosphere on similar spatial scales as seen in imaging spectrographs operating at optical and ultraviolet (or other) wavelengths.
This makes it possible to now directly measure the electron gas temperature in the mid-to upper-chromosphere at smaller spatial scales. Wedemeyer-Böhm et al. (2007) demonstrated that the contribution functions for several of ALMA's bands peak in distinct regions of the chromosphere, with increasing wavelengths sampling increasing heights. In a similar analysis, Loukitcheva et al. (2015) synthesized ALMA observations from a 3D radiation-MHD simulation of a quiet-sun network, finding that millimeter observations at a resolution of 1 arcsecond or better can provide a reasonable measurement of the electron temperature in the solar atmosphere. Thus, millimeter observations can provide a direct measurement of the temperatures at various, if ill-determined, heights in the solar chromosphere, though care should be taken to account for all relevant physics and observational biases, especially the effects of line-of-sight integration and limited spatial resolution on observed millimeter intensities, to ensure correct inferences (Martinez-Sykora et al. 2020;Eklund et al. 2021).
Temperatures in the chromosphere have traditionally been obtained indirectly through the inversion of spectral lines (del Toro Iniesta & Ruiz Cobo 2016). "Inversion" here refers to any mathematical problem that uses the value of an integral to determine the expression being integrated. In this case, the expression is the physical conditions within the solar atmosphere, and the effects thereof on radiation propagating outward, e.g. absorption, emission, scattering, polarization; the integral itself is a summation of these effects along the line-of-sight; and the value of that integral is the emergent intensity at a given wavelength, which is what is observed. These inversion techniques have already been applied with great success in the photosphere, with spectral lines such as the Fe I 6302Å line and inversion codes assuming local thermodynamic equilibrium (LTE) such as SIR (Ruiz Cobo & del Toro Iniesta 1992), and recent years have seen the development and application of new codes such as STiC (de la Cruz Rodríguez et al. 2019) and DeSIRe (Ruiz Cobo et al. 2022), which are capable of treating multiple lines in non-local thermodynamic equilibrium (non-LTE), particularly the Mg II h and k lines, the Ca II H and K lines, and the Ca II infrared triplet.
While spectral inversions have proved quite useful in inferring the properties of the solar atmosphere, they can be somewhat inconsistent in certain regimes: especially in non-LTE, inversions are an ill-constrained problem with several possible degeneracies and pitfalls. One way to remedy this problem is to include mm observations as an additional observational constraint in the inversions. This was first demonstrated by da Silva Santos et al. (2018), who used a slice from a 3D r-MHD simulation to compute synthetic observations, including both chromospheric spectral lines and millimeter continuum intensities, invert them with STiC in various combinations, and compare the results, finding that inclusion of millimeter diagnostics should greatly improve the reliability of inversions in recovering temperatures in the low-to mid-chromosphere. More recently, da Silva Santos et al. (2020) applied this technique to real observations of a plage in the Mg II h and k lines and the Band 6 continuum, finding that the inclusion of the millimeter continuum reveals the presence of pockets of hot and cool gas in the low chromosphere that are invisible to the UV lines.
In this work, we explore some of these questions by performing the first non-LTE inversions of spectral lines in the visible and near-infrared whose heights of formation significantly overlap with the formation region of the millimeter continua sampled by ALMA. We examine the preferred height scale for inversions including the solar chromosphere and the importance of including non-LTE hydrogen ionization into account in that process. We compare the temperatures derived from this inversion process with the brightness temperatures derived from the simultaneous high-resolution observations at 1.2 and 3.0 mm. We also perform tests where we include the ALMA continuum intensities directly in the inversion process to see if this additional information helps derive representative temperature profiles in the atmosphere. We examine the effects of the possible extended formation ranges of the millimeter continuum, especially for the Band 6 measurements, and the effect of uncertainties in the electron density in the inversion outputs.

OBSERVATIONS
We use observations obtained as part of coordinated program between ALMA and the Dunn Solar Telescope (DST, Dunn & Smartt 1991) on 2017 April 23. At the DST, the Interferometric Bidimensional Spectrometer (IBIS, an imaging spectrograph Cavallini 2006;Reardon & Cavallini 2008) observed the target region in several visible and near-IR spectral lines. In addition, the IRIS (De Pontieu et al. 2014) and Hinode (Kosugi et al. 2007) satellites were co-pointing during the time of these observations, though these datasets were not included in our analysis. Context images and magnetic field maps were available from SDO/AIA (Lemen et al. 2012) and SDO/HMI (Pesnell et al. 2012).
The observed target, shown in Figure 1, was the leading portion of NOAA active region 12651, which was stable and showed no significant flaring activity. This same region, with observations take one day prior, was the one analyzed by (da Silva Santos et al. 2020). The center of the field of view (located at E04, N11, or µ = 0.96, where µ is the cosine  Middle: HMI LOS magnetogram, scaled (linearly) between -500 to 500 G.; Right: IRIS SJI image at 2796Å, averaged over one minute. The field of view of ALMA Band 6 is shown as the blue circle and the ALMA Band 3 field of view is shown as the red circle. The rectangle indicates the region of interest for which inversions were performed.
of the angle between the surface normal and the line of sight) could be classified as magnetic plage, with some quieter areas present in the southern portion of the IBIS and ALMA Band 3 fields of view (FOV). Even though the target was close to an active region, the average unsigned flux in the ALMA field of view was only about 20 G, as measured from the HMI magnetograms. The leading spot of the active region was located about 50 east of the target center, outside the area covered by IBIS and ALMA, but the penumbra of the spot did slightly extend into the edges of the IBIS field of view.

DST/IBIS Observations
We observed the target with IBIS from 15:13 to 19:06 UT in conditions of good to excellent seeing. The observations were taken in repeated scans of several different spectral lines (including Ca II 8542Å Na I D1 5896Å, and H I 6563 A), in slightly different combinations or cadences. In this work, we use data taken during two series of scans acquired between 15:54-16:35 UT (hereinafter referred to as series 1554) and 17:25-18:12 UT (series 1725). The second series has been studied previously in Molnar et al. (2019Molnar et al. ( , 2021. The cadence for the scans through all the observed lines was 10 and 16 seconds in the two series, respectively. We perform our spectral inversions on a selected scan obtained during each series, in the time intervals 15:57:10-15:57:18 and 17:36:32-17:36:46. We display the time-average of line-minimum-intensity images taken from 15:54-16:01 UT in the left column of Figure 2. We applied linearity, dark correction, flat-fielding, and fringe-removal corrections to the IBIS data. In order to correct for optical and atmospheric image distortions, we used the nearest-in-time HMI continuum intensity images as a reference in order to precisely map first the IBIS whitelight-channel images, and then the co-temporal spectral data onto a regular, fixed spatial grid (with an image scale of 0.096 /pixel), effectively removing the bulk of the seeing distortions. 1

ALMA Observations and Processing
The ALMA data were obtained with the array in configuration C40-3, with a maximum baseline of 396 m (rather than the nominal 460 m). The data analyzed here were snapshots obtained in Band 6 (1.2-1.3 mm, 229-249 GHz, 0.8 beam size) and Band 3 (2.8-3.3 mm, 92-108 GHz, 1.9 beam size) at 15:57 and 17:36 UT, respectively. Even after the application of the standard interferometric calibration procedures using CASA (Petry & CASA Development Team 2012), the ALMA images are still heavily influenced by the phase disturbances due to water vapor variations in Earth's atmosphere, which causes phase errors across the array and small-scale distortions on the images (Shimojo   To counteract this effect, we apply the self-calibration technique, where an average (clean) image over some temporal interval is used to derive and remove the random phase variations in the individual images within that observing interval (as described in more detail in Chai et al. 2022). This allows us to minimize the atmospheric image distortions in a self-consistent manner that should not introduce significant artifacts.
The effective field of view of ALMA Band 3 (Band 6) images is only ∼ 60 (30 ) in diameter (Figure 2, right), which is smaller than the large-scale solar structures. To provide information on the more extended or background emission from the Sun, several single-dish observations that scan the full disk were taken nearly simultaneously with the interferometer array (White et al. 2017), albeit with a significantly lower spatial resolution. Through the feathering process provided in CASA, we combined the two datasets covering the full range of spatial scales. The full-disk image was normalized such that the brightness temperature in the central region of the disk (diameter 200 arcsec) in Bands 3 and 6 was set to 7300 K and 5900 K, as determined by White et al. (2017), providing an absolute temperature accuracy of about 2-5%.

INVERSIONS
The inversions were performed using the STockholm inversion Code (STiC, see de la Cruz Rodríguez et al. 2019), a regularized non-LTE inversion code based on the radiative transfer code RH, a multi-line, non-LTE spectral synthesis code (Uitenbroek 2001). The inversion code takes in the calibrated data, I obs (λ), and seeks a combination of plasma parameters, p, that reproduce the observed spectra; in this case, the fitted parameters are temperature, line-of-sight (LOS) velocity, and microturbulent velocity, at a number of discrete depths in the model atmosphere. For each iteration, the synthetic spectra, I syn , are compared to I obs using the merit function: Here, the first term is the chi-squared coming from the likelihood of the data given the model, and the second term is the sum of regularization terms for each model parameter, n is the number of degrees of freedom and is equal to (N − N p ), and σ i is the measurement uncertainty for intensity at λ i . Note that σ i can also be tuned to give certain wavelengths more weight than others; in our case, the line centers are given more weight in order to ensure a good fit at chromospheric heights. If the fit is poor, i.e. χ 2 is large, then the atmospheric parameters are adjusted according to the response of the emergent intensity, and the process repeats until convergence, upon which the code returns the best-fit atmospheric parameters and the corresponding synthetic spectra.
We inverted IBIS spectral scans of the Ca II and Na I lines together, in a similar manner to the multi-line inversions of Ca II and Na I performed by da Silva Santos et al. (2022), which required additional calibrations to the data. The spatially dependent instrumental blueshift was removed by interpolating and resampling onto a uniform wavelength grid, which is necessary in order for the code to properly take into account the broadening profile of the instrument. The wavelengths were weighted such that only the observed points contribute to the merit function. To allow accurate inference of LOS velocities, the mean spectra were aligned with the FTS ATLAS (Wallace et al. 1993(Wallace et al. , 1998 to provide an absolute wavelength calibration. Finally, the intensities were normalized to the disk center continuum.
An important aspect of the inversion process, especially in the non-LTE case, is the choice of the atomic models that incorporate radiative transitions other than the ones that are directly observed. The iterative treatment of the populations of all these atomic levels makes non-LTE inversions numerically demanding and requires a balance between accuracy and simplicity. We found that a 4-level sodium atom and a 6-level calcium atom provide a good compromise between computational speed and realism.

Height Scale Comparison
Inversion codes such as STiC seek a height-dependent parameterization of the solar atmosphere, which results in a model with a large number of free parameters; namely, if the atmospheric model parameterizes p physical properties over a range of N height steps, then the model has N × p free parameters. As a result, spectroscopic inversions are an ill-posed numerical problem with a high degree of degeneracy. To make the problem more tractable, inversion codes parameterize each atmospheric parameter in terms of a smaller number of nodes, specific heights in the atmosphere for which the numerical solutions are computed, and between which the atmospheric parameters are smoothly interpolated. In order to accurately infer the temperature stratifications, it would nominally be ideal to use a relatively large number of nodes to recover as much vertical detail as possible, assuming sufficient computational resources. Using many nodes, however, tends to lead to oscillatory solutions and poor convergence, which STiC circumvents through regularization, which penalizes oscillatory or divergent solutions. The spectroscopic inversion processes described above does not recover the atmosphere along a fixed geometrical height scale. Historically, inversions have always been performed using optical depth as the height scale, which has proved to work well for the photosphere where it gives a rather uniform vertical resolution. However, optical depth scale greatly compresses the chromospheric structures: while the solar atmosphere up to the transition region spans roughly six orders of magnitude in both τ and ξ, the chromospheric temperature rise spans two orders in ξ and less than one order in τ (de la Cruz Rodríguez et al. 2019, Figure 3); this can be a significant obstacle for high-fidelity inversions of chromospheric lines, especially when inverting for multiple chromospheric diagnostics simultaneously. The STiC code introduced the possibility of using column mass (ξ) as the height scale. This both simplifies the solution of the hydrostatic equilibrium and helps us to better resolve the chromospheric structures. For inversions seeking to probe the structures and gradients of the chromosphere, column mass provides a more natural height scale, as the  Column mass -> optical depth   inversion-derived temperature [K] Inversions vs. Ca Core Intensity Figure 5. Maps of the ratios between the logarithm of the minimum intensity in the Ca II 8542Å line core and inversionderived temperatures, using the two different height scales, at heights of log ξ = −3.9 (left) and log τ = −4.9 (middle). Displayed values are scaled to ±20% of the mean ratio for each case. Right: scatter plots of inversion-derived temperatures, using the ξ (red) or τ (blue) height scales, versus the Ca II 8542Å line core minimum intensity, plotted on a logarithmic scale. solar atmosphere is approximately hydrostatically stratified, and therefore, column mass is much more closely related to geometric height.
The nodes can be placed arbitrarily for either of these two height scales. To capture the chromospheric temperature gradients in the optical depth scale, however, would require using a large number of nodes in the regions above the (local) temperature minimum. Since the location and slope of the chromospheric temperature rise can vary significantly between profiles, it is difficult to have a optimal node placement in optical depth that consistently and accurately fits lines such as the calcium infrared triplet, which are sensitive to the location of this gradient. Because the chromosphere is not as compressed in column mass scale as in optical depth scale, we found that an equidistant or uniform node placement in column mass provided a much more robust inference of the atmospheric profiles.
We performed tests comparing inversions in optical depth and column mass, fitting the observed calcium and sodium lines together with simultaneous Band 3 continuum over the full plage region, in order to recover the temperature (T ), line-of-sight velocity (v los ), and microturbulent velocity (v turb ). For both the log ξ and log τ inversions, we determined [T, v los , v turb ] using [21,11,4] equidistant nodes covering the range of heights (−5.1 < log ξ < 0.9 for column mass, −6.0 < log ξ < 0.0 for optical depth) from the lower photosphere to the upper chromosphere. The regularization weights chosen such that the penalty term in the merit function is of a similar magnitude as the χ 2 term. The quality of the fits ( Figure 3) for the inversions in log τ was significantly worse overall than for log ξ: while the final χ 2 was overall similar for both (median χ 2 for log τ inversion was 1.35 and for log ξ inversion 1.44), the log τ inversions had ∼ 10 − 100× more anomalous pixels that did not converge to χ 2 ≤ 10, and the inferred chromospheric temperatures tended to be less smoothly stratified in height, especially in the central plage regions. This can be seen by comparing the histogram of the temperature profiles recovered from the inversions for the two height scales, as shown in Figure 4. To further highlight the differences, we apply a conversion between these profiles into their corresponding alternative height scale based on the assumption of hydrostatic equilibrium, as shown in the upper left and lower right panels. For the column mass inversions, we used the per-pixel optical depth scales as computed by STiC, while for the optical depth inversions, we obtained the column mass from the total gas pressure: ξ = p gas /g, where g is the acceleration due to gravity. This shows how the steep temperature rise into the chromosphere, which takes place smoothly between −4 ≤ log ξ ≤ −2, is compressed into a much smaller interval of −4.6 ≤ log τ ≤ −4. This makes the inversion in optical depth particularly sensitive to the placement of the nodes such that they properly sample this steep change (in this height scale) in atmospheric parameters.
Further, the recovered temperature profiles from the optical depth inversions showed significant oscillatory behavior in the upper photosphere, most notably in the interval −3 ≤ log τ ≤ −1). We also note that the inversions in optical depth tended to show a falling temperature in the chromosphere above the nominal optical-depth-unity height for the Ca II 8542Å line core (log τ ∼ −4.8), with some cases of anomalously low temperatures below 5000 K, while in the column mass inversions, the temperature profiles continued to rise above those heights. In addition, the column mass . Response functions of the Ca II 8542Å and Na I 5896Å lines, computed for the FALC model in log ξ (top) and log τ500 (bottom). The blue lines are the synthesized spectral profiles for each method. Note that the heights of the responses calculated in this standard model atmosphere are expected to differ from the heights inferred from the models derived through the inversion process.
inversions produce temperature structures that are better correlated with the Ca II 8542Å line core intensity ( Figure  5), with outliers primarily consisting of fibril structures and localized cold spots not apparent in the line core minimum image. The minimum intensity of this line is known to be sensitive to the temperatures of the middle chromosphere (see e.g. Cauzzi et al. 2009;Quintero Noda et al. 2016), so we can use it as a reasonable proxy to test our inversion results at the typical chromospheric heights of formation of the line core.
To further illustrate the differences between these two alternative height scales, we computed the response functions in both height scales for the FALC model (Fontenla et al. 1999), following the RF calculation methods from Quintero Noda et al. (2016). These response functions, shown in Figure 6, clearly show several key features discussed above: the photosphere is better sampled in optical depth, whereas the chromosphere has a greater extent in column mass; the sodium core is likewise more sensitive to temperatures over a broader range of optical depths; and the calcium line core response, which in log τ is localized to a range of ∼ 0.2 dex, extends in log ξ over almost a full scale height. Furthermore, in contrast to the RFs computed by Quintero Noda et al. (2016), the calcium line shows almost no sensitivity to temperature over the ∼ 1 dex range of the temperature minimum, corresponding roughly to −4.5 ≤ log τ ≤ −3.5 and −2 ≤ log ξ ≤ −1. Given our interest in studying the thermodynamic conditions in the temperature minimum and chromosphere, we chose to proceed with the analysis of the inversions done in column mass scale, as we consider them more reliable in this regime.

Electron Density Model
It may be essential to take non-LTE hydrogen ionization effects into account in order to properly model electron densities when attempting to accurately determine chromospheric temperatures. In the tenuous environment of the chromosphere, radiation is able to escape upwards, leaving fewer photons available to ionize hydrogen and thereby reducing the electron density for a given temperature. Therefore, when considering the non-LTE hydrogen ionization, the local electron density in the chromosphere is lower. This then reduces the coupling between the Ca II source function and the local temperature, which would then produce a lower intensity in the line core. To properly reproduce the observed line core intensities, the inversion process must therefore assume higher chromospheric temperatures. The authors of STiC modified the RH synthesis module such that it can treat hydrogen as an active non-LTE species and iteratively solve for the electron pressure (albeit with a significant increase in the computational time required for each inverted pixel).
We used the column mass inversions above, both with and without non-LTE hydrogen ionization, to quantify this effect. Figures 7 and 8 illustrate this effect, with non-LTE hydrogen ionization increasing the inferred chromospheric temperatures by ∼ 1000 K. While we do not expect the calcium and ALMA temperatures to be identical, the inferred temperatures from the calcium line should at least be close to those measured from ALMA Band 6, given their similar  expected heights of formation. Inclusion of these physics primarily alters electron densities in the chromosphere, though significant departures from LTE electron densities can be observed as far down as the upper photosphere and temperature minimum, around log ξ = −2 (log τ = −4.1). Unfortunately, including this correction to the electron densities introduces significant numerical instability, which manifests as both salt-and-pepper "inversion noise" and occasional failures of individual profiles to converge. We mitigated this issue by inverting with two cycles and strong regularization, and the atmospheric profiles of any pixels that failed to converge (∼ 1%) were replaced with the median of the surrounding pixels.

Final inversion configurations
The final configuration for all our inversions combining the Ca II and Na I lines used [21,11,4] nodes in [T, v los , v turb ], equally spaced in column mass between −5.1 ≤ log ξ ≤ 0.9, with regularization such that the penalty term in the merit function is typically of order unity. We used spectral scans taken simultaneously with the ALMA data at 15:57 (series 1554) and 17:36 UT (series 1725). For the series 1554 data, a simultaneous white-light continuum image is included to compensate for the poor sampling in the wings of the Na I line. The ALMA data, when included, are weighted similarly to the spectral data.
We performed the inversions of these datasets two times, once using only the IBIS spectral scans, and once including the cotemporal ALMA brightness temperature maps as inversion inputs, Band 6 in the case of the series 1554 scans and Band 3 with the series 1725 scans.
In addition to the observed ALMA data, we also took the opportunity to build on the work of Molnar et al. (2019), specifically in the nearly linear relationship between the H-α line width and the Band 3 brightness temperature. As all the data taken with IBIS also included H-α line scans, we were able to use a snapshot of the line width as a pseudo-band 3 image for the series 1554 data. The line-width snapshot was blurred to the resolution of ALMA Band 3 and converted to an approximate 3-mm brightness temperature using the relation in their work. We tested the effects of including this proxy data in the inversions, both with and without the real observed Band 6 data, in Section 4.4.

COMPARISONS BETWEEN SPECTROSCOPIC INVERSIONS AND ALMA
One primary aim of this work is to study the agreement between the temperatures inferred from spectroscopic inversions and the brightness temperatures observed with ALMA. With two largely independent measures of temperature, each with a different mechanism of formation, we have the means to not only test the fidelity of the temperatures recovered in the inversion process, but also to test our understanding of the formation of the millimeter continuum. This is particularly interesting because we are using the Na I D and Ca II 8542Å lines, whose heights of formation overlap significantly with the regions of the atmosphere contributing to the millimeter emission in the observed bands. Band 3 vs. Ca Inversion can also address the question of whether including observations in Band 3 or Band 6 (or both) in the inversion process might provide more realistic constraints when using these spectral lines.

ALMA observations and inversion-derived temperatures
We begin by comparing the ALMA Band 3 image, taken at 17:36 UT, with the inversion-derived temperatures from the cotemporal spectral line observations. Figure 9 shows the comparison of the inversion-derived temperature at log ξ = −3.8 (roughly corresponding to log τ = −5.0) and the ALMA band 3 brightness temperature. This specific height was chosen because the inverted temperatures show the best correlation with the core intensity (depth) of the Ca II 8542Å line profile, indicating this is the column mass layer corresponding, on average, to the peak of the response function (e.g. Quintero Noda et al. 2016;da Silva Santos et al. 2018). While the maps are similar on large spatial scales, with bright/hot plage in the center of the field and hotter fibrils extending outward into the internetwork region below, they differ significantly on smaller scales. Given the lower spatial resolution of the ALMA Band 3 data (approximately 2 ), it is not clear whether the small scale structures seen in the Ca-II-derived maps would similarly be present in the Band 3 continuum if higher resolution could be achieved.
Furthermore, it can be seen that the ALMA temperatures are systematically higher than the temperatures inferred from inversions. This is not unexpected, as previous works (e.g. Loukitcheva et al. 2015;Wedemeyer et al. 2016) show that the Band 3 (3.0 mm) radiation is typically expected to form at slightly higher chromospheric layers than the Ca II core, in a region where we expect the transition region temperature rise to have begun. We note that the difference in absolute temperatures is significantly larger than the expected errors in the absolute calibration (White et al. 2017) or the expected offsets between the different sub-bands (Wedemeyer et al. 2020). Nonetheless, the large-scale structures are fairly similar, with Pearson correlation coefficient R = 0.69.
We find that the correlation coefficient between the inverted temperatures and ALMA brightness temperature maps falls off rapidly at greater heights, for log ξ < −3.8, dropping nearly linearly from the peak to 0.2 at the upper bound of our inverted atmospheres (log ξ = −5.1). This indicates that caution should be taken when interpreting the inverted atmospheres as they become unreliable above the peak of the response function of the highest-forming line. This is expected, and here we confirm it by using the independent temperature measurement from ALMA.
In a similar manner, we compare the ALMA Band 6 image, taken at 15:57 UT, with the inversion-derived temperatures. Figure 10 shows the comparison between inverted temperatures at log ξ = −3.8 and the ALMA band 6 brightness temperatures. In this case, the temperatures for these two diagnostics are much closer to each other in terms of mean value (7000 K and 6900 K, respectively), but the correlation between them is significantly weaker (only 0.22), as seen in the large cloud of points in the scatterplot on the right side of the figure. There are no other values of ξ with a significantly improved correlation with the inversion temperature. While our average Band 6 temperature is similar to the value of 7000 K found by Jafarzadeh et al. (2019) for plage regions surrounding a sunspot, those authors found a much higher correlation coefficient with the radiation temperature of the upper chromospheric C II line (0.83). Similarly, Bastian et al. (2017) found temperatures of ∼7300 K for their definition of the plage areas in this same dataset, and a good correlation with the Mg II h line radiation temperature (0.71). However, we note that the correlations in these two studies were based on just a restricted portion of the field of view which tended to select regions of higher Band 6 temperatures, whereas our correlation is calculated over our full field and spanning a broader range of temperatures. A similar result was was found by Chintzoglou et al. (2021) who noted the sensitivity of these measured correlations to the definition of the plage mask. As discussed in Section 4.3, the regions of high Band 6 temperatures are indeed well reflected in the inversion-derived temperatures maps for the upper atmosphere, but the correspondence breaks down when the full span of observed temperatures is considered. The Band 6 contribution functions can vary over a significant range of physical heights, based on the local temperatures, resulting in a lack of overall correspondence between the brightness temperatures and the temperature structures at a single height.
Previous radiative transfer calculations using MHD simulations have shown that the heights of formation of the Band 6 (1.2 mm) emission may spread, even for the same column, from the temperature minimum region to heights above that of the Ca II 8542Å line core (Loukitcheva et al. 2015;da Silva Santos et al. 2018). However, recent calculations by (Martinez-Sykora et al. 2020), including the effects of nonequilibrium and time-dependent ionization of hydrogen and helium, found that the Band 6 radiation should form at a narrower range of heights very similar to that of the Band 3 radiation. Yet, in our case, while the ALMA temperature maps in Figure 9 and Figure 10 show some similar large-scale structuring following the plage magnetic fields, the details on smaller scales seems to be quite different and not entirely consistent with the idea that they are formed at nearly coincident heights in the atmosphere. Some of these differences may be due to the evolution of the chromospheric fine structure in the 100 minutes between the times when these two observations were acquired. However, we examined the magnitude of the changes seen in either band over the full one-hour duration of their respective observing blocks, and we find that the magnitude of the evolutionary changes over these timescales are smaller than the differences seen between the two bands. The distinct differences between the inversion correlations for these two bands in our dataset indicate there may be more significant variations in the formation of these two continua then expected.
The ALMA Band 6 temperature map shows some very cold pixels ( 3000 Kelvin), especially at the top of the field of view, that do not exist in the spectroscopic inversions (STiC enforces 3100 Kelvin as the lowest possible temperature because chemical equilibrium solvers cannot converge at lower temperatures). These cold regions are not present in the maps of inversion-derived temperatures at chromospheric heights, suggesting that the classical spectral lines in the visible and near-infrared might not be highly sensitive to such cold pockets in the atmosphere. This is supported by de la Cruz Rodríguez et al. (2012), who found a similar result for inversions of synthetic Ca II data, attributing this insensitivity to 3-D scattering effects. Indications of similar cold spots have been identified using observations of the CO 4.7 µm line (Stauffer et al. 2022), though it remains unclear at which heights the CO is forming.

Simultaneous inversions of spectroscopic data and ALMA Band 3 data
To test whether the inclusion of the millimeter continuum as an additional input in our inversion process leads to more reliable atmospheric profiles, we inverted the series 1725 spectroscopic data together with cotemporal ALMA  Figure 11. Effects of including Band 3 T b in the inversion of the series 1725 data. Top: inversion-derived temperatures using only IBIS spectroscopic data. Bottom: inversion-derived temperatures using both IBIS and ALMA Band 3 data. Temperature maps are shown for two heights, located at 0.6 dex above (log ξ = −4.4, center) and below (log ξ = −3.2 right) the Ca II 8542 A line core.
Band 3 data. The Band 3 continuum is expected to be formed in the upper chromosphere, therefore it is expected to act mostly as a constraint on the temperatures in the upper chromosphere, with less significant effects on the temperatures in lower regions to which the Ca II 8542Å line is sensitive. This was explored by (da Silva Santos et al. 2018), who compared inversions of synthetic Ca II 8542Å (and Fe I 6302Å) profiles with and without the ALMA Band 6 and Band 3 continua (their inversion setup i3 and i4). They found that without the Band 3 continuum, the derived temperatures at regions above the height of formation of the Ca II 8542Å line are poorly constrained and unrealistic. We find similar behavior in our inversions, as seen in Figure 11, where the central column shows the inversion-derived temperatures at log ξ = −4.6 (roughly corresponding to the log τ = −5.6 shown in da Silva Santos et al. (2018)). When only the spectral lines are inverted (upper row), the temperature profiles in the higher layers do not recover any significant temperature increase in the upper chromosphere. However, when the ALMA Band 3 temperatures are included, the inferred temperatures show the expected increase of temperatures in the upper atmosphere, ranging from approximately 7000-12000 K. As expected, the changes in the temperatures at lower heights, around the height of formation of the Ca II 8542Å line core, are less dramatic, but it does appear to increase the dynamic range of temperatures in these regions. The right column of Figure 11 shows the temperatures at log ξ = −3.2 (roughly log τ = −4.4) when inversions were done without and with the Band 3 brightness temperature. When the millimeter continuum is included, the average temperatures at this height appears to be lowered, from T ∼ 6500 K to 5800 K. Some of the smaller-scale structuring of the higher temperatures in the spectral-line-only inversions appear to be replicated when Band 3 is included. These inversions also introduce regions of significantly lower temperature at this height, that are not seen in the inversions without the Band 3 (e.g. x=-54, y=261), sometimes corresponding to regions of higher temperatures in the Band 3 temperatures. There is also a region centered on x=-51, y=240, where the lower temperatures seen at log ξ = −4.4 correspond to an increased temperature at log ξ = −3.2. Some of the changes seen at this lower height may be unrealistic or artifacts of the inversion technique, raising the question of whether the inclusion of the millimeter continuum observations may have mixed (height-dependent) effects on the reliability of the recovered atmospheres. A similar effect was seen in da Silva Santos et al. (2018) with synthetic inversions including both the Band 6 and Band 3 continua together with the Ca II 8542Å line, but we see evidence of similar issues even with only the higher-forming Band 3.

Simultaneous inversions of spectroscopic data and ALMA Band 6 data
We next explore the effects of including ALMA Band 6 observations in the inversions (instead of Band 3), using the series 1554 IBIS spectral scans together with cotemporal ALMA Band 6 data. This shorter wavelength is expected to typically form lower in the chromosphere, overlapping to a greater extent with the formation region for the Ca II 8542Å line and perhaps Na I 5896Å as well (Loukitcheva et al. 2015), though Chintzoglou et al. (2021) point out that in some dynamic phenomena the Band 6 emission might come from much higher in the atmosphere. Therefore, we expect this additional diagnostic to provide stronger, or more consistent, constraints for the inversion process.
Examining the results shown in Figure 12, we see that the most significant differences in the recovered temperatures are at the greater heights, log ξ = −4.4 (central column). In many cases, these inversions produce seemingly unphysical atmospheres, introducing extended regions of extremely cold (∼ 4000 K) gas in the upper-chromosphere (log ξ > −4.0), above the warmer (∼ 8000 K) gas at log ξ = −3.8 needed to match the Ca II 8542Å line core. This occurs both in cases where the Band 6 temperatures are low (such as the extended cool regions at the top or at the left of field of view) and where the temperatures are close to the average (for example, the extended log ξ = −4.4 cold region at x=-65 and 245 < y < 255). We note that the temperatures observed in Band 3, expected to be formed roughly at this height, doesn't show such low values, nor are these low temperatures consistent with those recovered from the inversions shown in Sec. 4.2. The inversion-derived temperatures at log ξ = −3.2 shown in the right column of Figure  12 show less significant differences compared to the spectral-line-only inversions. As with the inversions including Band 3, there does appear to be an puzzling inverse correlation in some locations, with lower temperatures at this height in correspondence to higher temperatures in the upper layers (e.g. x=-85, y=263 or -86 < x < 81 and 243 < y < 250). We note that STiC has a hard lower limit on the allowable temperatures of 3100 K. Even if the inversion code inferred that the cold Band 6 temperatures were coming from lower in the atmosphere, it might be difficult for it to properly fit atmosphere to some of the very low Band 6 temperatures (T < 4000 K). Looking back at Figure 4, we see that, even from spectroscopic-only inversions, a significant amount of pixels reach this hard-set minimum temperature around column masses log ξ = -1.0 to -2.0; however, this is more likely due to a lack of temperature sensitivity at these heights, combined with the tendency of the regularization to prefer smooth temperature profiles.
More problematically, there is no single atmospheric layer (or value of ξ) where the inverted temperatures on aggregate correlate well with the ALMA band 6 brightness temperatures. This can be seen in Figure 13, where we have partitioned the Band 6 brightness temperature map into three regions based on the observed Band 6 brightness temperatures. We define "hot" (T b > 7900 K), "tepid" (5900 < T b ≤ 7900 K), and "cold" (T b ≤ 5900 K) masks. We then applied those masks to the cubes of inversion-derived temperatures and for each case identified the column mass layer that had the best correlations with the Band 6 temperatures encompassed in that same mask. We found that the hottest regions seen in Band 6 (T b 7900 K) correlated best with temperatures at log ξ = −4.4, which falls slightly above the Ca II 8542Å line core (Pearson correlation coefficient, R = 0.4). The coldest regions in the brightness temperature maps instead seem to form lower down in the atmosphere, at log ξ = −2.9, just above the temperature minimum (R = 0.7). This indication of multi-height formation can be seen in the close visual similarity of the structures seen under the hot and cold masks in both the observed Band 6 brightness and the inverted temperatures. Interestingly, for the "tepid" regions in the middle column, while the best correlation is at log ξ = −3.9, there is less similarity between the structures in the two temperature maps (R = 0.1). This lack of correlation with a specific height in the atmosphere might indicate that the Band 6 T b for these intermediate temperatures are not generally representing a single or well-constrained region of the solar atmosphere, but might be mixing contributions from multiple heights. One way to determine from which layers the mm-continuum originates is to analyze the response functions (RFs) for the mm-continuum ( Figure 14). To do this, we sorted the pixels in order of increasing T b , then selected 100 pixels at regular intervals of 40 K, spanning from 2500-6500 K (labeled cool) and 7500-11500 K (labeled hot). Given the general level of uncertainty in the Band 6 temperature values, this gave a randomly distributed selection of points in these two different temperature regimes. For each pixel selected, we extracted the 1-D atmospheric models from the inversions of the series 1554 data, both with and without the inclusion of Band 6, and computed the RFs for both 1.2 and 3.0 mm millimeter continuum. The resultant RFs for both millimeter bands are shown in Figure 14, along with their respective temperature stratifications, with the columns comparing the two inversion schemes (with and without Band 6).
The most noticeable difference is the significant increase of the peak of the response functions (log ξ < −4.5) in the hot pixels for both the 1.2 and 3.0 mm continuum when the observed Band 6 temperatures are included in the inversion process (rightmost column). Without Band 6 as an inversion input, only a sampling of the 3.0 mm RF's and few of the Band 6 RF's extend to heights of log ξ < −4.0. The temperature profiles in the bottom row show that the inversions of the Ca II 8542Å line only capture the expected strong rise in temperatures at the top of the chromosphere if also provided with additional information about the presence of high temperatures as observed in the ALMA Band 6. It is interesting that these IBIS and Band 6 inversions show, for hot pixels in the range 7500-9500 K, a generally extended hot chromosphere from −2.8 > log ξ > −5.0, with the Band 6 RF's often extending over a large part of this range. Instead for the pixels with Band 6 temperatures above 9500 K the most common solution from the inversions indicates an initial temperature rise around log ξ = −3.0, followed by a temperature drop at log ξ ∼ −3.8 (the approximate height of formation of the Ca II 8542Å line), and then a very steep temperature rise above log ξ = −4.2. The Band 6 RFs in this case are very localized in height to this region of a steep temperature gradient. While this might represent real differences in the chromospheric conditions, the significant bifurcation in the class of temperature profiles might indicate that instead these are simply different solutions from the inversion process when trying reconcile the emissions from the spectral lines and millimeter continuum.
For the cool pixels in the left two columns, the biggest change in the RF's occur for the lowest temperature pixels (T b < 3500K). For these cases, when including Band 6 in the inversion process, the peak of the RF's are moved very low in the atmosphere, to log ξ = −1.5, the very heart of the temperature minimum in the atmosphere. While the inversions without Band 6 as an input reach similarly low temperatures in the middle atmosphere, when constrained by the presence of very low Band 6 T b , solutions must be found that keep the upper atmosphere relatively cool and compressed in order to minimize the millimeter continuum opacity in the upper atmosphere. Otherwise, the 1.2 mm radiation would be formed higher up in the atmosphere, where the Ca II 8542Å line places limits on the allowable minimum temperatures. Even with these atmospheric modifications, the inversions for these coldest pixels do not fully reproduce the extreme low temperatures seen in the Band 6 observations. For many of the other pixels, the RF's show similar average heights log ξ ∼ −3.5 for both the 1.2 and 3.0 mm radiation. However, the RF's are more uniform and sharply peaked for the inversions that include the Band 6 inversions, primarily because the solutions in this case again show a temperature decrease in the chromosphere above the height of formation of the Ca II 8542Å line (albeit at slightly lower column mass, greater height than for the hot pixels). We also note that for the IBIS and Band 6 inversions, the vertical extent of the photosphere also appears to increase, with the onset of temperatures below ∼ 5000 K occurring at log ξ ∼ −1.3 instead of log ξ ∼ −0.7 in the spectral line only inversions. As for the hot pixels, it is not clear whether the inclusion of the Band 6 brightness temperature together with the inversion of the Ca II 8542 A line leads to more realistic inferred temperature profiles.
These issues likely stem from the limitations of a hydrostatic, time-independent model atmosphere, in which the ionization fraction is calculated from the temperature and gas pressure. This model is unable to account for phenomena such as rapid adiabatic cooling, in which the cooling timescale is significantly shorter than the recombination timescale, which can lead to pockets of cold gas with an overabundance of free electrons. The opposite problem, in which the hotter gas of the mid-to upper-chromosphere is assumed to contribute too many free electrons, may also contribute to the poor fits. The proper distribution of free electrons along the atmospheric column is critical for determining the height sensitivity of the 1.2 mm radiation. This may lead to extended response functions, in which case the observed radiation temperatures become a weighted average across a range of heights that may extend from the temperature minimum to the upper-chromosphere. This temperature averaging, and its sensitivity to the time-dependent electron density, may make it critically challenging for inversion processes to properly interpret the observed T b . This may not be possible with current non-LTE inversion codes, and may require either an independent measurement of the electron density at various heights, or more likely, an inversion method that does not rely on time-independent hydrostatic equilibrium.
It should be noted that these shortcomings in our inversions are not necessarily due solely to the ALMA Band 6 observations. Recently, da Silva Santos et al. (2020) performed spectroscopic inversions on IRIS observations of a plage in the Mg II h and k lines, together with simultaneous observations in ALMA Band 6. They found that these observations synergize well together, with the Band 6 continuum effectively constraining temperatures in the low-to mid-chromosphere, slotting in neatly between the wings and core of the Mg II lines. Also noteworthy is that these joint IRIS and ALMA inversions obtained some extremely cold chromospheres in certain regions; however, in that case, the cold chromospheres were located in quiet regions corresponding to extremely low Band 6 temperatures, whereas some of our joint inversions find cold chromospheric regions even in locations with moderate Band 6 temperatures in the middle of the plage (Figure 12). It may be that inversions combining spectral lines that form at overlapping heights to the millimeter radiation, and whose information might conflict under the simplified assumptions of the inversion code, lead to solutions that are not representative of the true, dynamic solar chromosphere. Further work is needed to understand the effective source regions of the Band 6 radiation and how to properly integrate it with other spectral lines.

Simultaneous inversions with both ALMA bands
This dataset provides us the opportunity to investigate the value of employing simultaneous observations of both Band 3 and Band 6 wavelengths in the inversions with the goal of better constraining the chromospheric temperature profile. Obtaining such a dataset is not possible with the currently available observing modes with ALMA, but has been discussed as a possible development of the facility capabilities, in essence dedicating subsets of the array antennas  to different wavelength bands. In order to simulate such as dataset, we make use of the linear correlation between the ALMA Band 3 brightness temperatures and the Hα line widths (Molnar et al. 2019). We take the H-alpha width map obtained simultaneously with the actual Band 6 observation, smooth it with a 2 FWHM Gaussian smearing to match the ALMA resolution and apply their linear fit parameters to generate a corresponding Band 3 temperature map, shown in Figure 15. The simulated Band 3 image looks roughly similar to that observed 90 minutes later during series 1725 (compare to Figure 9), with changes consistent with the typical evolution on these time scales.
The two temperature maps show similarities in some of the spatial structuring, most notably in the regions where the brightness temperatures are either high (T 1.2mm b > 8000K) and low (T 1.2mm b < 5000K). In the regions with the intermediate temperatures, the correspondence between the two maps is much less pronounced. This might indicate that in hotter regions, the height of formation of the 1.2 mm continuum increases to more closely match that of the 3 mm wavelengths, and that conversely in colder regions the height of formation of the 1.2 mm interval might move lower in the atmosphere. The differences in the appearance in the intermediate regime might indicate that the heights of formation of the continuum at these two wavelengths might not always track as closely as suggested by (Martinez-Sykora et al. 2020).
In order to provide a direct comparison with the inversions incorporating the Band 6 data, using the same spectral line scans and solar structures, we used the proxy Band 3 temperature maps as an input into the inversion of the series 1554 data. The recovered temperatures at log ξ = -3.8, corresponding most closely to the height of formation of the core of the Ca II 8542Å line, is shown in the lower-left panel of Figure 16. Similarly to the results in Section 4.2, including the Band 3 radiation alone doesn't significantly change the temperatures in the mid-chromosphere, which remain determined primarily by the information from the Ca II 8542Å line.
The comparison for the temperatures recovered at this height with the spectral inversions incorporating only Band 6 is shown in the upper-right of Figure 16. In this case the method infers regions of very cold gas both in places with low Band 6 brightness temperature (e.g. at the top of the field), but also in the hotter central regions above the magnetic network. Both of these inversion-derived temperatures can be compared to the intensity in the core of the Ca II 8542Å line (upper-left) or to the temperatures derived at this same height when carrying out the inversions only the spectral lines (left panel of Figure 9).
Finally, we performed the same inversions using both the observed Band 6 and proxy Band 3 images as additional inputs into the inversions, with the derived temperatures at this same height shown in the bottom right of Figure 16. Compared to the inversions with just one of the ALMA bands, these joint inversions achieved overall worse fits and display similar issues over a wider physical region. In particular, the regions of extremely cold gas noted in section 4.3 now extend across much of the plage, while for the regions with particularly low Band 6 temperatures, many of the inversions were unable to reproduce the observed Band 6 temperatures.
To understand the nature of the solutions found with the different combinations, we examine the temperature profiles for a slice through the data cube (whose locations is indicated by the white dashed line in Figure 16. The temperatures inferred in the case of using only the spectral lines are shown in the upper left of Figure 17, while the other three panels show combinations that include one or the other ALMA bands, or both (lower right). The spectral-line-only inversions show a relatively well behaved rise in temperatures from T ∼ 4000 K in the upper photosphere to T ≥ 7000 K in the chromosphere. However, only with the inclusion of the Band 3 temperatures do we see a strong rise at the top of the inversion boundary, with T ≥ 11000 K at log ξ ≤ -5. Otherwise, the inclusion of the Band 3 information doesn't produce significant effects in the photosphere or low chromosphere.
However, when using only the Band 6 as an additional input into the spectral inversions, certain recovered atmospheric profiles become problematic, with very low temperatures T < 4000 K found in the mid chromosphere. It also changes the inferred conditions lower down, with high temperatures extending down into the nominal temperature minimum (251 < y < 255) or moving the temperature minimum region higher(265 < y < 269). This does produce some atmospheres that have a strong temperature rise at logξ = -5, but many some display surprisingly low temperatures at the top of the boundary.
When including both Band 3 and Band 6 as an input, the temperature rise at the top of the boundary is more uniform, but the regions of low temperatures in the middle chromosphere (logξ = -4) are even more extended (as was seen in the spatial maps in Figure 16. As a reminder, this portion of the slice lies above a region of relatively enhanced magnetic field and bright Ca II 8542Å line core intensity. Overall, these different combinations of the inversions reinforce what we had shown earlier. The addition of (proxy) Band 3 to the inversion of the spectral lines primarily enforces solutions showing a strong temperature rise above the height of formation of the core of the Ca II 8542Å line, as is physically reasonable. The inversions including Band 6 can be very unreliable in certain pixels, in a way that is not easy to predict. These problems are only exacerbated when including both ALMA wavelengths at the same time. This appears to create solutions, with pockets of very cold plasma in the middle chromosphere, that seem to be unrealistic and don't correspond to the signal from Ca II 8542Å expected to form in this same region.

DISCUSSION AND CONCLUSIONS
We have presented observations of a solar plage region observed simultaneously in two visible/near-IR spectral lines (Ca II 8542Å and Na I 5896Å) and in the millimeter continuum. We performed spectroscopic inversions on these data to recover the temperature profiles through the photosphere and chromosphere at each spatial pixel. For a typical inversion setup, when using a fixed-number of regularly space nodes, we found that using column mass for the height scale in the inversion, instead of optical depth, was more reliable in inferring the chromospheric physical parameters. This is because the column mass scale naturally spreads the chromospheric temperature rise over a greater range of nodes in this height scale, which can more easily and stably be reproduced in the inversion process. Even though the column mass scale tends compress the photosphere into a smaller range of column-mass values (−2 < log ξ < 1 versus −4 < log τ < 0), we found it was still reliable in fitting the wings of the lines and recovering the conditions in the lower atmosphere. In fact, our inversion-derived temperature profiles in the photosphere using the optical depth scale showed small-scale fluctuations with height that were not present in the column-mass inversions. Similarly, a dip in the temperature profiles in the upper chromosphere (at the limit of the Ca II 8542Å sensitivity was present in the optical-depth inversions, while in most cases the column-mass inversions showed a nearly monotonic temperature rise at those heights. We then compared the temperatures derived from the joint inversion of the Ca II 8542 and Na I 5896Å lines (using the column-mass height scale) with the simultaneous brightness temperatures measured by ALMA at 1.25 and 3.0 mm. These two, independent temperature measurements show some correlation, as expected. The Band 3 images exhibit similar large-scale structures as the temperature maps from the approximate height of formation of the Ca II 8542 line core. The correspondence breaks down at smaller scales around <∼ 5 (more than twice the ∼ 2 beam size). In addition, the observed T 3mm b is higher than the temperatures inferred for the Ca II 8542Å line core (∼ 9500K vs. ∼ 7000K); One explanation is that the 3-mm continuum forms, on average, somewhat higher in the chromosphere than the Ca II 8542Å line and is thus sensitive to the rising temperatures at the base of the transition region (similar to the case for Hα discussed by Molnar et al. 2019).
By contrast, the Band 6 temperature maps exhibit significant structural differences from the temperatures inferred from the inversion of the Ca II 8542 line, even if these two diagnostics could be expected to form in more similar heights in of chromosphere. In fact, there was no layer in the temperature cube recovered by the inversions that showed a good spatial correlation with the 1.2 mm continuum image. Curiously, the mean temperature over the field view (∼ 7000K) was very similar for the Band 6 and for the inversions at the layer best corresponding to the Ca II 8542 line core. This correspondence of the mean temperatures but bad spatial correlation may indicate that the temperatures observed at 1.2 mm might be an average over a large spread of heights in the atmosphere that show different spatial structuring, but scattered around a common average temperature. Previous works have already suggested that Band 6 could have a disjoint contribution function with components at various heights (Loukitcheva et al. 2015;Chintzoglou et al. 2021). The observed Band 6 continuum emission at any given resolution element thus might include contributions from multiple different plasma elements along the line of sight, leading to a local breakdown of the correlation with the Ca II 8542Å line, which is formed over a narrower range of heights. This may be especially problematic in the current situation, where the spatial resolution in the millimeter observations is lower than in the optical and IR lines. Alternatively, this might indicate an offset of unknown source in the absolute temperatures derived by the inversions or the ALMA calibration (Jafarzadeh et al. 2019). We also compared the appearance of the Band 6 and Band 3 observations of the same region taken only 90 minutes apart. The differences in the observed structures seem to be greater than what would be expected from solar evolution or due to the different spatial resolution, which isn't fully consistent with the result from non-equilibrium-ionization simulations that these two continuum emissions are generally formed at the same heights in the atmosphere (Martinez-Sykora et al. 2020).
We further investigated the relationship between ALMA observations and spectral line diagnostics by performing joint inversions, in which the observed millimeter T b were included as an additional constraint in the fitting process. Inclusion of ALMA Band 3 data provided a seemingly reliable constraint on the temperatures inferred in the upper chromosphere, ensuring that there was a strong temperature rise at the top of the atmospheric range. Without the inclusion of the Band 3 temperatures, the recovered atmospheric profiles tended to remain flat or even go down slightly about the primary height of formation of the Ca II 8542 line core. Otherwise, the inclusion of Band 3 in combination with the two spectral lines had little effect on the inferred temperature stratifications lower down.
Inclusion of Band 6 data, on the other hand, drastically alters the inferred atmospheres, introducing locations with temperatures at chromospheric heights as low as 4000 K above the plage. The inversion also found regions of lower temperatures also lower down in the atmosphere, but generally away from the strong magnetic concentrations. We spatially partitioned the T 1.2mm b temperature map into regions of high (T b > 7900K) or low (T b < 5900K) temperatures. We found good spatial correlations for the regions of hotter temperatures with the corresponding inversion-derived temperatures above the height of formation of the Ca II 8542 line (log ξ = −4.4)) and for the lower temperatures regions with the temperatures in those spatial locations but much deeper in the atmosphere (log ξ = −2.9), just above the temperature minimum). In the regions of intermediate temperatures, there was again no height that showed a strong correlation between these two temperature measures. This was shown to be consistent with the changes in the distribution of heights for the response functions we calculated for Band 6 (and to some extent Band 3), ranging from log ξ = −1.5 for the coldest pixels to log ξ < −4 for the hottest. Even attempting to perform an inversion with the two spectral lines and both Band 6 and Band 3 brightness temperatures (the estimate for the latter being derived from simultaneous Hα data using the relationship found by Molnar et al. (2019) didn't improve the reliability of the inversions. It appears that the inversion process has some difficulty in determining the proper height of formation of the 1.2 mm continuum in any given pixel. This is likely related to another issue we examined in this work, the need to include the option in the STiC code to treat non-LTE hydrogen ionization in the process of solving for the underlying atmospheric conditions (and in particular the electron pressure). The effect of including this mechanism is to reduce the hydrogen ionization, and hence the electron density, for a given temperature or height in the chromosphere. This reduced density also manifests itself in a reduced collisional coupling between the Ca II 8542 source function and the local temperature. To reproduce the observed line core intensities for this line, the inversion process must infer a higher local temperature. Therefore, when we included this non-LTE ionization in the inversions, we recovered higher temperatures in the chromosphere heights corresponding to the Ca II 8542 line, more closely matching the temperatures measured in the 3 mm continuum with ALMA. However, even including this effect, it seemed the inversion fitting process struggled to reconcile the temperatures measured in the 1.2 mm continuum in particular, with the temperatures needed to recreate the observed spectral lines forming at similar heights.
The underlying issue may be that both of these diagnostics are somewhat sensitive to the electron density profiles through the atmosphere. For the Ca II 8542 line this is through the collisional coupling described above, while for the millimeter continuum the density strongly controls the contribution height(s) of the continuum emission. With a strong gradient in temperature with height, uncertainties in the electron density can lead to significant confusion in the inversion. Uncertainties in understanding the electron density profile in the solar atmosphere are a known problem, as indicated by the discrepancies in this parameter for the FAL-like models based on spectral-line intensities and those models calculated to explicitly fit the limb intensities in the radio and millimeter wavelengths (Ewell et al. 1993;Selhorst et al. 2005Selhorst et al. , 2019. The latter show elevated electron densities extending several thousand kilometers higher in the atmosphere. This difficulty in understanding this component of the atmospheric structure is exacerbated by the fact that the calculations underlying the inversion and line fitting process do not strongly constrain the determination of the actual electron density profile, since variations can be accommodated by changes in the temperature profile in order to achieve a reasonable fit. Further, the electron density depends on complicated physical processes that control the hydrogen ionization (Rutten 2017), which are hard to fully incorporate in a computationally tractable manner. Further developments in this area could include a better treatment of the electron density in order to better reconcile disparate sources of temperature measurements in the chromosphere.