Energy Levels of Singly Ionized and Neutral Hafnium

Improved energy levels for singly ionized and neutral hafnium of both even and odd parity are determined from Fourier transform spectrometer data using a least-squares optimization procedure. Data from interferometric spectrometers provide much tighter control of systematic uncertainties of line position measurements than can be achieved using dispersive spectrometers. The strong optical and near-UV lines connecting these levels are most likely to be used in the determination of isotopic abundance patterns. Comparisons of new results to published ones strongly suggest that our energy levels have systematic uncertainties in the mK (1 mK = 0.001 cm−1) range or smaller, and that widely used tables of energy levels for ionized Hf have systematic errors of approximately 70 mK.


Introduction
It is desirable, for several reasons, to improve the National Institute of Standards and Technology (NIST) Atomic Spectra Database (ASD) energy level precision from 10 to 1 mK, where 1 mK = 0.001 cm −1 . This has already been done for some lighter elements (Kramida et al. 2020). Here, we emphasize astrophysical application. Hafnium (Hf, Z = 72) is a neutron capture element with isotopes from both the r(apid)-process and the s(low)-process nucleosynthesis. In the solar system, Hf is approximately a 50%/50% mixture of r-and s-process material (Kobayashi et al. 2020). Many details of neutroncapture (n-capture) nucleosynthesis are still uncertain, including the relative role of n-star mergers in the r-process when the universe was young. There has been work on isotope shifts of lines of Hf II (the second spectrum or singly ionized spectrum) and of Hf I (the first spectrum or neutral atom spectrum) (e.g., Anastassov et al. 1994;Zimmerman et al. 1994;Zhao et al. 1997;Bouazza & Wilson 2000;Nieminen et al. 2002;Bouazza 2012). At this time, energy levels of singly ionized and neutral Hf in the NIST-ASD are available to precision of 0.01 cm −1 or 10 mK with no uncertainties (Kramida et al. 2020). This is satisfactory for efficiently finding Hf I and Hf II lines. If no uncertainty value is available in ASD, it is usually safe to assume that the probable error is between 2.5 and ∼25 units in the last place. About 90% of data in ASD satisfies this assumption. 1 More accurate energy levels and some information on isotope shifts are both needed for isotopic abundance measurements on distant stars. The motivation for measuring more and better isotopic abundance patterns in metal-poor stars is to explore nucleosynthesis throughout the Galaxyʼs history.
Large telescope observations and observations from space with the Space Telescope Imaging Spectrograph (STIS; Kimble et al. 1998;Woodgate et al. 1998) on board the Hubble Space Telescope (HST) are gradually making it possible to observe more isotopic abundances (e.g., Roederer et al. 2020). Such observations require both high spectral resolution and high signal-to-noise ratios (S/Ns). At least some information on isotope shifts is essential for measuring isotopic abundances. Ground-based telescopes with a 30-50 m primary mirror diameter are in the design and construction phase. These facilities, such as the Thirty Meter Telescope, the Giant Magellan Telescope, and the (European) Extremely Large Telescope, will extend the ability to measure isotopic abundance patterns. The S/N of spectra from current telescopes is limited by quantum noise. An order-of-magnitude increase in primary mirror area enables observers to increase the number of metal-poor stars for study, or shorten integration times, or improve the spectral resolution, or improve the S/N, or some combination of the preceding improvements using various tradeoffs. It is possible, perhaps likely, that improved techniques for measuring isotopic abundance patterns will be developed after 30-50 m telescopes come online. Improved energy levels are needed for spectral synthesis and especially for blending problems. The current choice of Hf is motivated primarily by the existence of an independent published set of center-of-gravity (COG) transition wavenumber measurements based on FTS data.
A further motivation lies in the analysis of lab spectra to determine reasonably accurate atomic transition probabilities. Analysis is often challenged by spectral line blends that must be resolved in some fashion to use the radiative lifetime + branching fraction (BF) method (e.g., O'Brian & Lawler 1996). Branching ratios for lines from a common upper level, which sum to unity, are called BFs. BFs can be divided by radiative lifetimes from laser-induced fluorescence measurements and turned into reliable atomic transition probabilities. One of the simplest and most efficient methods for blend separation is sometimes called the COG method. In this method, the COG of the blended feature is compared to Ritz (energy level difference) values for the wavenumbers of the transitions. Accurate Ritz wavenumbers for the desired spectral line and for the blending spectral line are necessary, along with a highresolution and high-S/N lab spectrum from a Fourier transform spectrometer (FTS) or a similar spectrometer. A single equation with a single unknown is used to determine the fraction of the blended feature from the line of interest. We have successfully used the COG method of blend separation in two recent transition probability studies in Hf II (Den Hartog et al. 2021, hereafter D21) and in Fe II (Den Hartog et al. 2019).
One of the desirable features of FTS data is the exceptionally linear wavenumber scale. With a proper compensator and a phase correction, only a multiplicative "rubber ruler" correction near unity of the wavenumber scale based on observed Ar II lines is typically needed to remove the effect of slight misalignment between the lamp beam and calibration laser beam and the effect of the FTS aperture. These concepts are discussed more fully in Section 2, below. The linearity of interferometric or FTS data is important to the control of systematic uncertainties in COG transition wavenumbers.
Hafnium has six isotopes with measurable abundance in solar system material. The lightest isotope 174 Hf has a low and uncertain abundance of 0.0016 ± 0.0012 or [0.0016 (12) (Meija et al. 2016) is followed by the nuclear spin I of odd isotopes (Stone 2005). For stable even isotopes, I = 0 (Sansonetti & Martin 2005). The odd isotope abundances and I values are sufficiently large to produce hyperfine structure (HFS) and complex spectral line shapes. Examples of complex spectral line shapes for lines of Hf II can be found in Figure 1 in Lundqvist et al. (2006, hereafter L06) and in Figure 1 of D21. The even isotopes are in the central spike, while the pedestal is due to the HFS of the odd isotopes. Not surprisingly, the spectral line profiles for some lines of Hf I are similar, an example of which is shown in Figure 1 of this paper. Hafnium is not unusual as a heavy metal with multiple isotopes. Hafnium is unusual in that it is quite difficult to separate from zirconium (Zr), due to very similar chemistry. Both Hf and Zr are found in the same ore, with Hf being in the minority, at 1%-4% of the Zr concentration. This results in line-blending problems, because most Hf samples have some Zr contamination. Modern separation techniques can overcome this contamination problem, but a search for strong known Zr lines should be used with all older spectral data. The shape of a Hf II line along with a Zr line is clearly visible in Figure 1 of D21.
Fortunately, there are completely independent and published FTS measurements of COG wavenumber values for lines of Hf II in Table 5 of L06. "Completely independent measurements" in this context means that a different instrument was used, different software was used, and different personnel participated. This led us to choose Hf II and Hf I for our first study to improve energy levels. We start with Hf II because singly ionized Hf is the dominant ionization stage in the photospheres of F, G, and K stars of interest (e.g., Den Hartog et al. 2019; Roederer et al. 2020). The study of L06 was primarily aimed at transition probabilities not wavenumbers and energy levels. The authors of L06 do, however, provide COG wavelengths and wavenumbers and claim to have "improved wavelengths" in the title of their paper. Their tables also have energy levels listed to 1 mK, but there is no mention of the method of their derivation, nor do they make claims to have "improved energy levels." The L06 COG wavenumbers are reported to 1 mK, but there are no uncertainties given. Comparisons of independent measurements provide an important estimate of possible systematic errors.

FTS Data and Ar II Calibration
The 1 m FTS on the McMath Solar telescope was used to record all of the data we analyzed for improved energy levels. Hollow cathode lamps were used as sources. Table 1 summarizes information on the FTS data, which are available from the National Solar Observatory (NSO) electronic archive. 2 There is limited FTS header data in the electronic archives of the NSO. However, the high discharge currents indicate that the water-cooled demountable (home-made) hollow cathode lamp was used. These FTS data have been studied previously to measure branching fractions of Hf II lines (Lawler et al. 2007, hereafter L07, D21). Our familiarity with these spectra of relatively broad coverage and high S/N is valuable.
Research-grade FTS instruments can have an exceptionally linear wavenumber scale. This is the result of recording an interferogram while moving a mirror along the instrument axis with the use of a single-frequency reference laser to track the mirror position. There are at least two effects that lead to the need for a multiplicative or "rubber ruler" correction of the wavenumber scale. A very slight misalignment (on the order of 10 −3 -10 −4 radians) between the reference laser beam and the beam from the lamp is an important contributor to the multiplicative correction. This misalignment effect contributes a cosine or inverse cosine factor that is typically very close to unity. The other important "rubber ruler" correction is typically called an aperture correction. The etendue or (solid angle) times (cross sectional area) of the lamp beam is generally much larger than that of the reference laser beam. A larger range of cosine factors are necessary with a large etendue. Both of these effects can be removed with a single multiplicative correction very near unity.
Optical to UV research-grade FTS instruments have thick beam splitters of various materials depending on the desired wavelength coverage. Thick beam splitters provide the required fraction of a wavelength stability, but introduce dispersion into one beam of a Michelson interferometer. The compensator is an uncoated copy of the beam splitter and is used to introduce equivalent dispersion in the other beam of the interferometer. Residual phase errors, which could move part or all of the desired signal from the real part of the transform into the discarded imaginary part of the transform, are removed with a final phase correction. Zero path difference is the same for all wavelengths with a good phase correction.
The most widely used, and likely best, method for finding the multiplicative correction is to measure selected strong Ar II lines with accurately known wavenumbers and high S/Ns in the spectra of interest. The authors of L06 used Ar II lines with wavenumbers from Whaling et al. (1995). The preferred set of 28 Ar II lines for wavenumber calibration was proposed earlier by Learner & Thorne (1988). These Ar II lines are strong in spectra from hollow cathode lamps and typically have high S/ Ns, usually >100. The 28 lines of Ar II with wavenumbers from Whaling et al. (1995) are used to determine the multiplicative correction factor near unity for the FTS spectra of Table 1. Sansonetti (2007) has disputed the accuracy of the Ar I lines measured by Whaling et al. (2002), but endorsed the earlier Ar II measurements by Whaling et al. (1995) relative to precise wavelengths of CO molecular lines. Sansonetti (2007) at NIST concluded with a correction factor for putting the Ar I lines measured by Whaling et al. (2002) on the same absolute scale as the Ar II lines measured by Whaling et al. (1995).
The exceptionally linear wavenumber scale of FTS instruments immediately leads to the question of why FTS instruments are not used, except in some special cases, on large astronomical telescopes. Multiplex or photon noise is spread evenly throughout FTS spectra. Photon noise is not spread out using a large dispersive or echelle spectrometer, and thus echelle instruments are preferred for general purpose astrophysical instruments. Although a large echelle spectrometer coupled to a telescope may have a mosaic of multiple echelle gratings to match the etendue (or solid angle times area) of the telescope, the echelle spectrometer can, in principle, achieve the resolving power of a single echelle grating. It is very difficult and unnecessary to phase echelle gratings in a mosaic to make a single very wide diffraction grating. The diffraction-limited resolving power R of one echelle grating is high enough, similar to our 3 m focal length echelle spectrometer (Wood & Lawler 2012). The width W = 25.4 cm and blaze angle 63.5°at 30,000 cm −1 yield R = 30,000/(1/(2 × sin(63.5°) × 25.4)) = 1.36 × 10 6 . This estimate of the diffraction-limited resolving power at 3300 Å is from interferometry. It is based on the maximum path difference from the width of a single grating and the blaze angle. The diffraction-limited resolving power requires a very bright source to achieve a high S/N and very small pixel size if a  detector array is used. Large echelle spectrometers do not lack in resolving power, but they do not have the exceptional wavenumber linearity of FTS instruments. Systematic uncertainties of interferometric data are certainly smaller than those of older data recorded using dispersive spectrometers and film. The estimate of systematic uncertainty requires both skill and experience. One systematic that has been discussed, but not completely resolved, is due to the possible differences of illumination of the interferometer entrance aperture by Ar II lines and sputtered Hf lines. This systematic effect could have some dependence on the hollow cathode design and discharge current, the FTS design, and other parameters. The good agreement with COG measurements in L06 discussed below strongly suggests that this systematic effect is smaller than a mK for many levels. We have scrutinized our data for any indication of systematic differences between the six spectra and found them to be small. The offset between transition wavenumbers measured on each spectrum and the mean from all spectra, averaged over all lines, has been tabulated and found to be small for all spectra. The highest average offset is ∼9 × 10 −4 cm −1 on spectrum 2. Any systematic differences between spectra arising either from the alignment, spectral calibration, or from lamp operating parameters is on the level of a mK. Another source of small systematic uncertainty, 1 mK or less, is due to weak offdiagonal hyperfine structure which could be omitted from numerical integration if it has melted into background noise of the FTS data (e.g., Kramida et al. 2017). Integration limits in our work were chosen to be sufficiently wide as to include such structure unless another line was nearby. Final phase correction is also a source of small systematic uncertainty unless it is nearly perfect.
We have chosen this subset of Hf II lines because they have satisfactory (10-100 or more) S/N and appear clean or unblended in the archival spectra. They are transitions observable in a range of stellar spectra. Typically, the strongest lines from each upper level were chosen with an average between three and four lines per odd parity level. The average number of lines per even parity level was >5, but the decline was steeper as the even parity level energy increased. Additional work on near-infrared lines is certainly needed but is beyond the scope of this work. The creditability of this work is enhanced by comparison to independent COG measurements in L06. Table 2 includes COG values for 123 strong Hf II lines, with the multiplicative correction from Ar II lines included, from Spectra 1 through 6 of Table 1. Mean COG values along with standard deviations for the distribution and for the mean are included in Table 2. Mean is used when multiple spectra are included, while average is used when both multiple spectra and multiple lines are included. Also included in Table 2 are air wavelengths computed from our mean COG wavenumber using the standard index of air given by the five-constant, threeterm approximation (Equation (3)) of Peck & Reeder (1972) as well as the upper level, lower level, and Ritz wavenumber from the NIST ASD (Kramida et al. 2020). The final two columns are Ritz wavenumbers and their uncertainties from the optimization procedure described below.
L06 reported improved transition wavelengths in Hf II, but wavelengths are typically quoted in air, whereas the spectral output of an FTS is in vacuum wavenumbers. The energy levels of Table 6 of L06 have a 0.001 cm −1 precision, but are not described in any detail and may not have been improved in a self-consistent fashion as in the present work. We therefore make comparison to the vacuum COG wavenumbers from Table 5 of L06, which are measured directly from their Ar II calibrated FTS spectra. This table includes many weak lines from upper levels of interest. Some lines also have notes indicating that they are problematic, including 10 lines designated as self-absorbed (SA) and 11 lines designated as blends (bl). These 21 lines and a curious outlier at 23,620 cm −1 are dropped from further consideration. The mentioned outlier, originating from the upper level at 42,518 cm −1 , stands out in a comparison of old wavenumber measurements from the same upper level reported by Meggers & Scribner (1934). It is possible that the transition near 23,620 cm −1 is a misprint in L06. The remaining 174 lines from Table 5 of L06 are of interest in our study. Of the 174 lines, 52 overlap with our selected set of 123 strong lines of Hf II. The inclusion of many weaker lines from upper levels of interest in L06 serves to dilute the effect of Zr contamination on an optimized energy level determination, but increases the standard deviation of the (measured COGs-Ritz wavenumber). The average offset for the 52 lines in common between our measured COGs and those from L06 is less than 1 mK, and the standard deviation of the offset is between 3 and 4 mK. (We use mK in this text, for the readerʼs convenience, and cm −1 in all tables.) It is possible that some of the 174 lines of L06 had relatively low S/N and might better be culled for our energy level study, but no information is available on the S/N. Weaker branches from a selected upper level may have larger uncertainties unless multiple spectra were recorded at higher discharge current and analyzed to yield a similar S/N. All of the COGs of L06 are reported to 1 mK. A decision on which additional lines to cull is thus too subjective, and we have retained the full set of 174.

Least-squares Refinement of Energy Levels from Inversion of Large Matrices
Kramida (2011) describes the determination of least-squares optimized energy levels from a set of transition wavelengths or wavenumbers and provides a program (LOPT) for doing such an optimization. The LOPT method requires the inversion of a matrix of dimension N × N followed by operation on a vector with N elements. The number N is the total number of energy levels to be determined of both even parity and odd parity. Because N is typically larger than 50 and could in principle range up to 1000 or more, the simplest matrix library inversion routines are not sufficiently accurate. Although we studied this reference closely, a coauthor of the present study preferred to write his own code using familiar software instead of using LOPT. Alexander Kramida, acknowledged below, provided many tests of our mean COGs and guidance of our work. He also ran LOPT with our data and found that LOPT returned energy levels equal to ours within ±1 in the least significant digit. Our work represents an independent test and verification of the LOPT program. We are using the Cholesky decomposition method. Kramida (2011) presents formulas and discussions of the quantities D1, D2, and D3, which are all uncertainties of leastsquares adjusted energy levels. The uncertainty D1 best represents the uncertainty of an optimized energy level with respect to other levels connected via a strong optical or UV transition. Thus, D1 is a minimum uncertainty. D2 is the level uncertainty with respect to the ground level at 0.000 cm −1 . D3 is the level uncertainty with respect to the lowest fixed level,  . b Wavelength in air is calculated from the mean observed center-of-gravity wavenumber and the standard index of air from Equation (3) of Peck & Reeder (1972). c Upper-and lower-level energies and Ritz wavenumber used to identify transitions are taken from NIST ASD (Kramida et al. 2020).
(This table is available in its entirety in machine-readable form.) but not the ground level. D3 is not applicable in our case, but is used when some excited levels are weakly connected and have a different spin than the ground level. We have included evaluation of D1 and D2 statistical uncertainties in our energy level tables.

Methodology
We followed the formulation of the optimization as a linear model described in an authoritative text by Graybill (1961) and by Radziemski et al. (1972): 1. Inputs: Transition wavenumbers, y = {y j |j = 1, L , p}, their variances s j 2 , and a p × N array x consisting of zeros and plus/minus ones. (Here, p is the number of transition inputs and N is the number of levels to be determined.) For example, if˜b b = -  y j k ℓ where the tilde denotes the optimal values of transition wavenumbers and levels, then x j,k = 1 = − x j,ℓ and x j,a = 0 for a ≠ k, ℓ. Bold face denotes column vectors and matrices, superscript T denotes transpose (row vectors or transposed arrays).

Outputs: Optimal energy levels
and the numerical value of s , which is a parameter used in the optimization of the level variances, not the actual energy levels, that is a measure of the goodness of fit, and an array S −1 of correlations.
We define two sets of random variables, where E[q] indicates mean or "expectation value" of the quantity q with respect to the probability distribution function f. The (log of) the probability distribution function is given by in which V is a p × p matrix of transition weights. Here, f ln is optimized with respect to β i and σ: The solutions are the maximal likelihood estimates b and s  of the variables β and σ. These are given by the equation for s  2 , whereỹ i is the Ritz wavenumber for the transition computed from the optimized levels, and the equations for the optimized levels that define the N × N matrix S (called W in Kramida 2011), There are two differences between our method and that of the LOPT program of Kramida (2011). The first difference is that we use as a basis the space of transition wavenumbers, not levels, so that our vector space is larger. The second difference is the additional optimization variable σ 2 . A simple scaling of the weight matrix, which is what σ 2 provides, has no effect on the optimized energies. We have confirmed this statement. Radziemski et al. (1972) clearly state: "A proportional change in all the weights will not affect the estimated a j and b j (optimized energies) at all, simply changing the value of the parameter σ 2 . Individual changes do matter, hence it is important that the ratios of the weights be correct." We can compute variances and covariances of levels, in the standard way, as the average deviation (squared) from the optimized (mean) values of the levels, expressed as compactly as expectations with respect to the probability distribution f above (the mean of a set of β i values is b  i , the optimized values): Acknowledging that the estimation of uncertainties of the input wavenumbers is based on standard deviation of the mean of several measurements, we set σ 2 = 1 in this formula, if the fit is such that s   1, and report a conservative error bar (equivalent to Equation (26) of Kramida for the uncertainty in the absence of systematic shifts): We find typical values of s  ranging from 0.68 to 0.91, indicating that the fit is good, but the D2 values in our tables are calculated using σ 2 = 1. The uncertainty in a Ritz Again, it is important for a reader to recall that the use of interferometric spectra instead of the older dispersive spectra has the primary advantage of reducing dominant systematic uncertainties and statistical uncertainties. Statistical uncertainties are quite small in this paper. Systematic uncertainty, even using interferometric data, is likely to be at the mK level.
Our coding was done in R using the Rmpfr arbitrary precision library to ensure that the inversion would be robust with numerical errors much smaller that the input data precision. Precision was set to a value such that S −1 · S has summed absolute values of off-diagonal elements no larger than 10 −19 .

Energy Level Optimization for Singly Ionized Hafnium
The 123 transition COGs determined using Hf II lines from Spectra 1 through 6 yield a set of least-squares optimized energies for 24 levels of even parity and 35 levels of odd parity. The ground level is 0.000 cm −1 by definition, and thus only 58 energy levels are determined. The weighting factor for mean transition COGs is the inverse square of the sum of the standard deviation of the mean unweighted COG measurements from Spectra 1 through Spectra 6 and an additive constant. The simplest estimate of the uncertainty of the means would be to use the inverse square of the standard deviation without an additive constant. However, we found that, due to the small number of spectra, five or six, contributing to the mean wavenumber, the standard deviation was unrealistically small for some lines, resulting in those lines being too heavily weighted. The underestimate of the wavenumber uncertainties is confirmed by calculating the residual sum of squares (RSS) divided by the number of degrees of freedom (N df ): RSS/N df for the data set. 3 This statistical measure (discussed in detail in Kramida 2011) should be near unity, but without the additive constant is 1.6. Using the additive constant in the uncertainty brings RSS/N df down to 0.9 indicating the uncertainties to be reasonably estimated.
This additive constant of 0.5 mK applied to the standard deviation of the wavenumber mean is a lower limit for the wavenumber uncertainty. It might be considered a dark uncertainty, but it is due to the small number of spectra included in our mean. Thus, it is not completely "dark," and the question arises as to whether simple addition or addition in quadrature should be used to combine the constant and the standard deviation of the mean. We have used simple addition, which is more conservative, but the effect of switching to quadrature is small. In all cases, the change caused by switching from addition to quadrature lowers the D 2i uncertainty by at most a factor of 1/3, and for most levels considerably less. It is important for the reader to recall that the use of interferometric spectra instead of the older dispersive spectra has the primary advantage of reducing systematic uncertainties-not statistical uncertainties. Statistical uncertainties are quite small in this paper. Systematic uncertainties are difficult to estimate in most cases, but the comparison to independent results from L06 is critical for estimating systematic uncertainties. The systematic uncertainties are in many cases 1 mK or less, and in the worst cases will be only a few mK for our refined energy levels.
The optimized levels are tested by tabulating offsets between measured COG mean values and Ritz transition wavenumbers calculated from differences of optimized energy levels. The comparison is summarized in row 1 of Table 3. The "mean COG-Ritz" comparison for any transition that solely determines an energy level is dropped from the comparison, because it contributes zero offset by definition. The average offset for 115 remaining comparisons is much less than 0.1 mK = 0.0001 cm −1 , and a standard deviation of the offset is ∼2 mK for the 115 lines.
The 174 transition COGs from L06 are analyzed in the exact same fashion except that, due to lack of information on the S/ N, the COG measurements are all weighted equally. This analysis yields least-squares optimized energies for 27 levels of even parity (not counting the defined ground energy level at 0.000 cm −1 ) and 18 levels of odd parity. These are tested by tabulating offsets between measured mean COG values and Ritz transition wavenumbers calculated from differences of optimized energy levels. The comparison is summarized in row 3 of Table 3. The average offset for 171 comparisons is much less than 0.1 mK, and the standard deviation is ∼13 mK. Any transition that solely determines an energy level is dropped in computing the standard deviation, because it contributes zero offset by definition. The relatively large standard deviation is due to some combination of lines with poor S/Ns equally weighted with stronger lines, and possibly some contribution from Zr blends. Row 2 of Table 3 reminds the reader that there are 52 transition COGs in common between our line list and the 174 transition COGs from L06. The average offset for these 52 transition COGs is −0.9 mK and the standard deviation of the Table 3 Comparison of Mean Measured COG Wavenumbers from Spectra 1 through 5 or 6 of UW to Ritz Wavenumbers and to Wavenumbers from Notes. a The comparison is for 115 measured mean COGs from spectra 1 through 5 or 6 and Ritz wavenumbers from this study. (In this comparison and in c and d below, any transition that solely determines an energy level is dropped because it contributes zero offset by definition.) The matrix weighting of each mean measured COG is 1/(σ mean + 0.0005 cm −1 ) 2 . See text for discussion. b The comparison is for 52 transition wavenumbers in common between this study and that of 3 RSS and N df are given by Kramida (2011). In our notation, offset is 3.3 mK. This small average offset indicates that we are reproducing measurements from L06 sufficiently accurately. The last row of Table 3 is devoted to Hf I and will be discussed below.
Tables 4 and 5 compare energy levels of even parity and odd parity, respectively, from the NIST ASD with the least-squares optimized values from our COG measurements of 123 strong transitions of Hf II and with the least-squares optimized values from COG measurements of 174 transitions of Hf II as reported in L06. We made COG measurements of lines from additional levels including 42,770.56, 45,643.25 cm −1 , and numerous other higher levels, but the S/Ns did not yield a sufficiently small standard deviation of the mean COG for inclusion in this study. Although the spectra of Table 1 list broad UV and near-IR coverage, again the S/N of Hf II lines at shorter UV and longer IR wavelengths did not yield a sufficiently small standard deviation of the mean COG for inclusion in this study.
The problem of Zr contamination in Hf samples is not new. Consider the Hf II paper by Meggers & Scribner (1934). They discussed on page 636 a highly purified oxide sample that they acquired and used. Later on the same page, Meggers & Scribner concluded that the sample, while better than earlier samples, "contained considerable amounts of zirconium and nickel." The selection of 123 strong lines of Hf II in our analysis above should reduce effects from Zr blends. However, we still have concerns, and have checked in the MIT wavelength table (Phelps 1982) for Zr lines within 0.5 cm −1 of our selected Hf II lines as suggested above. Two of our 123 strong lines of Hf II are of particular concern, including the line near 27,789 cm −1 from the upper level near 43,044 cm −1 and Note. a D1 and D2 are statistical uncertainties relative to connecting levels and relative to the ground level, respectively, in this Table and Tables 5, 7, and 8. These statistical uncertainties are discussed in the text and fully defined in Kramida (2011). Although occasionally D1 > D2, the level uncertainty with respect to the ground level should be the maximum of D1 or D2. Systematic uncertainties are more difficult to evaluate, but are expected to be smaller using interferometric data than older dispersive spectrometer data.
the line near 26,539 cm −1 from the upper level near 46,674 cm −1 . These transitions solely determined the level energies for these two upper levels in the study using our COGs, whereas the transitions were labeled as blends in L06 and thus not included in our 174 line comparison. The offset for these two instances between the energy levels of Table 5 determined from 174 line COGs of L06 in the column labeled "Levels from L06 COGs" and the energy levels from our COGs in the column labeled "Levels from UW COGs" provides an estimate, likely high, for the effect of a Zr blend in our COGs. The inclusion of . b Wavelength in air is calculated from the mean observed center-of-gravity wavenumber and the standard index of air from Equation (3) of Peck & Reeder (1972). c Upper-and lower-level energies and Ritz wavenumber used to identify transitions are taken from NIST ASD (Kramida et al. 2020).
(This table is available in its entirety in machine-readable form.) weaker lines with equal weight in Table 5 of L06 should dilute the effect of Zr blends but increase the standard deviation of (measured mean unweighted COG-Ritz wavenumber). The offsets are 7 mK and 6.4 mK for the 43,044 cm −1 and 46,674 cm −1 levels, respectively. A simple average of the energy level listed in the "Levels from L06 COGs" column with the energy level listed in the "Levels from UW COGs" column might be appropriate for those two levels.
Isotope shifts are dominated by mass shifts at the top of the periodic table and by nuclear field shifts at the bottom of the periodic table, if an s-electron is involved in the transition. Although isotope shifts are small near the middle of the periodic table, including those of Hf, the existence of independent published COG measurements from FTS data for Hf is a huge advantage for demonstrating the reproducibility of such measurements.

Energy Level Optimization for Neutral Hf
We measured COGs for 123 lines of Hf I as listed in Table 6, which includes the same information for the Hf I COGs as Table 2 includes for the Hf II COGs. These lines have satisfactory (10 to 100 or more) S/N in Spectra 1 through 5 or 6 and appeared clean or unblended. Improved energy levels for nine even-parity levels of neutral Hf are given in Table 7, and 60 odd-parity levels of neutral Hf are given in Table 8. The exact same least-squares optimization of energy levels used for singly ionized Hf is used for neutral Hf with weighting based on the inverse square of the sum of the standard deviation of the measured mean unweighted COG measurements and an additive constant. As in the Hf II analysis, the additive constant of 0.0005 cm −1 applied to the standard deviation of the mean wavenumber is needed to reduce the variation of the weighting factors due to the relatively small number of spectra, 5 or 6, contributing to the mean wavenumber. The additive constant is a lower limit for the wavenumber uncertainty. As with Hf II, we compare the statistical measure RSS/N df for the cases with and without the additive constant. Without the additive constant, this ratio is 3.4, and with the constant added, it reduces to 1.1, indicating the uncertainties are now reasonably estimated. The optimized levels are tested by comparing the mean COGs to the Ritz wavenumbers determined from the improved energy levels for 96 cases. As with singly ionized Hf, the "observed-Ritz" comparison for any transition that solely determines an energy level is dropped in computing the standard deviation, because it contributes zero offset by definition. The average offset, reported in the fourth row of Table 3, is much less than 0.1 mK = 0.0001 cm −1 , and the standard deviation is ∼1.3 mK.
Our concern about possible Zr blends could not be reduced by an independent set of COG measurements from FTS data using a range of strong and weaker transitions as was the case for lines of Hf II. We have thus culled 12 Hf I lines of interest due to Zr I or Zr II lines in the M.I.T. wavelength tables (Phelps 1982). The 123 lines of Table 6 remained after culling 12 from the initial 135 lines.
The last digit of NIST ASD energy levels for singly ionized and neutral Hf is 10 mK. Our refined energy level determinations in Tables 4, 5, 7, 8 have average offsets ± standard deviation from NIST ASD values of −40 mK ± 40 mK, −70 mK ± 30 mK, −10 mK ± 10 mK, −30 mK ± 20 mK, respectively. The average offset of −70 mK is significantly larger than a few times the last digit, although it is within the ASD probable error goal of 25 times the last digit. The good agreement we found with mean COG wavenumbers for Hf II reported by L06 indicates that the older wavelengths used to determine NIST ASD energy levels (Meggers & Scribner 1934) have a systematic error, which is not uncommon in wavelength measurements from the 1930s. Typical isotope shifts, which depend primarily on the electron configurations of the upper and lower level of a transition, are comparable to, and often smaller than, the ∼70 mK systematic offset.

Summary
We have produced improved energy levels from selected lines of Hf II and Hf I using FTS data from hollow cathode discharges. A total of 24 even-parity and 35 odd-parity levels of singly ionized Hf are included, not counting the ground level defined at 0.000 cm −1 . A total of nine even-parity and 60 oddparity levels of neutral Hf are included, not counting the defined ground level. These energy levels include low-lying levels of both even parity and odd parity that are most likely to be used in isotopic abundance studies. These data are intended to supplement the more comprehensive energy level data in the