Open Access
24 December 2022 Underdetermined polarimetric measurements for Mueller extrapolations
Author Affiliations +
Abstract

Polarized light-matter interactions are mathematically described by the Mueller matrix (MM)-valued polarized bidirectional reflectance distribution function (pBRDF). A pBRDF is parameterized by 16 degrees of freedom that depend upon scattering geometry. A triple degenerate (TD) MM assumption reduces the degrees of freedom to eight: one for reflectance, six for non-depolarizing properties, and one for depolarization. When the non-depolarizing dominant process is known or assumed (e.g., Fresnel reflection), the degrees of freedom are further reduced to two. For a given material, if the TD model is appropriate and the dominant non-depolarizing process is known, then these two degrees of freedom can be estimated from as few as two polarimetric measurements. Thus, the MM can be extrapolated from a reduced number of measurements. The primary contribution of this work is the development and demonstration of a linear estimator for an MM’s dominant eigenvalue (i.e., single depolarization parameter) that requires fewer measurements than a full MM reconstruction. MM extrapolations from single snapshot acquisitions with a Sony Triton 5.0MP polarization camera are performed at 30 acquisition geometries and two wavelengths on an ensemble of LEGO bricks treated to have varying surface roughness. These extrapolated MMs are compared with MMs reconstructed from a complete dual rotating retarder Mueller imaging polarimeter. The flux error mean and mode are 11.06% and 1.03%, respectively, despite a 10 × reduction in the number of polarimetric measurements.

1.

Introduction

Polarized light-matter interactions are mathematically described by Mueller matrix (MM)-valued functions called polarized bidirectional reflectance distribution functions (pBRDF) corresponding to 16 degrees of freedom at each scattering geometry. Although analytical models exist, building an experimental model with no prior information requires a minimum of 16 polarimetric measurements at each combination of input and output scattering geometries. Using more than 16 polarimetric measurements creates an overdetermined system and a pseudoinverse solution that becomes more robust to noise as more linearly dependent measurements are taken.1,2

Cloude coherency matrix eigenanalysis is a standard technique for analyzing MMs.39 An MM with three equal coherency eigenvalues, a triple degeneracy (TD), has only eight degrees of freedom.10 When the dominant non-depolarizing process is known, this reduces further to two degrees of freedom. These degrees of freedom correspond to the average reflectance and the dominant coherency eigenvalue. Based on empirical observation, the non-depolarizing process is assumed to be a Fresnel reflection, even in the case in which subsurface scattering is stronger.

The objective of this work is to use prior knowledge about a material to extrapolate its MM from a small quantity of measurements. Our strategy is to relate polarimetric measurements to the two degrees of a freedom of a TD MM model when the dominant process is known, rather than to relate measurements to the MM directly. The primary contribution of this work is a linear estimator for an MM’s dominant coherency eigenvalue that requires as few as two polarimetric measurements. To the authors’ knowledge, this is the first method for extrapolating depolarizing MMs with rank four coherency matrices from fewer than 16 measurements. We make use of a single illumination polarization state and a commercial division of a focal plane (DoFP) linear Stokes camera, meaning that four measurements can be performed in a single snapshot acquisition. Experimental results using measurements taken in a snapshot configuration with a Sony Triton 5.0MP polarization camera are presented and compared with a complete dual rotating retarder (DRR) Mueller polarimeter.

2.

Background

The bidirectional reflectance distribution function (BRDF) of a material is a radiometric property defined as the ratio of the differential output radiance to the differential input irradiance.11 The term bidirectional refers to the dependence on both the incident and observation directions. These are specified with the angles θi and ϕi, which are the zenith and azimuth angles of the incident direction, respectively, and θo and ϕo, which are the zenith and azimuth angles of the observation direction, respectively. θi and θo are also referred to as the incident and scattered angles, respectively.

Torrance and Sparrow12 introduced the theory of microfacets to explain off-specular scattering observed in most BRDFs. According to the theory, off-specular reflection from a surface with normal n^ is the result of small, randomly-oriented mirror-like microfacets with varying surface normals h^ combined with a diffuse component. The light that impinges on a microfacet obeys the law of reflection, so its behavior is angle-dependent. For a given pair of incident and exitant directions, ω^i and ω^o, respectively, the amount of reflected light depends on the average of the Fresnel reflectance coefficients, the probability of the microfacet normal satisfying the law of reflection, and the probability of interaction between adjacent microfacets. Light incident on one microfacet may be blocked by a neighbor, preventing full illumination (called shadowing), or light reflected from a microfacet may be blocked by a neighbor from reaching the observer (called masking). The diffuse component is an angle-independent term that describes both light scattered by multiple microfacets before reaching the observer and light that transmits into and then out of the material. The microfacet BRDF model is characterized by the angle-dependent distribution on microfacet orientations, shadowing-masking functions, and the relative weight of specular to diffuse.

A scalar BRDF is a single-valued function of four variables, so a pBRDF is a MM-valued function of four variables. A MM is a 16-element matrix that describes the linear transformation of a polarization state due to light-matter interaction and is given as

Eq. (1)

M=M00m=M00(1m01m02m03m10m11m12m13m20m21m22m23m30m31m32m33).

Having 16 elements means that there are 16 degrees of freedom: one for average reflectance, three for diattenuation, three for retardance, and nine for depolarization.13 Unfortunately, the retardance and depolarization degrees of freedom do not correspond to individual Mueller elements and are instead coupled among multiple Mueller elements. Depolarization is the random modulation of polarization in time, space, and/or wavelength faster than the particular detector in use can resolve.14 The average reflectance is M00, which can be factored out to produce a normalized MM that describes polarimetry separately from radiometry. In this work, an MM designated with a lowercase letter is normalized as in Eq. (1). A challenge inherent to pBRDFs is that polarization properties are defined in the plane transverse to the ray direction. When considering all possible combinations of input and output directions, keeping track of the basis vectors describing the transverse planes becomes paramount.

