Diagnostics from Polarization of Scattered Optical Light from Galactic Infrared Cirrus

We propose polarization of scattered optical light from intermediate Galactic latitude infrared cirrus as a new diagnostic to constrain models of interstellar dust and the anisotropic interstellar radiation field (aISRF). For single scattering by a sphere, with Mie scattering phase functions for intensity and polarized intensity for a dust model at a given wavelength (Sloan r and g bands), and with models of anisotropic illumination from the entire sky (represented in HEALPix), we develop the formalism for calculating useful summary parameters for an integrated flux nebula (IFN): the average of the phase function weighted by the illumination, polarization angle (ψ), and polarization fraction (p). To demonstrate the diagnostic discrimination of polarization from scattered light, we report on the effects of different anisotropic illumination models and different dust models on the summary parameters for the Spider IFN. The summary parameters are also sensitive to the IFN location, as we illustrate using FRaNKIE illumination models. For assessing the viability of dust and aISRF models, we find that observations of ψ and p of scattered light are indeed powerful new diagnostics to complement joint modeling of the intensity of scattered light (related to the average phase function) and the intensity of thermal dust emission. However, optically thin IFNs that can be modeled using single scattering are faint and p is not large, as it could be with Rayleigh scattering, and so these observations need to be carried out with care and precision. Results for the Draco nebula compared to the Spider illustrate the challenge.


INTRODUCTION
Polarization of scattered optical light from intermediate Galactic latitude infrared cirrus is investigated as a new diagnostic to further constrain models of interstellar dust and the anisotropic interstellar radiation field (aISRF).Sandage (1976) identified the nature of these faint optical diffuse cloud structures, which are called integrated flux nebulae (IFN) because they are illuminated not by the light of single star but by the ISRF of the Galaxy.Zhang et al. (2023) demonstrated the ben-efits of joint modeling of the intensity of scattered light and thermal emission from dust, focusing for illustration on an optically thin IFN and infrared cirrus called the Spider.
However, for the available data the models were not perfect and revealed some degeneracy between effects of the properties of the dust and of the aIRSF.Hence our exploration of the potential benefits of additional independent constraints from complementary observations, namely the polarization of light scattered by the IFN.
In Section 2, we describe how to calculate the intensity and polarization of light scattered by a sphere resulting from illumination from a single incident direction, introducing the scattering geometry, Stokes parameters, and the phase functions for intensity and polarized intensity, illustrated for various dust models.
We present some models of the aISRF in Section 3 and describe how the illumination of the IFN from over the entire sky is used to calculate summary parameters for the intensity and polarization of scattered light seen by an observer at the Sun (Earth).
In Section 4 we provide illustrative results for an IFN at a fixed Galactic location, 1 varying the illumination model for a fixed dust model, and vice versa.
The summary intensity and polarization parameters also depend on the Galactic location, as demonstrated in Section 5.
We end with a discussion, summary, and conclusions in Section 6.

SINGLE SCATTERING BY A SPHERE
In this paper we concentrate on optically thin IFNs, because the 3D distribution of dust within infrared cirrus structures is not known.Lacking knowledge of the shape and alignment of dust in IFNs, a further simplifying assumption is that the dust particles are spherical.The scattering geometry is described in Figure 1 taken from Figure 4.1 right in Martin (1978).The scattering angle θ is the angle between the forward direction of incoming light (incident interstellar radiation field) and of for linear polarization to the orientation ψ of the electric vector, measured from North through East.For electric vector orientation along the vertical axis (North-South; ψ = 0), Q is positive and U = 0, for which we use +Q as shorthand.For electric vector rotated counterclockwise to ψ = 45 • , we have correspondingly +U (Q = 0).Negative and positive of each Stokes parameter correspond to electric vector orientations separated by 90 • .At 180 • the orientation is the same as for 0 • , and so Q is positive and U = 0.The actual polarization, e.g., the red line, is a pseudovector characterized by a length P and orientation (not direction) ψ defined in the half plane, with ψ ranging over [0 • , 180 • ] (see Equation 5).The signs of Q and U in each 45 • sector are indicated by the signs on the boundary black lines, stemming from the "projections" Q = P cos(2ψ) and U = P sin(2ψ).For an observer of scattered light originating as unpolarized light single-scattered by a sphere, "North" is taken to be along 2 in the scattering plane (Figure 1), and the scattered light is described entirely by I and Q in that frame of reference (U = 0).
the outgoing scattered light.The scattering plane is defined to contain these two vectors.There are two principal orientations of the electric vector, perpendicular and parallel to the scattering plane and denoted 1 and 2, respectively.The incident intensity is assumed to be unpolarized, and so components 1 and 2 are equal.However, upon scattering they are not equal and by symmetry the light has become linearly polarized in an orientation either perpendicular or parallel to the scattering plane.

