Estimating index of refraction for specular reflectors using passive polarimetric hyperspectral radiance measurements

. Results of a method of estimating index of refraction from passive, polarimetric hyperspectral imaging radiance measurements are presented. As off-nadir viewing hyperspectral imaging platforms gain prominence, estimating index of refraction, which is invariant to viewing angle, may prove advantageous to estimating the emissivity, which is not. Results show that index of refraction can be retrieved to within 8% rms error for fused silica and sapphire glass targets, while simultaneously estimating object temperature. The accuracy and self-consistency of this technique for estimating index of refraction are shown to compare favorably to the maximum smoothness temperature – emissivity separation algorithm. Additionally, the results show that atmospheric downwelling radiance can also be accurately estimated, to within the noise of the instrument, concurrently with index of refraction. © The Authors


Introduction
One of the primary challenges in material classification and ID is dealing with target variability, that is, the same material producing different signatures based on scene conditions, illumination, viewing angle, etc. Target variability is often dealt with using a subspace, which is a group of signatures, to describe the target.As these subspaces become larger, however, it becomes more difficult to distinguish between spectrally similar materials.Most material ID work done in the longwave infrared (LWIR) relies on some form of temperature-emissivity separation (TES) algorithm.Generally, the retrieved spectral emissivity is treated as an invariant quantity unique to a material, i.e., one signature per target.Emissivity, however, can vary with viewing angle 1 and surface morphology, 2 making it necessary to employ a subspace to fully describe a target.This paper will deal primarily with the variability due to viewing angle and how estimating index of refraction, in place of emissivity, may be advantageous in certain scenarios.
To date, most hyperspectral imaging (HSI) platforms have been nadir or near-nadir viewing, so the target variability due to changing viewing angles was usually not an important factor.With additional emphasis being placed on sensing within A2/AD environments, off-nadir viewing sensors and algorithms to process data from these platforms are becoming increasingly important.When dealing with these oblique viewing geometries, the target variability due to the viewing angle of the sensor relative to the target's surface normal can no longer be ignored.
Other efforts to remotely estimate index of refraction are available in the literature.Thilak et al. [3][4][5] used polarization to simultaneously estimate index of refraction and surface normal angle.Their work only focused on broadband measurements, however, limiting its utility in material identification.Huynh et al. 6 examined using multispectral measurements to estimate index of refraction, using a model to reduce the number of parameters.While using a model to emulate the spectral variation of the refractive index is something leveraged in this work, their model did not account for materials with complex indices of refraction.Fetrow et al. 7 used polarimetric hyperspectral imaging (P-HSI) measurements to estimate index of refraction.Their method, however, estimated index of refraction independently at each wavelength leading to very noisy retrievals with fit uncertainties exceeding 1 (index of refraction is unitless) at some wavelengths.Additionally, they assumed a priori knowledge of the target surface temperature and the downwelling radiance.
The goal of the research presented here is to solve for index of refraction in a way that is both accurate and invariant to scene conditions without the need for control over the illumination source.Currently, this work is limited to specular, likely manmade, targets but future work will examine fitting semidiffuse targets as well.A fitting routine is developed using spectrally resolved polarimetric radiance measurements, often from multiple viewing geometries, to solve for index of refraction, surface temperature, and downwelling radiance.The fit models the spectral variation in index of refraction using physics-based constraints to reduce the number of parameters, making a customarily underdetermined problem into an overdetermined one.While previous work has demonstrated this capability in a laboratory setting, 8 this paper will demonstrate how this can be done outdoors with a more realistic and spectrally structured downwelling radiance profile.This paper will first introduce the theory underlying this research, then discuss how index of refraction is estimated, and, finally, present some results as well as a comparison with a TES algorithm.