One of the first models for a pBRDF is an extension of the popular microfacet model by Priest and Germer.15 In this model, the Fresnel equations are not combined to an average term; instead, the Fresnel reflection MM is used with careful consideration for the coordinate systems before and after scattering. This MM has the form

Eq. (2)

M(ω^i,ω^o,n0,n1)=R(αo)FR(n0,n1,θd)R(αi),
where R(α) are Mueller rotation matrices of the form

Eq. (3)

R(α)=(10000cos(2α)sin(2α)00sin(2α)cos(2α)00001),
where αi and αo are the rotations from the polarization state generator (PSG) coordinates to local microfacet coordinates and from local microfacet coordinates to the polarization state analyzer (PSA) coordinates, respectively, and FR(n0,n1,θd) is the Fresnel reflection MM of the form

Eq. (4)

FR(n0,n1,θd)=12(|rp|2+|rs|2|rp|2|rs|200|rp|2|rs|2|rp|2+|rs|200002Re(rp*rs)2Im(rp*rs)002Im(rp*rs)2Re(rp*rs)),
where rs and rp are the well-known Fresnel amplitude reflection coefficients that depend on the refractive indices n0 and n1 and the angle of incidence onto the microfacet θd. It is worth noting that typically αi is considered to be the rotation from the coordinate system defined by n^ and ω^i to local microfacet coordinates, and αo is the rotation from the microfacet coordinates to the coordinate system defined by n^ and ω^o. These rotations are needed for the forward problem of computer graphics rendering from a pBRDF, but because the goal of this work is to recover MMs, the input and output coordinate systems are different.

Other similar pBRDFs based on the microfacet model make use of some form of depolarizing component. The most basic approach is to set the diffuse term to an ideal depolarizer. The ideal depolarizer MM has the form

Eq. (5)

mID=(1000000000000000),
and has the property of converting any incident polarization state to unpolarized light.

Baek et al. proposed another pBRDF that interprets the diffuse term as a depolarizer.16 In this model, the ideal depolarizer is sandwiched between two Fresnel transmission matrices. The Fresnel transmission matrix has a similar form as in Eq. (4), but with the amplitude transmission coefficients replacing the reflection coefficients. The rationale for this diffuse term is that some fraction of the incident light undergoes Fresnel transmission into the material, undergoes subsurface scattering that completely depolarizes the light, and then is transmitted out. A pBRDF model by Kondo et al.17 is a sum of a specular and two diffuse terms.

As stated above, a MM has 16 degrees of freedom and, therefore, must be constrained by at least 16 linearly independent measurements. The most common architecture for a complete MM polarimeter is a DRR polarimeter, which takes a sequence of measurements with different retarder positions.18,19 Taking 16 or more sequential measurements at sufficiently fine scattering geometry sampling to build an empirical pBRDF is a time-intensive process. This motivates the development of snapshot polarimeters. Technologies for snapshot acquisition of MMs include encoding polarization to wavelength,2023 channeled polarimetry using polarization gratings,24,25 structured illumination combined with a DoFP full Stokes camera,26 and splitting the beam to enable multiple DoFP cameras.27 Hybrid approaches also exist that use rotating retarders for illumination but decrease the number of acquisitions using a DoFP camera.28

These technologies are fast, but remain complex. Using partial polarimetry simplifies the measurement requirements. Partial polarimetric systems take fewer than 16 linearly independent measurements, making them underdetermined and incapable of reconstructing the full MM. Recent work by Gonzalez et al.29 demonstrated analysis and decomposition of 3×4 partial MMs measured using four polarized illumination states and a snapshot linear Stokes camera. A partial polarimeter that only makes use of linear illumination and analyzer states can at most reconstruct the upper 3×3 portion of an MM. Swami et al. showed that, for a non-depolarizing MM, symmetry arguments can be applied to the linear partial MM to obtain the full 4×4 matrix.30 Ossikovski and Arteaga showed symmetry arguments for obtaining a full 4×4 MM from 12 elements in which a row or column is missing31 or from nine elements in which a row and column are missing.32 In the 12-element case, it is possible to recover a depolarizing MM but only if it obeys certain symmetry constraints and has only two non-zero coherency eigenvalues.

The pBRDF model explored in this work was derived by Li and Kupinski10 and demonstrated in polarimetric computer graphics renderings by Omer and Kupinski.32 This model is based on the coherency matrix introduced by Cloude.4,5 The coherency matrix and its eigen decomposition (sometimes called the spectral decomposition) are standard techniques for analyzing MMs.39 Using the spectral decomposition, a depolarizing MM is rewritten as a convex sum of up to four non-depolarizing MMs as

Eq. (6)

m=n=03ξnm^n,
where ξn are the Cloude coherency matrix eigenvalues normalized so that n=03ξn=1 and m^n are the non-depolarizing (indicated with the hat ·^) MMs that also have the property that 14n=03m^n=mID. Li and Kupinski showed that, when the smaller three eigenvalues of the Cloude coherency matrix are equal, the MM has the form

Eq. (7)

M=4M003[(ξ014)m^0+(1ξ0)mID],
where m^0 is the dominant non-depolarizing MM and ξ0 now controls the relative weight between this dominant process and the ideal depolarizer. This MM is referred to as being triply degenerate (TD) because the last three eigenvalues are identical. The 16 degrees of freedom of a general MM are reduced to eight: one for M00, one for ξ0, three for the diattenuation orientation and magnitude of m^0, and three for the retardance orientation and magnitude of m^0. Both m^0 and ξ0 are functions of the scattering geometry: θi, ϕi, θo, and ϕo, though this dependence is omitted for brevity. The decomposition of an MM into a non-depolarizing component and an ideal depolarizer has been discussed in the literature, such as in the textbook by Brosseau in which an MM is written as a sum of a non-depolarizing MM that depends on the input Stokes vector and an ideal depolarizer.33 For the TD model, the non-depolarizing MM does not depend on the input Stokes vector and the relative contribution of each component is determined by the unique coherency eigenvalue.