Stokes parameters
The linear polarization of light can be described using the Stokes parameters, I, Q, and U .The total intensity is described by I and the linearly polarized intensity by the two parameters Q and U .2As described more fully in Figure 2, the polarization (red line) corresponding to the parameter pair (Q, U ) is a pseudovector with orientation defined in a half plane, range [0 • , 180 • ].For an observer of the scattered radiation,"North" is taken to be along 2 in the scattering plane.

Polarization of scattered light
Because the incident light is unpolarized, only Stokes parameter I i is non-zero.On scattering by a spherical dust grain, the radiation acquires polarization, with electric vector either parallel or perpendicular to the scattering plane (sp) and so for the scattering plane reference system U sp = 0.The intensities of the incoming light (the illumination) and scattered light are described by the Stokes parameters (I i , 0, 0) and (I s , Q sp , 0), respectively.The transformation maxtrix that acts on the incident (I i , 0) column vector to polarize the scattered light can be written (Martin 1978): where C s is the scattering cross section and M 1 and M 2 are calculated for incident electric vector perpendicular and parallel to the scattering plane, respectively.
Only the first column of the matrix is relevant because Q i = 0.In the calculations below, both the normalization factor of the illumination and the factor C s /4π are ignored, because they affect I s and Q sp of the scattered light equally.Zhang et al. (2023) discuss how the product of these factors can be related to the brightness of the thermal dust emission.
A familiarizing case re the matrix elements is Rayleigh scattering, where the size of the sphere is much smaller than the wavelength. 3In that case M 1 /2 = 3/4 and M 2 /2 = 3/4 × cos 2 θ.
For our applications, M 1 and M 2 are calculated using Mie theory from closed-form infinite series that are functions of the scattering angle θ.The dust models have size distributions for a given composition and up to two compositions.Elements M 1 and M 2 depend on size and composition, and their combined effect is accomplished through weighting by corresponding cross section C s , which depends on size and composition as well.
Note that the scattered light distribution depends only on θ and is axisymmetric around the direction of the incident radiation.The Mie calculations are carried out on a grid of θ within which values of the matrix elements at any θ can be interpolated as needed.The grid used is every degree from 0 • to 180 • inclusive (181 points).
The scattering angle used for the interpolants below was calculated as the arccos of the dot product on the unit sphere between the incident and scattered light vectors, where in our applications the scattered light of interest is that scattered toward Earth.We used cubic spline interpolation and comment on this choice below in Section 3.5.

Phase functions for intensity and polarized intensity
The phase function for intensity I ph describes the intensity distribution as a function of θ, relative to the average of the scattered light intensity integrated over the sphere (with this normalization, the integral of I ph over the sphere is 4π).In terms of the above, There is a corresponding phase function for polarized intensity, in the specified reference frame (I ph does not depend on the reference frame).Both I s = I i I ph and Q sp = I i Q ph scale in the same way with I i .

Phase functions for representative dust models
For purposes of illustration, we consider the DustEM dust model (Compiègne et al. 2011) as used by Zhang et al. (2023) and Astrodust (Draine & Hensley 2021).The DustEM model has two separate populations of grains, "LamC" and "aSil" (large amorphous carbon and amorphous silicate grains), whereas Astrodust has composite grains of amorphous carbon and silicate material; we adopt the Astrodust default model with porosity 0.2 and no iron inclusions (Hensley & Draine 2023).Composite grains such as these have been discussed before (Mathis 1996;Dwek 1997).We also consider two models by Weingartner & Draine (2001) (WD) that like DustEM are composed of large amorphous carbon and amorphous silicate grains.These are historically viable dust models constrained by many observational considerations (e.g., Hensley & Draine 2021), for which we are seeking further discrimination including their response to different aISRFs (Section 3).
For these and similar dust models the adopted size distributions are constrained by observations such as the Even in the event of color blindness, the curves in this and following figures are tied to their correct identity through the inset internal labels.
wavelength dependence of interstellar extinction.Most match a ratio of total to selective extinction R V ≈ 3.1 characteristic of the diffuse interstellar medium, but one (WD5.5)has size distributions reproducing R V ≈ 5.5, which is more characteristic of dense regions not encountered at intermediate and high Galactic latitudes.Results from WD5.5 can be compared to those from WD3.1 to show the differential effects from the changing size distribution (larger on average).