Theory
In the LWIR, the at-sensor radiance comes from three sources: downwelling radiance reflecting off the target, radiance emitted by the target, and radiance emitted by the atmosphere along the line of sight.Presuming that the atmosphere is unpolarized in the LWIR, 9 the total radiance in each polarization state arriving at the sensor can be expressed as E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 1 ; 6 3 ; 6 4 1 where τ a and L a are the atmospheric transmission and path radiance, respectively, ϵ denotes the emissivity, ρ represents the reflectivity, BðT e Þ represents the Planckian blackbody radiance at the temperature of the material's surface, and L d is the downwelling radiance, incident on the object.All of these are spectral quantities, denoted in wavenumbers, ν.Subscripts s or p indicate the polarization state.
Although dealing with rough surfaces is the focus of future work, this research deals exclusively with specular reflectors.As such, a material's reflectivity-and thus, for opaque materials, their emissivity-can be related to index of refraction through the Fresnel equations 10 E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 2 ; 6 3 ; 4 5 7 where Ñ ¼ n þ iκ represents the complex index of refraction and θ i represents the angle of incidence.
Incorporating the Fresnel expression for reflectivity, Eq. (1) has 5N þ 2 unknown parameters, where N is the number of spectral points: τ a ðνÞ, L d ðνÞ, L a ðνÞ, nðνÞ, κðνÞ, T e , and θ i .Measuring P-HSI radiance values gives 2N (as discussed below) unique measurements, indicating that this problem is highly underdetermined.τ a ðνÞ, L d ðνÞ, and L a ðνÞ, however, can be either modeled or calculated using in-scene atmospheric correction techniques, significantly reducing the number of parameters.Furthermore, nðνÞ and κðνÞ are fundamentally linked via the Kramers-Kronig relationship.Finally, the spectral variation in nðνÞ and κðνÞ can often be described using a small number of parameters through models, such as the Lorentz oscillator model. 11Incorporating these constraints potentially convert a highly underdetermined problem into an overdetermined one.
It is common to express polarimetric information in terms of a four-element Stokes vector.For remote sensing applications, the fourth element, which represents the amount of circular polarization, is almost always ignored. 12 From the full P-HSI dataset, a linear Stokes vector is calculated at each spectral point and for each pixel.The first three elements of the Stokes radiance vector are expressed as E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 3 ; 3 2 6 ; 7 5 2 Note that, while there are three measurements at each spectral point, S 1 and S 2 are highly correlated, meaning they do not represent unique measurements at each spectral point.ϕ represents the azimuthal angle of what the sensor considers vertical polarization relative to the plane of reflectance.This is often referred to as the angle of polarization (AoP) and adds another parameter that needs to be determined.AoP can potentially give information about the azimuthal angle of a target relative to the sensor, but measurements are generally very noisy and are not utilized in this work. 6It is therefore beneficial to remove this ϕ dependence from the fit parameters.Instead, the total linear polarization, , is used in place of S 1 and S 2 when fitting index of refraction.