When the dominant process is believed to be Fresnel reflection, the angular dependence of m^0 is known and the degrees of freedom are reduced to two. Based on knowing a priori the dominant coherent process m^0, Li and Kupinski proposed a method to measure ξ0 using only two linearly independent polarimetric measurements.10 Our method would, in principle, work with only two linearly independent polarimetric measurements as well. However, with a commercially available linear Stokes camera, such as the Sony Triton 5.0MP polarization camera used in this work, simultaneous acquisition of four (three linearly independent) polarimetric measurements is possible. The MM extrapolations presented make use of all four measurements. Jarecki and Kupinski34 presented initial results for low albedo measurements. In this work, high albedo measurements and additional analysis are included.

3.

Methods

3.1.

Polarimetry

The MM of a particular light-matter interaction is measured by illuminating the sample with a known PSG and then measuring the scattered light through a known PSA as

Eq. (8)

pj=ajTMgj=wjM16×1,
where M is the object MM, gj is the j’th PSG state, aj is the j’th PSA state, and pj is the j’th irradiance measurement at the detector. This measurement equation is rewritten in the right side of Eq. (8) as an inner product of two vectors, where wj=ajgjT and M16×1 is the object MM elements ordered in a vector form. Conventional full polarimetry requires that at least J=16 linearly-independent measurements be performed with different illumination and analyzer states to constrain the 16 degrees of freedom of M. Most polarimeters have J>16 to create an overdetermined system for mitigating the effects of noise.1,2 The MM estimate is then reconstructed via the pseudoinverse as

Eq. (9)

M˜16×1=W+P,
where W is called the polarimetric measurement matrix and has rows wj, P is the vector of flux measurements with elements pj, and the tilde on M˜ indicates the quantity is an estimate.

3.2.

Triple Degenerate Model

Using the Cloude coherency eigenanalysis of MMs, Li and Kupinski derived the TD model in which an MM consists of a weighted sum of a dominant non-depolarizing MM and an ideal depolarizer, and the relative weights are controlled by a single depolarization parameter as shown in Eq. (7).10 Under the assumptions of the TD model, the degrees of freedom in an MM are reduced from 16 to eight: one for M00, three for the diattenuation and three for the retardance of dominant non-depolarizing MM m^0, and one for the dominant coherency eigenvalue ξ0 that controls the relative weights in the model.

The partial polarimetric method presented in this work relies on a priori knowledge of the dominant process m^0 and how it varies over the measured scattering geometries. Based on empirical observations, Fresnel reflection due to microfacets is used as the dominant process. The coordinate system transformations needed to calculate Fresnel reflection over the field of view of a measurement are described in Appendix A.

For a TD MM, the contribution of mID outweighs the contribution of m^0 when ξ0<0.625. In this regime, subsurface scattering dominates the interaction, but the single largest coherent process can still be Fresnel reflection.

3.3.

Mueller Matrix Extrapolation

A noise-free model for flux measurements P of a TD MM are written as a linear system

Eq. (10)

P=ΦTα.

Here, Φ is a matrix with rows that are the measurement matrix W applied to the dominant process m^0 and ideal depolarizer mID from the TD model as

Eq. (11)

Φ=(p0TpIDT)=([Wm^0]T[WmID]T),
and the elements of α are the weights in the TD model as

Eq. (12)

α=(α0αID)=4M003(ξ0141ξ0).

It is worth reiterating that a benefit of the TD model is that the relative weights are controlled by a single depolarization parameter ξ0 rather than varying independently. An estimate of the coefficients α˜ is recovered with the Moore–Penrose pseudoinverse of ΦT as

Eq. (13)

α˜=[ΦT]+P=(Wm^0WmID)+P,
where P here is a vector of noisy flux measurements. Solving the system in Eq. (12) for the model parameters, estimates for ξ˜0 and M˜00 are

Eq. (14)

ξ˜0=14+α˜0/α˜ID1+α˜0/α˜ID,M˜00=3α˜04ξ˜01,
where α˜0 and α˜ID are the elements of α˜. ξ˜0 is the parameter of interest because it determines the fractional contributions of the dominant coherent process and the ideal depolarizer. This fractional contribution adjusts the depolarization of the MM that changes with scattering geometry, albedo, and surface texture.35 The estimated average reflectance M˜00 is ignored in the results because the polarimeters used were not calibrated to produce data in absolute radiometric units. Instead, the MM results will be normalized.

3.4.

Polarimeters

The ground-truth MM images were taken using a full Mueller imaging polarimeter called the RGB950.36 The RGB950 [shown in Fig. 1(a)] is a DRR polarimeter that reconstructs an MM image from a sequence of 40 polarimetric measurements at different retarder positions. Data were taken at two wavelengths: 662±11.17  nm (red) and 451±9.83  nm (blue). Thirty scattering geometries shown in Table 1 were measured using a rotation stage for the sample and a goniometric swing arm for the camera.

Fig. 1

(a) The full Mueller imaging polarimeter referred to as the RGB950. It takes 40 polarimetric measurements. (b) The goniometer and PSG of the RGB950 but with the PSA replaced by a commercial off-the-shelf linear Stokes camera. The Stokes camera takes 4 polarimetric measurements in a single snapshot with a static PSG state.

OE_61_12_123104_f001.png

Table 1

A total of 30 acquisition geometries specified on-axis where ϕi and ϕo are 0 deg. For each angle between the sample surface normal and source, θi, measurements are performed for six angles between the surface normal and the camera, θo. The scattering geometries across the field of view of an image will have θi, ϕi, θo, and ϕo that deviate from these on-axis values. Some acquisition geometries are omitted from analysis because exposure issues with the linear Stokes camera produced non-physical MM extrapolations.

θiθo,1θo,2θo,3θo,4θo,5θo,6
−10 deg10 deg20 deg30 deg40 deg50 deg60 deg
−25 deg15 deg25 deg35 deg45 deg55 deg65 deg
−40 deg20 deg30 deg40 deg50 deg60 deg70 deg
−55 deg25 deg35 deg45 deg55 deg65 deg75 deg
−70 deg30 deg40 deg50 deg60 deg70 deg80 deg