Phase functions for intensity, I ph
Phase functions I ph for these dust models are shown in Figure 3, calculated for the Sloan g (λ 0 ∼ 0.48) and Sloan r (λ 0 ∼ 0.62) bands used by the Dragonfly telescope array (Zhang et al. 2023).These can be contrasted to the phase function for Rayleigh scattering, which is close to isotropic (I ph = 1).Forward scattering is heavily favored.Relatively little light is scattered "sideways" and there is a slight back-scattering peak.
The strong tendency toward forward scattering can be quantified using the asymmetry parameter, G, which is the average cos θ weighted by I ph over the sphere.This can be evaluated in closed form from Mie theory and av- eraged over the size distribution.As seen by the entries in Table 1 in Section 4 below,4 G is close to 0.6 in the optical for these dust models, increasing with the ratio of the grain size to the wavelength.For the Rayleigh phase function G = 0, just as for isotropic scattering.
We also tabulate the albedo, ϖ, which is the ratio of the scattering cross section to the extinction cross section (itself the sum of the scattering cross section and absorption cross section).

Phase functions for polarized intensity, Q ph
Similarly, Figure 4 shows the polarized intensity phase functions for the dust models, with the Rayleigh scattering function for comparison.All functions are zero in the forward and back-scattering directions.Rayleigh scattering, with electric vector perpendicular to the scattering plane (negative Q ph ), produces strong polarization in the "sideways" direction.Overall, the amplitude for the dust models is much smaller.Peak negative Q ph for the dust models occurs in the forward-scattering hemi- sphere.In the back-scattering hemisphere Q ph changes sign (electric vector parallel to the scattering plane) and has a low amplitude positive peak.The behavior is more complex for WD5.5 in the g band, where the ratio of size to wavelength is largest.

ANISOTROPY OF THE INCIDENT ISRF
Even with anisotropic phase functions, the scattered radiation is isotropic and unpolarized unless the incident radiation making up the ISRF is also anisotropic.For illustrative purposes, we use the same three aISRF models considered by Zhang et al. (2023).

Sandage
One illumination model that we tested was proposed in Sandage (1976), in which an IFN is illuminated by the entire Galactic plane below, taken to be a slab of uniform brightness and extending to infinity.This produces an axisymmetric illuminating intensity proportional to 1/ cos Θ, where Θ is measured from the perpendicular dropped from the IFN to the plane.Different asymmetry can be explored by choosing Θ max < 90 • , thus illumination within a cone from a truncated slab.The total illumination integrated over the lower hemisphere is divergent as − log cos Θ max and so Sandage arbitrarily limited Θ max to be 89 • .This is the value used in Figure 5, top, where over the lower hemisphere the intensity ranges in relative units from unity to about 57.5 To capture such small-scale non-uniformity in a HEALpix representation (Section 3.5), an appropriately large N side is needed.6A value N side = 256 (786,432 pixels) is sufficient for these models, based on convergence of the computed p and ψ and of the integral of the intensity.

Uniform slab with internal extinction (USIE)
A second model is like the Sandage model, based on illumination from below by a uniform slab, but with the distinction that it has internal extinction characterized by optical depth τ ∞ from the midplane vertically to the slab edge.The intensity is then proportional to 1 − exp(−2τ ∞ / cos Θ). 7 In the limit τ ∞ → 0, it becomes the Sandage model for the same Θ max < 90 • .
A scale for τ ∞ is required and we adopted theirs (Mezger et al. 1982;Mathis et al. 1983), using the CCM extinction curve (Cardelli et al. 1989) for R V = 3.1 to obtain values of τ ∞ = 0.139 and 0.191 for the r and g bands, respectively.The r filter illumination displayed in Figure 5, middle panel, uses Θ max = 90 • , which can be computed because the internal extinction removes the singularity.The intensity over the lower hemisphere ranges in relative units from 0.24 to unity, much more uniform (isotropic) than the Sandage model.Again a value N side = 256 is sufficient.With increasing τ ∞ (as in the g band) the illumination becomes even more isotropic.

FRaNKIE
Third, we used the suite of aISRF models from the FRaNKIE code, a tool for calculating multi-wavelength interstellar emission in galaxies (Porter & Strong 2005;Porter et al. 2015Porter et al. , 2017) ) as part of the GALPROP code for cosmic-ray transport and diffuse emission production.These FRaNKIE models are based on realistic distributions of the sources of the radiation, including in the Galactic bulge, radiative transfer and dust, and the resulting illumination distribution changes with wavelength and location of the IFN (Galactocentric distance and height above the plane).From the grid available on the GALPROP web site site,8 we downloaded those most appropriate to application to the Spider IFN and the Dragonfly bands: Galactocentric distance 8.5 kpc, height 0.2 kpc, and filters at 0.6 and 0.45 µm.We also explored results for the illumination at different heights. 9he illumination displayed in Figure 5, bottom, for the Spider IFN near the r band, is seen to be more realistic than the other two models.The IFN is illuminated from all directions, not just from below, and there is an azimuthal asymmetry (with Galactic longitude) because the illumination from directions in the inner Galaxy is higher.