Solving for Index of Refraction
The Lorentz oscillator model is commonly used to describe index of refraction of materials, which treats the electrons in a material as harmonic oscillators being driven by incident radiation.Different materials produce different resonant frequencies, which determine the dielectric constant as a function of frequency.The functional form of the Lorentz oscillator is 13 E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 4 ; 3 2 6 ; 3 7 9 where ε ∞ is the dielectric constant as frequency goes to infinity, ν p is the plasma frequency, f j is the strength, ν j is the frequency center, and Ξ j is the damping constant of the j'th oscillator.
One of the drawbacks of this model is that it is difficult to replicate index of refraction of amorphous materials. 14These materials tend to have broader more slowly varying refractive indices that are difficult to describe using the functional form of the Lorentz oscillator.A number of different models exist to better describe amorphous materials, see Refs.14-16 for some examples.Many of these have been tested for this research but were found to be either too slow for the purposes of this research, taking several minutes to converge to a solution using the optimization technique detailed below, or limited in the number of materials they could accurately model. 17Instead, as part of this work, a method was developed to solve for the imaginary component of index of refraction at a few equally spaced points (or knots) and then use interpolation to resample onto the spectral axis of the measured values.The Kramers-Kronig relationship is then used to solve for the real component of index of refraction.More information on this model can be found in Ref. 8.
While this fitting model works well for amorphous materials, results show that it does not perform as well as the Lorentz oscillator model for crystalline materials.For this work, the model used to fit index of refraction is selected by the user.An error metric based on the measured and modeled S 0 and P is defined as follows: E Q -T A R G E T ; t e m p : i n t r a l i n k -; e 0 0 5 ; 6 3 ; 6 8 6 EðxÞ ¼ where x is a vector of the fit parameters.Included in x are parameters describing index of refraction model, a variable for object surface temperature as well as four variables describing the atmosphere.An initial estimate and bounds for the surface temperature can either be user-supplied or calculated based on the measured radiance and extreme cases in the atmospheric lookup table.For all results presented here, temperature bounds of 300 AE 15 K were used.The initial estimate for the index parameters set κ to 1 across the band when using the interpolation method.When using the Lorentz oscillator model, the initial guess is to have oscillators of equal strength spaced equally across a spectral region extending AE200 cm −1 beyond the measured spectral band.Because there may be several local minima in the search space, the Nelder-Mead search algorithm is used to minimize the error metric given by Eq. ( 5).
The four atmospheric parameters are temperature at the ground, temperature lapse rate, water concentration at the ground, and an ozone scaling factor.The temperature profile with altitude is defined by the ground temperature and lapse rate up to either an altitude of 16 km or the temperature of 216.7 K. Above that altitude, the standard model values are used.The water vapor profile is defined using the equation for water vapor scaling in Ref. 18.The ozone scaling factor is a unitless scalar, which multiplies the ozone concentration at each altitude.Temperature varied from 0°C to 35°C, water vapor concentration from 1000 to 10,000 ppmv, and ozone scaling from 0.5 to 20 in 20 increments.Temperature lapse rates of 4, 5, 6, 6.5, and 7 K∕km were used as well.
Adjusting these parameters, a four-dimensional lookup table of atmospheric downwelling values was generated.For now, since testing was done at relatively short ranges, it is assumed that τ a ¼ 1 and L a ¼ 0, but this same lookup table could be used to determine a transmission and path radiance or in-scene atmospheric compensation methods could be applied.There is also an option to manually input an atmospheric profile.This work is not intended to develop a new means of atmospheric compensation.Instead, this provides an example of how the downwelling radiance can be estimated concurrently with the refractive index.For most of the results presented, except where otherwise stated, it is assumed that the atmospheric downwelling is known a priori.

Instrumentation
Data were collected using a Telops longwave hypercam imaging Fourier transform spectrometer with a linear polarizer mounted at the entrance pupil.Spectral resolution can be adjusted anywhere between 0.25 and 150 cm −1 .The camera used a 320 × 256 mercury-cadmium-telluride focal plane array, with sufficient responsivity from 875 to 1250 cm −1 .
Radiometric calibration was performed using two onboard blackbodies.A spectral axis correction was done by scaling the wavenumber axis, so atmospheric features (namely the peaks of water lines) matched their known locations.After this, a simple linear interpolation was used to resample the data on to a common spectral axis.Polarimetric information was collected by measuring a scene through the linear polarizer set to a number of different angles.For this work, the modified Pickering method was used for data collection, so the scene was measured with the polarizer set to 0, 90, 45, and 135 deg. 12It has been shown by Holder 19 that, under the assumption of an ideal polarizer, a two-point blackbody radiometric calibration, as outlined in Ref. 20, is sufficient to compensate for polarimetric effects of the instrument.The polarizer used is nearly ideal, with a peak extinction ratio of 350:1, so it is assumed that the polarimetric effects of the instrument are properly compensated for using the radiometric calibration.