The linear partial polarimetric experiment was performed using a Sony Triton 5.0MP polarization camera, shown in Fig. 1(b). This camera has an array of micropolarizers in front of the detector elements, so four polarimetric measurements (three of which are linearly independent) are taken simultaneously at the cost of spatial resolution. The rank three measurement matrix of this system is underdetermined for full Mueller polarimetry but is overdetermined for recovering the two unknown degrees of freedom in the TD model. Designing polarimeters to be overdetermined systems is a well-established practice in the literature.1,2,37 The PSG of the RGB950 was used to keep the illumination consistent, but only horizontal linear polarization was used with the Stokes camera measurements. In principle, unpolarized illumination would also work. In this work, the calibration and operation software provided by the manufacturer are used. This could lead to exaggerated errors in the extrapolation results.

3.5.

Samples

The samples measured in this work are a collection of blue LEGO bricks shown in Fig. 2. These are the same LEGOs used by Li and Kupinski and are a collection of individual objects with the same material properties and albedo but with varying textures.10 The surface roughness of each brick was measured using a white light interferometer. The roughness averages (Ra) in microns for each brick are top: 0.49, 0.56, 3.35, middle: 3.55, 2.62, 0.35, and bottom: 1.68, 1.26, 6.32. Because the bricks are blue, the different wavelengths represent different albedo cases: 662 nm illumination is low albedo and 451 nm illumination is high albedo. Umov’s effect states that the amount of depolarization is expected to trend positively with albedo, so these albedo cases also represent cases with different amounts of depolarization.38

Fig. 2

A tower of blue plastic LEGO bricks, each being sanded with a different grit of sandpaper. This represents a group of objects with similar properties and albedos for a given wavelength but with different surface textures. The Ra in microns for each brick are top: 0.49, 0.56, 3.35, middle: 3.55, 2.62, 0.35, and bottom: 1.68, 1.26, 6.32.

OE_61_12_123104_f002.png

Fig. 3