Rotation to a common reference frame
The polarization pseudovector (Q sp , 0) of the scattered radiation is relative to the scattering plane.For a constant scattering direction, toward the Earth, the scattering angle and orientation of the scattering plane change with the direction of incident light.Therefore, when considering more than a single direction of incident light, it is necessary to rotate the pseudovector on the sky about the vector joining the IFN to Earth to a common reference frame in which to sum the contributions to the individual Stokes components of (Q s , U s ).We chose the common reference frame to be the meridional plane containing the IFN, the Galactic North pole, and Earth, which we refer to as the sky reference frame with polarization in the Galactic convention.The vector joining the IFN to Earth is a common vector of the two planes, and the required amount of rotation of the column pseudovector is by twice the dihedral angle ϕ between the planes.The rotation matrix, of which only the first column is relevant because U sp = 0 in the scattering plane reference frame, is where ϕ is measured from North through East in the meridional Galactic frame.The scattered light intensity I s is of course independent of the reference frame.
In vector form, the dihedral angle is obtained via the dot product of the unit normals to the two planes and the normal for each plane is obtained from the cross product of the IFN to Earth vector and another vector in that plane.This can be reduced to dot products among the three vectors, appropriately normalized and coded to take into account edge cases: forward and back scattering and coplanar.

Using HEALPix over the whole sky
In practice, the whole sky illuminating intensity for FRaNKIE is given on a HEALPix sphere (Górski et al. 2005), with N side = 32 (12,288 pixels). 10We use HEALPix for the other illustrative aISRFs as well.HEALPix divides a sphere into equal-area pixels.Each indexed pixel corresponds to a single incident light vector as seen by dust in the IFN, with an associated radiation field intensity.We perform calculations of the scattered light intensity and polarized intensity pixel by pixel over the whole sky using the method outlined above, with the scattering direction always toward Earth.Summing the scattered light contributions to I, Q, and U in the common reference frame, weighted by the varying incident intensity over the whole sphere, results in a single Stokes parameterization (I s,tot , Q s,tot , U s,tot ) describing the scattered light seen at Earth.
For efficiency, the steps in the calculation can be grouped according the goals at hand.Here, the quantities of interest can be produced using vectors on the HEALPix indexed sky.For the scattering geometry, vectors ⃗ θ and ⃗ ϕ and then ⃗ cos 2ϕ and ⃗ sin 2ϕ can be (pre-)computed for a given N side .Then for different dust models and bands, their vectors ⃗ I ph and ⃗ Q ph can be produced by interpolation to the given ⃗ θ.Different incident aISRFs are represented by their illumination vectors ⃗ I i .The total scattered light I s,tot is obtained via the dot (scalar) product ⃗ I i • ⃗ I ph .The vector describing the polarized intensity in the scattering plane frame of reference ⃗ Q sp is obtained by the Hadamard (or Schur or element-wise) product ⃗ I i • ⃗ Q ph (like the dot product without the summation).The totals Q s,tot and U s,tot in the common reference frame are obtained by dot products ⃗ Q sp • ⃗ cos 2ϕ and ⃗ Q sp • ⃗ sin 2ϕ, respectively.To follow the steps during code development, we examined many intermediate vectors using Mollweide projection.
Simple checks of the accuracy of the interpolation and integration are available.For example, the numerical integral of the phase function over the sphere (as the sum of the elements of ⃗ I ph divided by the number of pixels, i.e., their average) can be compared to the expectation value unity.Similarly, G can be evaluated numerically (as the dot product ⃗ cos θ • ⃗ I ph divided by the number of pixels) for comparison to that computed directly from Mie theory.For example, for the typical phase functions used here these numerical integrals show agreement to within one in the fourth figure, using N side = 32 and linear interpolation in a uniform grid based on 181 angles θ.Cubic spline interpolation yields an order of magnitude improvement even for G, and so we adopt that.

Summary intensity and polarization parameters
One useful diagnostic factor is I, a relative intensity obtained by dividing I s,tot by the total illumination (the sum of the elements in ⃗ I i ).This is simply the average of the phase function I ph weighted by the illumination.The value of I is unity for an isotropic radiation field and/or a phase function with G = 0.A value of I < 1 means that the observer is in a position less favorable than average relative to the IFN, given the asymmetries of the aISRF and/or I ph , i.e., the IFN is dimmer as discussed by Zhang et al. (2023).
Two parameters that describe polarization independent of the scale of the illuminating intensity and the scattering cross section are the polarization fraction p, also called the degree of polarization, and the polarization angle ψ (see Figure 2), also called the position angle or orientation of the electric vector.The conventional formulae are where the two-parameter arctangent ensures that for the orientation ψ the correct sector in the half-plane is obtained (Figure 2).We use these formulae on (I s,tot , Q s,tot , U s,tot ) to obtain the summary parameters reported below.