Results and Discussion
To test this method of estimating index of refraction, a fused silica wafer and sapphire glass window were observed at different times of day from a number of viewing angles.For the results presented here, data were collected at 1 cm −1 spectral resolution.Downwelling radiance was estimated using a gold mirror placed in the scene and tilted to the same angle as the targets.

Fitting with Known Downwelling Radiance
Results of index of refraction retrieval on the fused silica wafer are shown in Fig. 1.For this, retrieval viewing angles of 30 deg, 50 deg, and 70 deg relative to the target's surface normal are used in conjunction with one another.A 96 pixel window corresponding to the wafer in each image is used, and the fit is performed independently on each pixel.The rms error between the true and retrieved is 0.134 or 7.0% in n and 0.098 or 7.9% in κ.The spectral angles between the retrieved and ellipsometry-measured spectra are 3.24 and 4.40 deg in n and κ, respectively.Furthermore, the fit is fairly consistent pixel-to-pixel.The rms standard deviation across all 96 fits is 0.036 (1.9%) for the real component and 0.038 (3.1%) in the imaginary component.
While it is encouraging that the fitted index matches the expected reasonably well, it is also valuable to consider how well the retrieved parameters actually model the measured spectra.This can help determine how much of the discrepancy between the retrieved and true index of refraction is due to systematic errors in the fit as opposed to systematic errors in the instrument and calibration used to collect the data.A comparison between the fitted S 0 and P spectra and the measured for all three viewing angles is presented in Fig. 2.
The rms error between the measured and retrieved spectra is 205 nW cm 2 srcm −1 in S 0 and 129 nW cm 2 srcm −1 in P. Both of these are well within the noise equivalent spectral radiance (256 nW cm 2 srcm −1 ) of the instrument at this spectral resolution. 19ith the exception of the ozone feature region between 1000 and 1100 cm −1 at 70 deg, the fitted values lie within the error bounds of the measurement at all points.An error in the "known" downwelling may be to blame for this error, but this requires further examination.
A piece of sapphire glass was also examined as part of this experiment.While fused silica and sapphire appear very similar in the visible spectrum, their spectral signatures in the LWIR and material properties are very different.While fused silica is amorphous, the sapphire glass is a more crystalline solid, meaning it is well described by the Lorentz oscillator model.Additionally, sapphire is strongly birefringent.This means that the material has an index of refraction that varies with the polarization state of the light incident on the surface.To account for this, two indices of refraction are solved for, one each for the s-and p-pol radiance.The results are compared with the model for the ordinary and extraordinary rays of sapphire used in the JA Woollam ellipsometer software. 21gain, viewing angles of 30 deg, 50 deg, and 70 deg are used, but the sapphire window is larger than the fused silica wafer, so a 400 pixel window is used.The results of this fit, performed independently on each pixel, are shown in Fig. 3.
Again, this fit is quite accurate, with rms errors comparable to the fused silica wafer, indicating that index of refraction can be retrieved under atmospheric downwelling, even for a c-axis birefringent material.The rms error in the ordinary ray index is 0.028 (2.8%) in the real component and 0.006 (4.5%) in the imaginary component, while the rms error in the extraordinary index is 0.028 (2.7%) and 0.004 (8.2%).The spectral angle between retrieved and expected is 0.79 deg, 2.3 deg, 1.4 deg, and 3.9 deg for n o , κ o , n e , and κ e , respectively.As with the fused silica wafer, the fit is also fairly consistent pixel-to-pixel with an rms standard deviation across all pixels of 0.007 and 0.008 in n and κ for the o-ray and 0.007 and 0.002 for the e-ray.
Figure 4 shows the retrieved S 0 and P compared with the measured values.The rms error between the measured and retrieved S 0 and P is 252 nW cm 2 srcm −1 and 243 nW cm 2 srcm −1 , respectively.Both are close to, but still less than, the NESR of the instrument.It is apparent from Fig. 4 that this discrepancy is primarily due to an unexplained artifact in the measured value between 875 and 925 cm −1 .Note that the 2-sigma pixel-to-pixel variability in this spectral region, denoted by Fig. 2 Measured (blue) compared with retrieved (black) S 0 and P for a fused silica wafer assuming the atmospheric downwelling is known a priori.The shaded regions represent plus/minus two standard deviations across all pixels.
Fig. 1 Retrieved (blue) and true (green) index of refraction for a fused silica wafer when not fitting the atmospheric parameters.Atmospheric downwelling radiance is instead estimated from the gold mirror placed in the scene.The solid blue line represents the median retrieval across all pixels, and the shaded blue region represents plus/minus two standard deviations.The green line is taken from the ellipsometry measurements.
the shaded area, greatly exceeds the NESR of the instrument as well.This could indicate a contaminant existing in varying amounts on different parts of the sapphire window.