Comparison of the MM image results at 451 nm of (a) full reconstruction with 40 polarimetric measurements to (b) MM image results of partial polarimetric extrapolation (Video 1) and at 662 nm of (c) full reconstruction and (d) extrapolation (Video 2). The still images are for the geometry θi=25  deg, θo=25  deg (Video 1, MP4, 3 MB [URL: https://doi.org/10.1117/1.OE.61.12.123104.s1]; Video 2, MP4, 4 MB [URL: https://doi.org/10.1117/1.OE.61.12.123104.s2).

OE_61_12_123104_f003.png

4.

Results

Figure 3(a) shows a MM image at 451 nm calculated from conventional Mueller polarimetry with 40 measurements compared to Fig. 3(b) which shows an MM image extrapolated from 4 measurements. Figures 3(c) and 3(d) show this same comparison at 662 nm. Figures 3(a) and 3(b) are a still from Video 1 and Figs. 3(c) and 3(d) are a still from Video 2. Videos 1 and 2 compare full MM reconstructions to extrapolations at multiple geometries. Depolarization can be qualitatively observed by comparing the magnitude of M00 with other matrix elements: regions of the image where all matrix elements have a smaller magnitude than M00 have larger depolarization. Depolarization is expected to be stronger for the high albedo case of blue bricks under blue illumination, and this can be seen by the relatively lower magnitudes across the field of view in both the reconstruction and extrapolation. Likewise, the expectation of lower depolarization is met for the low albedo case of the blue bricks under red illumination. The trend of increased depolarization with surface roughness is also captured by the extrapolation.

4.1.

Error in Dominant Eigenvalue Estimate

Figure 4 shows the estimated values of ξ˜0 for the smoothest and roughest textured bricks. Larger values correspond to a larger estimated contribution of the dominant process, or equivalently lower depolarization. For both brick textures, estimates of ξ˜0 are larger in the low albedo case of 662 nm illumination than in the high albedo case of 451 nm illumination. This is in agreement with expectations from Umov’s effect in which depolarization trends positively with albedo. Furthermore, for both wavelengths, estimates of ξ˜0 are larger for the smooth brick and smaller for the rough brick. This trend matches the expectation of a rougher texture resulting in higher depolarization.

Fig. 4

The estimate ξ˜0 calculated with Eq. (14) versus acquisition geometry at 451 nm for (a) the smoothest brick and (b) the roughest brick, and at 662 nm for (c) the smoothest brick and (d) the roughest brick. Geometries at which the dynamic range of the linear Stokes camera caused non-physical MM extrapolations are omitted.

OE_61_12_123104_f004.png

Figure 5 shows the difference between the true value of ξ0 from the complete MM reconstruction and the estimated ξ˜0 from the linear Stokes measurements. Positive-valued differences correspond to an overestimation of the dominant non-depolarizing process or equivalently an underestimation of the amount of depolarization.

Fig. 5

The estimate ξ˜0 calculated with Eq. (14) minus the true ξ0 versus acquisition geometry at 451 nm for (a) the smoothest brick and (b) the roughest brick, and at 662 nm for (c) the smoothest brick and (d) the roughest brick. When the difference is positive valued, the contribution of the dominant non-depolarizing term is overestimated. Geometries at which the dynamic range of the linear Stokes camera caused non-physical MM extrapolations are omitted.

OE_61_12_123104_f005.png

Fig. 6

Flux vectors at 451 nm of (a) the smoothest brick and (b) the roughest brick (Video 3) and at 662 nm of (c) smoothest brick and (d) the roughest brick (Video 4). Flux vectors are calculated by averaging the normalized MM on a 2×2  pixel region of interest and then applying W40 to simulate what the RGB950 would measure. The error bars are ±1 standard deviation in the region of interest. The flux vectors shown here are for the full MM reconstructions (blue), the nearest TD approximation of the reconstruction (orange), and the extrapolated MMs (green) at θi=25deg, θo=25deg. Geometries at which the dynamic range of the linear Stokes camera caused non-physical MM extrapolations are omitted from the videos (Video 3, MP4, 9 MB [URL: https://doi.org/10.1117/1.OE.61.12.123104.s3; Video 4, MP4, 9 MB [URL: https://doi.org/10.1117/1.OE.61.12.123104.s4).

OE_61_12_123104_f006.png

Extrapolations at 451 nm tended to overestimate ξ˜0 at more geometries than at 662 nm. 662 nm is the low-albedo case, where Umov’s effect indicates that depolarization is lower, so it is possible that the method is most successful for low-depolarization cases.

4.2.

Simulated Flux Vectors

To compare the MMs, the measurement matrix of the RGB950 W40 is applied to a 2×2  pixel average of the normalized extrapolated MM and the full reconstructed MM to simulate the flux measurements that the RGB950 would take. Additionally, the measurement noise is indicated by the standard deviation error bar on each flux measurement. The resulting flux vectors for two bricks multiple geometries are shown in Video 3 at 451 nm and in Video 4 at 662 nm. Figures 6(a) and 6(b) are a still from Video 3 and Figs. 6(c) and 6(d) are a still from Video 4. Also shown in Fig. 6 is the nearest TD approximation of the reconstructed MM. This is calculated by setting the three smallest coherency eigenvalues to 1ξ03. This TD approximation has the exact correct dominant process m^0 and represents the best possible extrapolation based on a TD model. Flux vectors in which the TD approximation and extrapolation show a similar deviation from the full reconstruction could indicate that the TD assumption is not valid. However, Videos 3 and 4 do not consistently show disagreement between the reconstruction and TD approximation where the extrapolation deviates from the reconstruction. The measurements at 662 nm (the low albedo case) exhibit larger error bars for both the smooth and rough brick compared with measurements at 451 nm (high albedo). This matches the expectation because, despite larger polarization modulation for low albedo per Umov’s effect, the overall amount of light is lower.38 The largest realizations of measurement noise occur for the rough brick at scattering geometries near those with non-physical results but are not yet themselves non-physical. However, for all other measurements, the disagreement between the extrapolations and the full reconstruction is larger than the error bars. This means that errors in the extrapolation are more likely the result of discrepancies in the assumed dominant process.

To compare MMs with a single-valued metric, the flux error ϵ is defined as

Eq. (15)

ϵ=j40|pjp˜j|j40pj,
where pj are the elements of the flux vector simulated by applying the RGB950 measurement matrix W40 to the ground truth MM (i.e., the full reconstruction) and p˜j are the elements of the flux vector simulated by applying W40 to the MM being tested. This can be interpreted as adding up all of the discrepancies and normalizing by the total expected flux. A flux error ϵ=0 would mean that the two MMs yield the same RGB950 measurements. This physical interpretation is the motivation for choosing Eq. (15) as our figure of merit instead of a sum of squared differences between two MMs. Furthermore, small disagreements in multiple off-diagonal MM elements could yield a small squared difference in MM elements but be an appreciable retardance difference.

Figure 7 shows histograms of flux errors calculated from the same flux vectors as in Fig. 6 but also includes the other textures. Figure 7(a) is the histogram of flux errors between the full reconstruction MMs and the MMs extrapolated from linear Stokes images using an assumed Fresnel reflection dominant process. The sources of error are measurement noise, the assumed dominant process, and the assumption of a TD eigenspectrum. The mean is 11.65%, and the mode is 1.03%. Figure 7(b) shows the histogram of flux errors between the full reconstruction MMs and those same MMs truncated to have a TD eigenspectrum. The process of TD truncation, explained above, preserves the exact dominant process and is not a new noise-realization, so the only source of error is the difference in eigenspectrum. The mean is 3.61%, and the mode is 0.54%. The narrower distribution in Fig. 7(b) compared with that in Fig. 7(a) indicates that the assumption of a TD eigenspectrum is not the largest source of error.

Fig. 7

Histograms of flux error, as defined in Eq. (15), of (a) the extrapolated MMs based on assuming the dominant process is Fresnel reflection and (b) the nearest TD truncations of the true MMs relative to the true MMs. Each histogram contains the data over all textures, all geometries, and both wavelengths. The narrower distribution in (b) compared with (a) indicates that the assumption of a TD eigenspectrum is not the largest source of error. The mean of (a) is 11.06%, and the mode is 1.03%. The mean of (b) is 3.61%, and the mode is 0.54%.

OE_61_12_123104_f007.png

For both bricks and both wavelengths, Fig. 8 shows that the larger flux errors tend to occur for the larger incident and scattering angles. The maximum flux error is 0.42 and occurs for the rough brick under 662-nm illumination at θi=60  deg, θo=65  deg, despite the maximum error in ξ0 occurring for the smooth brick at 451 nm at θi=40  deg, θo=60  deg.

Fig. 8

Flux vector error ϵ, as defined in Eq. (15), of the LEGO brick MMs extrapolated from linear Stokes images and MMs reconstructed by the RGB950 plotted versus acquisition geometry for (a) and (c) the smoothest brick and (b) and (e) the roughest brick. (a) and (b) The high albedo case and (c) and (d) the low albedo case. Each flux vector was calculated by averaging the normalized MM on a 2×2  pixel region and then applying W40 to simulate what the RGB950 would measure. Each point corresponds to an acquisition geometry according to Table 1. Geometries at which the dynamic range of the linear Stokes camera caused non-physical MM extrapolations are omitted.

OE_61_12_123104_f008.png

Table 2 gives the flux error ϵ averaged over the acquisition geometry for each brick and at both wavelengths. Both the overall minimum and maximum flux errors occur for 662 nm illumination on the 0.56 and 3.35  μm Ra bricks, respectively. Because these extrema do not correspond to the smoothest or roughest textures, it is likely that texture is not the dominant contributing factor to the error. Averaging over texture, the high-albedo case has an error of 10.50% and the low-albedo case has an error of 11.65%.

Table 2

The flux error ϵ, as defined in Eq. (15), averaged over acquisition geometry for each brick. Acquisition geometries that produced non-physical MM extrapolations due to the dynamic range of the linear Stokes camera are omitted.

Brick Ra (μm)0.350.490.561.261.682.623.353.556.32
ϵ at 451 nm (%)11.5011.408.968.337.7911.0613.2812.289.93
ϵ at 662 nm (%)10.428.427.2312.828.7512.9515.5515.1614.89

5.

Conclusion

For materials described by the simple triple-degenerate polarized light scattering model, this work proposes a new and simplified way to measure the MM. Typical Mueller polarimetry requires 16 or more polarimetric measurements to reconstruct an MM, whereas the TD MM model allows for extrapolation from as few as two measurements when the dominant process is known a priori. Existing methods for extrapolating full MMs from partial polarimetry require non-depolarizing MMs30,31 or, at most, a rank-two coherency matrix.39 To the authors’ knowledge, this work is the only method for extrapolating full-rank depolarizing MMs. Additionally, this method is compatible with existing DoFP polarimeter technology and, therefore, can be made into a snapshot polarimeter.

To demonstrate the method, extrapolations at different geometries of LEGO bricks with varying roughness are performed with a commercial linear Stokes camera and compared with the full MM polarimeter reconstructions. The depolarization, which varies with surface roughness, is apparent even on visual inspection of the diagonal elements of the extrapolated MMs. Over varying textures, geometries, and albedos, the partial polarimetric extrapolations achieve flux error mean and mode of 11.06% and 1.03%, respectively, despite a 10× reduction in the number of polarimetric measurements.

6.

Appendix A: Calculation of Rotated Fresnel Reflection Matrix

From the vectors ω^i and ω^o, the Fresnel reflection matrix at a given scattering geometry can be found using conventions from polarization ray trace calculus and then converting to an MM formalism.

6.1.

Unrotated Fresnel Reflection

Based on the microfacet assumption, there are subresolution facets that cause specular reflection from ω^i to ω^o. The surface normal of such a microfacet is called the halfway vector h^ and is calculated as

Eq. (16)

h^=ω^oω^i|ω^oω^i|.

The angle of incidence onto the facet is called the difference angle θd and is calculated as

Eq. (17)

θd=arccos(h^·(ω^i)).

From the angle of incidence, the Fresnel reflection coefficients are

Eq. (18)

rs(θd)=cos(θd)n1sin2(θd)n2cos(θd)+n1sin2(θd)n2andrp(θd)=ncos(θd)1sin2(θd)n2ncos(θd)+1sin2(θd)n2,
where n is an estimate of the index of refraction for the measured material and the index of the incident material is assumed to be 1. The MM for Fresnel reflection is given in Eq. (4).

6.2.

Coordinate System of Microfacet

The input eigenpolarization basis Σ^iΠ^i for the microfacet with normal h^ in the transverse plane to ω^i is calculated as

Eq. (19)

Σ^i=h^×ω^iandΠ^i=ω^i×Σ^i.

The output eigenpolarization basis Σ^oΠ^o for the microfacet with normal h^ in the plane transverse to ω^o is calculated as

Eq. (20)

Σ^o=h^×ω^oandΠ^o=ω^o×Σ^o.

Greek letters are used to differentiate the microfacet s-p basis coordinates from the macrosurface s-p basis coordinates conventionally used in polarization ray tracing.

6.3.

Coordinate Systems of PSG and PSA

When light is reflected off a surface, the polarization states also undergo a geometric transform independent of the material properties. The polarization state is initially parameterized in the coordinate system of the PSG and ends in the coordinate system of the PSA. Using a point source illumination model and a pinhole camera model, the polarization states are best described using a double pole basis.

Local basis vectors in a double-pole system are determined by an antipole direction, a rotation matrix from the antipole to the propagation direction, and a reference basis vector. With the object centered at the origin, the antipole directions of the input and output, a^i and a^o, respectively, are calculated as

Eq. (21)

a^i=s|s|anda^o=c|c|,
where s is the coordinate of the point source and c is the coordinate of the camera pinhole. A 3×3 rotation matrix by angle ϕ about the axis r^=(rx,ry,rz) is calculated as

Eq. (22)

R3(r^,ϕ)=[rx2(1cosϕ)+cosϕrxry(1cosϕ)rzsinϕrxrz(1cosϕ)+rysinϕryrx(1cosϕ)+rzsinϕry2(1cosϕ)+cosϕryrz(1cosϕ)rxsinϕrzrx(1cosϕ)rysinϕrzry(1cosϕ)+rxsinϕrz2(1cosϕ)+cosϕ.].

The basis vectors for the input direction ω^i are calculated using rotation axis r^i and rotation angle ϕi

Eq. (23)

r^i=ω^i×a^iandϕi=arccos(ω^i·a^i),
where basis vectors for the outward direction ω^o are calculated using rotation axis r^o and rotation angle ϕo

Eq. (24)

r^o=ω^o×a^oandϕo=arccos(ω^o·a^o).

Assuming the global vertical is y^=(0,1,0), then the local vertical polarization in the PSG and PSA coordinate systems, y^i,loc and y^o,loc, respectively, are calculated as

Eq. (25)

y^i,loc=R3(r^i,θi)y^andy^o,loc=R3(r^o,θo)y^.

6.4.

Rotated Fresnel Mueller Matrix

The rotation matrix in Mueller calculus has the form

Eq. (26)

R(α)=[10000cos2αsin2α00sin2αcos2α00001].

The angle of rotation in the transverse plane αi from the coordinates to the microfacet basis is calculated as

Eq. (27)

αi=2arccos(y^i,loc·Π^i)sign((y^i,loc×Π^i)·ω^i),
and the angle of rotation in the transverse plane from the microfacet basis to the PSA coordinates is calculated as

Eq. (28)

αo=2arccos(Π^o·y^o,loc)sign((Π^o×y^o,loc)·ω^o).

The rotated Fresnel reflection MM is then calculated as

Eq. (29)

M=R(αo)FR(θd)R(αi).

Acknowledgment

The authors declare that they have no conflicts of interest.

Code, Data, and Materials Availability

Mueller matrix measurement data from the RGB950 and linear Stokes images from the Sony Triton 5.0MP polarization camera can be found at https://doi.org/10.25422/azu.data.19763086.

References

1. 

M. H. Smith, “Optimization of a dual-rotating-retarder Mueller matrix polarimeter,” Appl. Opt., 41 2488 –2493 https://doi.org/10.1364/AO.41.002488 APOPAI 0003-6935 (2002). Google Scholar

2. 

K. M. Twietmeyer and R. A. Chipman, “Optimization of Mueller matrix polarimeters in the presence of error sources,” Opt. Express, 16 11589 –11603 https://doi.org/10.1364/OE.16.011589 OPEXFF 1094-4087 (2008). Google Scholar

3. 

S. R. Cloude, “Conditions for the physical realisability of matrix operators in polarimetry,” Proc. SPIE, 1166 177 –187 https://doi.org/10.1117/12.962889 PSISDG 0277-786X (1990). Google Scholar

4. 

S. R. Cloude and E. Pottier, “Concept of polarization entropy in optical scattering,” Opt. Eng., 34 (6), 1599 –1610 https://doi.org/10.1117/12.202062 (1995). Google Scholar

5. 

S. R. Cloude, “Group theory and polarisation algebra,” Optik, 75 26 –36 OTIKAJ 0030-4026 (1986). Google Scholar

6. 

Y.-Q. Jin and S. Cloude, “Numerical eigenanalysis of the coherency matrix for a layer of random nonsphericalscatterers,” IEEE Trans. Geosci. Remote Sens., 32 (6), 1179 –1185 https://doi.org/10.1109/36.338366 IGRSD2 0196-2892 (1994). Google Scholar

7. 

M. Kupinski et al., “Polarimetric measurement utility for pre-cancer detection from uterine cervix specimens,” Biomed. Opt. Express, 9 5691 –5702 https://doi.org/10.1364/BOE.9.005691 BOEICL 2156-7085 (2018). Google Scholar

8. 

C. J. R. Sheppard, “Parameterization of the Mueller matrix,” J. Opt. Soc. Am. A, 33 2323 –2332 https://doi.org/10.1364/JOSAA.33.002323 JOAOD6 0740-3232 (2016). Google Scholar

9. 

C. J. R. Sheppard, A. L. Gratiet and A. Diaspro, “Factorization of the coherency matrix of polarization optics,” J. Opt. Soc. Am. A, 35 586 –590 https://doi.org/10.1364/JOSAA.35.000586 JOAOD6 0740-3232 (2018). Google Scholar

10. 

L. Li and M. Kupinski, “Merit functions and measurement schemes for single parameter depolarization models,” Opt. Express, 29 18382 –18407 https://doi.org/10.1364/OE.425295 OPEXFF 1094-4087 (2021). Google Scholar

11. 

F. O. Bartell, E. L. Dereniak and W. L. Wolfe, “The theory and measurement of bidirectional reflectance distribution function (BRDF) and bidirectional transmittance distribution function (BTDF),” Proc. SPIE, 0257 154 –160 https://doi.org/10.1117/12.959611 PSISDG 0277-786X (1981). Google Scholar

12. 

K. E. Torrance and E. M. Sparrow, “Theory for off-specular reflection from roughened surfaces*,” J. Opt. Soc. Am., 57 1105 –1114 https://doi.org/10.1364/JOSA.57.001105 JOSAAH 0030-3941 (1967). Google Scholar

13. 

R. Chipman, W. Lam and G. Young, Polarized Light and Optical Systems, CRC Press( (2018). Google Scholar

14. 

B. H. Billings, “A monochromatic depolarizer*†,” J. Opt. Soc. Am., 41 966 –975 https://doi.org/10.1364/JOSA.41.000966 JOSAAH 0030-3941 (1951). Google Scholar

15. 

R. Priest and T. Germer, “Polarimetric BRDF in the microfacet model: Theory and measurements,” in Proc. Military Sens. Symp. (MSS) Specialty Group Meeting on Passive Sens., (2000). Google Scholar

16. 

S.-H. Baek et al., “Simultaneous acquisition of polarimetric SVBRDF and normals,” ACM Trans. Graph., 37 1 –15 https://doi.org/10.1145/3272127.3275018 ATGRDF 0730-0301 (2018). Google Scholar

17. 

Y. Kondo et al., “Accurate polarimetric BRDF for real polarization scene rendering,” Lect. Notes Comput. Sci., 12364 220 –236 https://doi.org/10.1007/978-3-030-58529-7_14 LNCSD9 0302-9743 (2020). Google Scholar

18. 

R. M. A. Azzam, “Photopolarimetric measurement of the Mueller matrix by fourier analysis of a single detected signal,” Opt. Lett., 2 148 https://doi.org/10.1364/OL.2.000148 OPLEDP 0146-9592 (1978). Google Scholar

19. 

D. H. Goldstein, “Mueller matrix dual-rotating retarder polarimeter,” Appl. Opt., 31 6676 –6683 https://doi.org/10.1364/AO.31.006676 APOPAI 0003-6935 (1992). Google Scholar

20. 

N. Hagen, K. Oka and E. L. Dereniak, “Snapshot Mueller matrix spectropolarimeter,” Opt. Lett., 32 2100 https://doi.org/10.1364/OL.32.002100 OPLEDP 0146-9592 (2007). Google Scholar

21. 

M. Dubreuil et al., “Snapshot Mueller matrix polarimeter by wavelength polarization coding,” Opt. Express, 15 13660 –13668 https://doi.org/10.1364/OE.15.013660 OPEXFF 1094-4087 (2007). Google Scholar

22. 

P. Lemaillet, S. Rivet and B. L. Jeune, “Optimization of a snapshot Mueller matrix polarimeter,” Opt. Lett., 33 144 –146 https://doi.org/10.1364/OL.33.000144 OPLEDP 0146-9592 (2008). Google Scholar

23. 

A. S. Alenin and J. S. Tyo, “Task-specific snapshot Mueller matrix channeled spectropolarimeter optimization,” Proc. SPIE, 8364 836402 https://doi.org/10.1117/12.921911 PSISDG 0277-786X (2012). Google Scholar

24. 

M. W. Kudenov et al., “Snapshot imaging Mueller matrix polarimeter using polarization gratings,” Opt. Lett., 37 1367 –1369 https://doi.org/10.1364/OL.37.001367 OPLEDP 0146-9592 (2012). Google Scholar

25. 

M. W. Kudenov et al., “Snapshot imaging Mueller matrix instrument,” Proc. SPIE, 8897 88970S https://doi.org/10.1117/12.2028546 PSISDG 0277-786X (2013). Google Scholar

26. 

A. Zaidi et al., “Towards compact and snapshot channeled Mueller matrix imaging,” Opt. Lett., 47 722 https://doi.org/10.1364/OL.446755 OPLEDP 0146-9592 (2022). Google Scholar

27. 

T. Huang et al., “Fast Mueller matrix microscope based on dual DoFP polarimeters,” Opt. Lett., 46 1676 –1679 https://doi.org/10.1364/OL.421394 OPLEDP 0146-9592 (2021). Google Scholar

28. 

I. J. Vaughn et al., “A portable imaging Mueller matrix polarimeter based on a spatio-temporal modulation approach: theory and implementation,” Proc. SPIE, 9613 961312 https://doi.org/10.1117/12.2188675 PSISDG 0277-786X (2015). Google Scholar

29. 

M. Gonzalez et al., “Introduction of a 3 4 Mueller matrix decomposition method,” J. Phys. D Appl. Phys., 54 424005 https://doi.org/10.1088/1361-6463/ac1622 (2021). Google Scholar

30. 

M. Swami, H. Patel and P. Gupta, “Conversion of 3 × 3 Mueller matrix to 4 × 4 Mueller matrix for non-depolarizing samples,” Opt. Commun., 286 18 –22 https://doi.org/10.1016/j.optcom.2012.08.094 OPCOB8 0030-4018 (2013). Google Scholar

31. 

R. Ossikovski and O. Arteaga, “Complete Mueller matrix from a partial polarimetry experiment: the nine-element case,” J. Opt. Soc. Am. A, 36 403 –415 https://doi.org/10.1364/JOSAA.36.000403 JOAOD6 0740-3232 (2019). Google Scholar

32. 

K. Omer and M. Kupinski, “Compression, interpolation, and importance sampling for polarized brdf models,” Opt. Express, 30 25734 –25752 https://doi.org/10.1364/OE.455126 OPEXFF 1094-4087 (2022). Google Scholar

33. 

C. Brosseau, Fundamentals of Polarized Light: A Statistical Optics Approach, Wiley, New York (1998). Google Scholar

34. 

Q. Jarecki and M. Kupinski, “Extrapolating Mueller matrices from linear Stokes images,” Proc. SPIE, 12112 77 –89 https://doi.org/10.1117/12.2619055 PSISDG 0277-786X (2022). Google Scholar

35. 

M. K. Kupinski et al., “Angle of linear polarization images of outdoor scenes,” Opt. Eng., 58 (8), 082419 https://doi.org/10.1117/1.OE.58.8.082419 (2019). Google Scholar

36. 

J. M. López-Téllez et al., “Broadband extended source imaging Mueller-matrix polarimeter,” Opt. Lett., 44 1544 –1547 https://doi.org/10.1364/OL.44.001544 OPLEDP 0146-9592 (2019). Google Scholar

37. 

I. J. Vaughn and B. G. Hoover, “Noise reduction in a laser polarimeter based on discrete waveplate rotations,” Opt. Express, 16 2091 –2108 https://doi.org/10.1364/OE.16.002091 OPEXFF 1094-4087 (2008). Google Scholar

38. 

V. N. Umow, “Chromatische depolarisation durch lichtzerstreuung,” Physikalische Zeitschrift, 6 (20), 674 –676 https://doi.org/10.1134/S1024856021060269 PHZTAO 0369-982X (1905). Google Scholar

39. 

O. Arteaga and R. Ossikovski, “Complete Mueller matrix from a partial polarimetry experiment: the 12-element case,” J. Opt. Soc. Am. A, 36 416 –427 https://doi.org/10.1364/JOSAA.36.000416 JOAOD6 0740-3232 (2019). Google Scholar

Biography

Quinn Jarecki is a third-year PhD candidate at the Wyant College of Optical Sciences at the University of Arizona. He received his BS degree in optical science and BM in music performance from the University of Arizona in 2019. His research interests include full and partial Mueller polarimetric imaging and polarized light scattering.

Meredith Kupinski is an assistant professor at the Wyant College of Optical Sciences at the University of Arizona. She received her MS degree and PhD in optical science from the University of Arizona in 2003 and 2008, respectively. Her research interests include task-relevant metrics for imaging system design, estimation/detection theory, and stochastic system analysis and information quantification.

CC BY: © The Authors. Published by SPIE under a Creative Commons Attribution 4.0 International License. Distribution or reproduction of this work in whole or in part requires full attribution of the original publication, including its DOI.
Quinn Jarecki and Meredith K. Kupinski "Underdetermined polarimetric measurements for Mueller extrapolations," Optical Engineering 61(12), 123104 (24 December 2022). https://doi.org/10.1117/1.OE.61.12.123104
Received: 24 June 2022; Accepted: 30 November 2022; Published: 24 December 2022
Lens.org Logo
CITATIONS
Cited by 2 scholarly publications.
Advertisement
Advertisement
KEYWORDS
Polarimetry

Depolarization

Cameras

Polarization

Reflection

Light sources and illumination

Matrices

Back to Top