RESULTS FOR THE SPIDER IFN
To demonstrate the diagnostic discrimination of polarization from scattered light, we report on the effects of different illumination and different dust models on the summary intensity and polarization parameters I, ψ, and p for the location of the Spider IFN.

Varying incident radiation fields
With the dust model held fixed, the summary parameters were calculated for the three different aISRFs.For the Sandage and USIE models, because of axisymmetry the results cannot depend on ℓ and only the location parameter b is required (40 • for Spider).Using a FRaNKIE model, both ℓ and b need to be specified, as well as the height above the Galactic plane (which implies a distance).These examples illustrate the additional and complementary discriminatory power of polarization observations in assessing proposed aISRFs.

Sandage and USIE models of illumination as a function of Θmax for fixed b
In the Sandage model, the (half) opening angle of the cone within which there is incident radiation is Θ max , and so this can be used to explore different geometries of the radiation field.The dependence of I, ψ, and p on Θ max for the DustEM model in the r and g bands is shown in Figure 6.
The case Θ max = 0 requires no integration over the sphere.There is a single incident direction and a single scattering angle, θ = 90 • + b.This corresponds to the backscattering lobe of I ph and I is simply I ph at that θ (130 • for Spider) and so less than 1.The value of p is simply |Q ph /I ph | at that θ.The scattering plane is equal to the meridional plane, and so the polarization arises simply from Q ph , which for the r band is slightly negative according to Figure 4, top panel, and so ψ = 90 • .However, for the g band, the polarization arises from the positive range of Q ph , and so ψ flips to 0 • or equivalently 180 • .Note that an IFN location at a higher b would also increase the scattering angle and induce the flip.
For increasing Θ max , more scattering angles and dihedral angles come into play.These are the same for the two bands and all dust models.For the Sandage model, the illumination is also the same and so what determines the shape of each curve is the dust model, for which the dependence of I ph and Q ph on θ is unique (Figures 3 and  4).For this dust model at g, there is a lot of cancellation among contributions to Q tot and U tot , causing a rapid change in ψ along with a minimum in p, before settling down to a behavior with Θ max like seen in the r band.
For Θ max approaching 90 • , the illumination is from a more and more extended slab filling the entire hemisphere.Still, all scattering angles are greater than 90 • and so I is less than unity.We have computed the summary parameters only to Θ max = 89 • because of the rapid increase of the illumination that cannot be resolved for the finite N side = 256 of the HEALpix representation.1.
Figure 6 also shows the corresponding results for the USIE illumination.As expected, they are identical to the Sandage results at Θ max = 0 • .Furthermore, the values of τ ∞ are not so large as to produce a radical qualitative difference in the relative illumination with latitude at intermediate Θ max , and so the net polarization is not much different than that produced by the Sandage model.However, at large Θ max , the internal extinction limits the excessive growth of the illumination and the results depart from those of the Sandage model.Even at Θ max = 90 • , I, ψ, and p can be calculated (e.g., the DustEM entries in Table 1).

FRaNKIE
Still for the DustEM models (r and g), we specify the illumination to be a FRaNKIE model.Now both ℓ and b need to be specified, because the radiation field is not axisymmetric.A height needs to be specified as well to select among the available FRaNKIE models described in Section 3.3.Focusing again on the Spider IFN, we adopt location (135 • , 40 • ) and height z = 0.2 kpc above the plane.
The summary parameters are tabulated in the DustEM rows and columns 7 -9 in Table 1.Compared to the results for the full hemisphere USIE illumination model (Θ max = 90 • ), shown in Figure 6 and in columns 4 -6 in Table 1, the value of I is higher, ψ is slightly lower, and p is significantly reduced.Again, this shows the possibility of discriminating between different aIS-RFs.

Results for various dust models with the same aISRF
For the same USIE illumination with Θ max = 90 • , and the Spider IFN location (b = 40 • ), results for several different dust models are compared in columns 4 -6 of Table 1.This illustrates the additional discriminatory power of scattered light polarization observations in assessing dust models, adding to the constraints from joint modeling of the intensity of thermal dust emission and scattered light intensity described by Zhang et al. (2023).
Table 1 shows, in columns 7 -9, that this dust model discrimination is also the case with the FRaNKIE illumination used in Section 4.1.2.