Fitting Downwelling Radiance
It is also interesting to see how the fit performs while simultaneously fitting index of refraction and downwelling radiance.Because this testing was done over short ranges, atmospheric transmission and path radiance can be ignored but, in principle, could be solved for in the same way.The results of the fit on the fused silica wafer are shown in Fig. 5.
The rms error between the true and retrieved is 0.146 and 0.178 or 7.7% and 14.3% in n and κ, respectively.The spectral angle difference between the retrieved and ellipsometrymeasured spectra is 2.6 deg and 7.3 deg in n and κ, respectively.As expected, this is slightly larger than when the downwelling radiance is known a priori.It is still within the pixel-to-pixel variability bounds, indicating that this difference is still on the order of the sensor noise, which is encouraging.The pixel-to-pixel variability when fitting both index and downwelling radiance does more than double, however, from 0.036 and 0.038 to 0.073 and 0.087 in n and κ, respectively.
Figure 6 shows the measured and retrieved S 0 and P. Bear in mind that only 24 parameters, 21 fitted and 3 known (the viewing angles), are being used to describe both n and κ, each with 583 measured points, as well as object temperature, and downwelling radiance.Again, however, despite using so few parameters, the spectra can be accurately modeled to well within the noise level of the instrument.This illustrates a way that the often highly correlated nature of hyperspectral information can be exploited to solve for additional parameters describing scene conditions beyond just object temperature.
The RMS error between the fitted and measured S 0 and P is 210 and 150 nW cm 2 srcm −1 , respectively.This is only slightly worse than the 205 and 129 nW cm 2 srcm −1 values obtained with  the known atmosphere, indicating that the fit does reasonably well in describing the atmospheric downwelling.With the exception of between 1050 and 1100 cm −1 , the retrieved values also lie well within the pixel-to-pixel variability margins.One potential source for the discrepancy around 1075 cm −1 could be incorrectly estimating the amount of CO 2 in the atmosphere.When generating the lookup table, the standard MODTRAN5 value of 380 ppmv was used for CO 2 concentration.In reality, the CO 2 concentration is closer to 400 ppmv and may be even higher locally since this test was conducted in a relatively populated area.