RESULTS FOR DUSTEM FOR VARYING IFN LOCATION
The intensity and polarization parameters are sensitive to the location of the IFN being observed in scattered light.Holding constant the dust model (DustEM) and the illumination (USIE with Θ max = 90 • ), we explored the sensitivity to b that arises through Q ph (there is no dependence on Galactic longitude).This can be calculated over the full range 0 • to 90 • .The results are shown in Figure 7.The summary intensity and polarization parameters in the r and g bands follow similar trends, differing only because of slight differences in I ph (Figure 3) and Q ph (Figure 4).One might question the physical basis for the model at b near 0 • .It is the only case where I is greater than unity because some radiation is scattered in the forward scattering lobe of I ph .There Q ph becomes small, suppressing the resulting p.As b increases, so do scattering angles θ, and I decreases while p increases.
At b near 90 • , the IFN is almost directly above the slab with scattering angles in the range [90 • , 180 • ].Cancellations among finite contributions from different illumination directions result in both Q s,tot and U s,tot near zero, so that p approaches zero.However, the distribution of Q sp , while symmetrical about the direction ℓ of the IFN, is displaced along it.Furthermore, sin 2ϕ has the same sign in the hemisphere centered on that ℓ, while cos 2ϕ does not, so that |U s,tot | is orders of magnitude greater than |Q s,tot |.For this radiation field ψ approaches 45 • (entirely net +U in Figure 2).FRaNKIE illumination at Galactocentric distance 8.5 kpc is available at discrete heights z [kpc] above the Galactic plane, ranging from 0 to 30 kpc.Though it instructive to examine all of these, we take the models at z = [0.2,0.5, 1] kpc as most likely to be relevant for IFNs.
For one set of illustrative polarization calculations we hold b constant, so that increasing z corresponds to increasing distance from the Sun (d = 0.31(z/0.2) sin(40)/ sin(b) kpc).We adopt (ℓ, b) = (135 • , 40 • ) and the dust model DustEM so that we overlap the previous calculations in Section 4.1.2and Table 1 (right).The results are tabulated in Table 2, columns 3 -5.With increasing height and distance, I decreases while p increases and there is a modest increase in ψ.Because the scattering geometry remains the same (θ, ϕ, I ph , Q ph ), these systematic changes originate in the greater anisotropy of the illumination with respect to Galactic latitude.
For another set of illustrative calculations, we can pair a corresponding b with z via b = arctan(tan(40)z/0.2),which corresponds to IFNs whose distance from the Sun projected in the Galactic plane is held constant.We again fix ℓ = 135 • and tabulate the results in Table 2, columns 7 -9.Now there is a larger decrease of I with height, as b increases.There is a substantial change in ψ accompanied by a minimum in p.These arise from the additional effects of scattering geometry.Generally, I, ϕ, and p are diagnostic of the aISRF and scattering geometry, for a given dust model (in this case DustEM).Unlike the other illumination models, the FRaNKIE models are anisotropic in longitude (as well as in latitude).This results in the variations of I, ψ, and p with ℓ shown in Figure 8, for z = 0.2 kpc and b = 40 • and for DustEM (solid curves).The values at ℓ = 135 • were tabulated in Tables 1 and 2.
In the r band, the I curves for both dust models are not perfectly symmetrical about ℓ = 0 • .There is an even greater asymmetry in g.The asymmetry reflects the fact that the FRaNKIE illumination model is not symmetrical about ℓ = 0 • .Nor in reality is the actual aISRF in the Milky Way expected to be symmetrical about ℓ = 0 • , because of the distributions of stars (bar, spiral arms) and interstellar dust.
Recall that when the IFN is in a location relative to observer less favorable than average then I< 1, because of the coupling of the θ dependence of the phase function I ph and the illumination.The least favorable location is for IFNs toward the anticenter (around ℓ = ±180 • ), where backscattering is predominant.Locations in the direction of the Galactic bulge (around ℓ = 0 • ) are somewhat more favorable because of contributions from the forward scattering lobe.
The ψ and p curves are also not symmetrical about ℓ = 0 • .The polarization asymmetry is also greater in g.The polarization asymmetry now reflects the coupling of the more complex θ dependence of the polarization phase function Q ph and the illumination asymmetry, with the latter having a dominant role.11

For the Astrodust model
We carried out the same calculations for the Astrodust model (dashed curves in Figure 8).The variation of p with ℓ is smaller for Astrodust than for DustEM.This is because of more cancellation of contributions from the positive and negative lobes of the Astrodust polarization phase function, Q ph (see Figure 4, top).As well, the similarity of the polarization parameter variations of Astrodust in the r band and of DustEM in the g band is related to the similarity of their corresponding Q ph in these respective bands.

DISCUSSION AND CONCLUSIONS
Considering optically thin IFNs, Zhang et al. (2023) demonstrated that observations of thermal dust emission combined with any dust model and aISRF model could be used to predict the expected intensity of scattered light.Comparing this prediction with the observed intensity of scattered light tests combinations of dust models and aISRF models and when applied to the Spider IFN demonstrated that further refinement of the models was needed to get a closer match of predictions and observations.Furthermore, it was recognized that it would be difficult to pin down deficiencies in the dust models separately from those in the aISRF models.Such joint modeling is related to the scattered light parameter I, which we have shown above is also dependent on the location of the IFN in the Galaxy.
We therefore undertook to predict the polarization properties of the scattered light as well, namely polarization angle ψ and polarization fraction p, to assess whether they would be useful in providing complementary diagnostics of dust models and aISRF models to help break the degeneracies, including for IFNs in different locations.After developing the formalism for these calculations, we demonstrated systematically through many examples spanning different combinations of dust models, aISRFs, and Galactic locations that the potential for discrimination is there.
However, it should be appreciated from these calculations that such observations are challenging, with optically thin IFNs being only a few percent of the brightness of the night sky.According to our results across many models, the polarization fraction is likely to be only a few percent and not more than 10 %.This is The vertical line corresponds to calculations for the Spider IFN in columns 4 -6 and DustEM rows in Table 1.
much less than one might intuitively hope, based on familiarity with Rayleigh scattering in the daytime sky, where the magnitude of Q ph is much more favorable and the illumination is very anisotropic.This is a reminder that moonlight is to be avoided and that other sources of illumination of the night sky could induce an unwanted polarization contamination as well.Nevertheless, this contamination would be smooth across a field, unlike the intensity and polarized intensities of the scattered light that are modulated by the column density of  dust in the IFN.The spatially different behavior would enable component separation.
Interstellar polarization produced by aligned dust grains in the foreground could contaminate the polarization as well.It is unlikely to be oriented with the same ψ.If this contamination were not removed, and the wavelength dependence of p of the scattered light were different than that of interstellar polarization, then the net ψ would change (even more) with wavelength.For IFNs at intermediate to high Galactic latitudes, the interstellar polarization would likely be smaller than p for the scattered light; it could be mapped for foreground stars whose distances are known from Gaia (Gaia Collaboration et al. 2016Collaboration et al. , 2023) ) and its Stokes parameters subtracted.
By comparison to simple diffuse Galactic light models, Chellew et al. (2022) found evidence for an excess at 0.65 µm that could arise from the "extended red emission" (ERE) phenomenon (Witt et al. 2008).The ERE is thought to be excited by near-UV photons and emitted isotropically, not subject to I ph .The ERE seems unlikely to be polarized.Its admixture in an IFN would therefore frustrate the simplest modeling in the r band.
We have shown that location of the IFN is important.One promising IFN to contrast to the Spider is the Draco nebula discussed by Zhang et al. (2023).The Draco field is centered on (91 • , 38 • ) and is at a distance of about 600 pc, thus at a height z about 0.37 kpc above the Galactic plane.Being at a different Galactic location, the scattering geometry would provide a distinctly different coupling to the anisotropy of the aISRF.In Table 3 we have previewed the magnitude of the differences that might be expected, using the DustEM model and FRaNKIE illumination at various heights.The first row for each band reproduces the previous results for the Spider from Tables 1 (right) and 2. The second row is for the Draco coordinates and the same FRaNKIE illumination for 0.2 kpc, which is probably too low; because b is similar to that for Spider, these results are close to those shown in Figure 8.The fourth row is for the Draco location but FRaNKIE illumination for 0.5 kpc, which is probably too high.In the third row, we simply averaged the two FRaNKIE files to provide an estimate for 0.35 kpc (note that the simple average of the values of I for 0.2 and 0.5 kpc is very close to the tabulated value, the average for ψ is only about 0. • 045 degrees lower, and the average for p is higher by a factor less than 1.025).The Draco values are certainly distinguishable from those for the Spider.
In conclusion, we have found that observations of I, ψ, and p of scattered light, along with the intensity of thermal dust emission, are a powerful combination with which to assess the viability of dust models and models of the aISRF, but they need to be carried out with care and precision.
Looking at the various tables and figures, one can see that depending on the combinations of aISRF, dust model, and location either r or g could offer more discrimination.Quite often g might be preferred.On the other hand, both the brightness and the polarization fraction can be lower at g, making observations more challenging; this would worsen to the blue.Moving further to the red than r, the brightness is again lower (Zhang et al. 2023).Overall, using a combination of both r and g seems optimal, rather than devoting observational resources to just one.Both Planck thermal dust polarization and interstellar polarization of starlight indicate that in reality dust particles are aspherical and aligned with respect to the ambient magnetic field (Planck Collaboration XII 2020) though the details in any particular IFN are presently difficult to prescribe.Though the nature of the asphericity is uncertain, the scattered light properties might be tackled in a next approximation with calculations for aligned spheroids, where I ph and Q ph would depend on the direction of the incident radiation with respect to the instantaneous orientation of the spheroid.There would be averaging for spinning/tumbling spheroids over a range of shapes and sizes and illumination by the allsky aISRF and one might reasonably expect (though it has not been demonstrated) that the net results would be close to those obtained for spheres and that trends in results for the r band relative to the g band would be preserved.
We can expect substantial improvement both in our understanding of how the aIRSF changes with position in the Galaxy, for example thanks to constraining data now available from Gaia (Gaia Collaboration et al. 2016, 2023), and in dust models.These will have a salutary impact on the testing previewed here, limiting the parameter space for model solutions consistent with the data.