Comparison with Temperature-Emissivity Separation
To better put these results into context, a comparison was done against the maximum smoothness TES algorithm. 22To compare, the fitted index of refraction for each pixel was forward modeled to an emissivity using the Fresnel equations and Kirchhoff's law.The true emissivity, as a function of viewing angle, was obtained by forward modeling the ellipsometry measurements in the same way.Because TES requires a downwelling radiance spectrum be provided a priori, it is compared against the fixed atmosphere fit from Sec. 5.1.Additionally, because TES does not utilize polarization, a sensor utilizing TES would likely not have a polarizer.The effect of the polarizer is to reduce the amount of light by a factor of ∼2.To model this polarimetric data, as if there was no polarizer, a two-pixel average was taken when employing the TES algorithm.Figure 7 shows the results of this comparison.The presence of atmospheric features in the residuals indicates that the downwelling estimate was likely not perfect, which could again be due to interpolation errors around sharp atmospheric features after performing the spectral axis correction.Either way, this effect will be present for both techniques.The rms error between retrieved and expected emissivity is 0.0187, 0.0390, and 0.0370 for the TES algorithm at 30 deg, 50 deg, and 70 deg, respectively.For the index retrieval method, the rms error is 0.026, 0.029, and  0.025 at these same three viewing angles.The spectral angle between the retrieved and true emissivities is 0.79 deg, 1.2 deg, and 2.7 deg for the TES algorithm and 1.1 deg, 1.4 deg, and 1.8 deg for the index retrieval.While the two techniques yield similar accuracies in this case, the biggest difference is in the fit-to-fit variability.The rms standard deviation across all retrievals is 0.01 at all three angles for the maximum smoothness TES algorithm while it is less than half that for the index retrieval, 0.004 at all three angles.This illustrates one of the key advantages of estimating index of refraction; because the index is invariant to viewing angle, multiple viewing geometries can be used in conjunction with one another to constrain the optimization.These multiple viewing angles could be obtained either by considering multiple sensors or persistently surveilling a target with a single sensor in motion.When considering this in the context of a low-cost penetrating sensor platform, which is likely to have less than ideal noise characteristics, this ability to exploit several images may significantly improve performance.
The biggest issue, in the view of the authors, with remotely estimating index of refraction is that there is no closed-form expression for index of refraction in terms of radiance.This means that nonlinear search algorithms, which are generally much slower, are necessary.Additionally, there are many local minima in the search space that can be difficult to avoid without sufficient constraining information.This effect has been observed in other datasets related to this work. 17Still, the fact that index of refraction is invariant to viewing geometry, as well as other factors, such as particle size, may overcome these disadvantages in certain scenarios.

Conclusions
This work demonstrates that index of refraction, for specular reflectors, can be accurately retrieved from passive P-HSI radiance measurements.Because index of refraction is invariant with viewing geometry for the vast majority of materials, data from multiple look angles can be used to constrain the fit, effectively producing an averaging effect.The index of refraction was accurately estimated to within 0.14 rms error of 4.5-deg spectral angle error for a fused silica wafer and to within 0.03 rms error and 4-deg spectral angle error for both the ordinary and extraordinary ray index of a sapphire window.Additionally, the ability to simultaneously estimate downwelling radiance without significantly affecting fit performance was demonstrated.Forward modeling these results to an emissivity displays similar accuracy and increased selfconsistency when compared to those estimated via the maximum smoothness TES algorithm.Currently, this research is still in its infancy and is intended to serve as a proof of concept.In the future, we hope to further refine this work by applying it to more realistic test cases, incorporating a polarimetric bidirectional reflectance distribution function model to account for rough surfaces and continuing to explore new ways to model and fit the index.

Fig. 3
Fig. 3 Retrieved (blue) and true (green) index of refraction for both the o-(plot a) and e-ray (plot b) of a sapphire glass window.The solid blue line represents the median retrieval across all pixels, and the shaded blue region represents plus/minus two standard deviations.The truth values were taken from existing models used in the JA Woollam IR-VASE ellipsometry software package for both the o-and e-ray indices of refraction for sapphire.21

Fig. 4
Fig.4Measured (blue) compared with retrieved (black) S 0 and P for a sapphire window when birefringence is accounted for in the fit.The shaded regions represent plus/minus two standard deviation across all pixels.

Fig. 6
Fig. 6 Measured (blue) compared with retrieved (black) S 0 and P for a fused silica wafer while simultaneously fitting atmospheric downwelling.The shaded regions represent plus/minus two standard deviations across all pixels.

Fig. 5
Fig.5Retrieved (blue) and true (green) index of refraction for a fused silica wafer while simultaneously fitting the atmospheric downwelling radiance.The solid blue line represents the median retrieval across all pixels, and the shaded blue region represents plus/minus two standard deviations.The green line is taken from the ellipsometry measurements.