Figure 1 .
Figure1.The scattering geometry, including the scattering angle θ, the scattering plane outlined in dashed lines, the two principal orientations of the electric vector, 1 and 2, and the direction of propagation, 3.

Figure 2 .
Figure2.Relationship of the Stokes parameters Q and U for linear polarization to the orientation ψ of the electric vector, measured from North through East.For electric vector orientation along the vertical axis (North-South; ψ = 0), Q is positive and U = 0, for which we use +Q as shorthand.For electric vector rotated counterclockwise to ψ = 45 • , we have correspondingly +U (Q = 0).Negative and positive of each Stokes parameter correspond to electric vector orientations separated by 90 • .At 180 • the orientation is the same as for 0 • , and so Q is positive and U = 0.The actual polarization, e.g., the red line, is a pseudovector characterized by a length P and orientation (not direction) ψ defined in the half plane, with ψ ranging over [0 • , 180 • ] (see Equation5).The signs of Q and U in each 45 • sector are indicated by the signs on the boundary black lines, stemming from the "projections" Q = P cos(2ψ) and U = P sin(2ψ).For an observer of scattered light originating as unpolarized light single-scattered by a sphere, "North" is taken to be along 2 in the scattering plane (Figure1), and the scattered light is described entirely by I and Q in that frame of reference (U = 0).

Figure 3 .
Figure 3. Phase functions I ph for various dust models in the r and g bands.The phase functions for Rayleigh scattering and isotropic scattering (unity) are included for comparison.Top: DustEM and Astrodust.Bottom: WD3.1 and WD5.5.Even in the event of color blindness, the curves in this and following figures are tied to their correct identity through the inset internal labels.

Figure 4 .
Figure 4. Polarized intensity phase functions Q ph for various dust models in the r and g bands.The function for Rayleigh scattering (scaled down by 1/3) is included for comparison.Top: DustEM and Astrodust.Bottom: WD3.1 and WD5.5.

5. 1 .
Location as a function of b, for USIE illumination with Θ max = 90 •

5. 2 .
Location as a function of height, for FRaNKIE illumination, ℓ = 135 • , and either constant b = 40 • or changing b

5. 3 .
Location as a function of ℓ, for FRaNKIE illumination at z = 0.2 kpc and b = 40 •

Figure 7 .
Figure 7. Effect of IFN location in Galactic latitude b (ℓ = 135 • ) on summary scattered light parameters for USIE illumination (Θmax = 90 • ) and DustEM grain model in the r and g bands.Top: relative intensity.Middle: polarization angle in degrees.Bottom: polarization fraction in percent.The vertical line corresponds to calculations for the Spider IFN in columns 4 -6 and DustEM rows in Table1.

Figure 8 .
Figure 8.Effect of IFN location in Galactic longitude ℓ (b = 40 • ) on summary scattered light parameters for FRaNKIE illumination (z = 0.2 kpc) and two dust models in the r and g bands (DustEM and Astrodust).Top: relative intensity.Middle: polarization angle in degrees.Bottom: polarization fraction in percent.The vertical line corresponds to calculations for the Spider IFN in columns 7 -9 and DustEM and Astrodust rows inTable 1 and the z = 0.2 kpc, b = 40 • rows in Table 2 for DustEM.

Table 2 .
Summary scattered light parameters for FRaNKIE illumination at various heights z above the Galactic plane, ℓ = 135 • , and DustEM model Table 1 and the z = 0.2 kpc, b = 40 • rows in Table 2 for DustEM.

Table 3 .
Summary scattered light parameters for DustEM model and FRaNKIE illumination for different heights for Spider and Draco locations