Interstellar Polarization Survey. IV. Characterizing the Magnetic Field Strength and Turbulent Dispersion Using Optical Starlight Polarization in the Diffuse Interstellar Medium

Angular dispersion functions are typically used to estimate the fluctuations in polarization angle around the mean magnetic field orientation in dense regions, such as molecular clouds. The technique provides accurate turbulent-to-regular magnetic field ratios, 〈Bt2〉1/2/Bpos , which are often underestimated by the classic Davis–Chandrasekhar–Fermi method. We assess the technique's suitability to characterize the turbulent and regular plane-of-sky magnetic field in low-density structures of the nearby interstellar medium (ISM), particularly when the turbulence outer scale, δ, is smaller than the smallest scale observed, ℓ min. We use optical polarization maps of three intermediate-latitude fields (∣b∣ ≳ 7.°5) with dimensions of 0.°3 × 0.°3, sourced from the Interstellar Polarization Survey–General ISM catalog. We decomposed the H i emission detected by the Galactic All-Sky Survey within our fields to estimate the multiphase ISM properties associated with the structure coupled to the magnetic field. We produced maps of the plane-of-sky magnetic field strength (B pos), mass density (ρ), and turbulent velocity dispersion (σ v,turb). In the regions with well-defined structures at d < 400 pc, the average B – ranges from ∼3 μG to ∼9 μ G, depending on the method and physical properties. In the region where structures extend up to 1000 pc, B pos varies from ∼1 to ∼3 μ G. The results agree with previous estimations in the local, diffuse ISM. Finally, optical starlight polarization and thermal dust polarization at 353 GHz consistently reveal a highly regular plane-of-sky magnetic field orientation unfazed by diffuse dust structures observed at 12 μm.


INTRODUCTION
The Galactic magnetic field (GMF) interacts with dust, gas, and cosmic rays, influencing the interstellar medium (ISM) from the diffuse medium to molecular clouds and star-formation regions (e.g., see Crutcher 2012;Haverkorn 2015;Vázquez-Semadeni 2015, and references therein).Starlight polarization produced by the dichroic extinction of background stars reveals the orientation of the Galactic magnetic field (GMF) projected on the plane of the sky (Hiltner 1949;Davis & Greenstein 1951).Moreover, linear polarization observations are useful for estimating the average plane-of-sky magnetic field strength with the Davis-Chandrasekhar-Fermi method (Davis 1951;Chandrasekhar & Fermi 1953) or any analogous approach.The DCF method estimates the magnetic field strength based on the dispersion of the polarization angle under the assumption that this dispersion measures the turbulence in the medium.In cases where this assumption is not true, a structure function (SF) or angular dispersion function (ADF) can be used to estimate the dispersion around the large-scale field and, subsequently, the magnetic field strength.
The second-order SF of the polarization angle, ⟨∆θ 2 (ℓ)⟩, measures the behavior of the dispersion about the mean magnetic field direction as a function of angular or spatial distance between two measurements.The SF method was first explored by Serkowski (1958Serkowski ( , 1962) ) and then by Kobulnicky et al. (1994) to study the polarizing properties of the ISM in open clusters.Later, Falceta-Gonçalves et al. (2008) used the SFs along with simulations of polarization maps to characterize the turbulent cascade and the magnetic field.Hildebrand et al. (2009) defined the ADF, ⟨∆θ 2 (ℓ)⟩1/2 , to measure the turbulent and large-scale contribution to the angular dispersion in molecular clouds.The technique (hereafter referred to as ADF-DCF) finds the magnetic field strength by avoiding inaccurate estimates of the field dispersion, σ θ , based on simple Gaussian fits to the polarization angle distribution.However, the ADF-DCF method is valid only at scale lengths larger than the turbulent correlation length, δ, or the observed resolution scale.
The ADF-DCF method and its variations (Houde et al. 2009) have been extensively used with starlight and submillimeter polarization maps of discrete and localized dense regions, e.g., nearby molecular clouds (e.g., see Hildebrand et al. 2009;Franco et al. 2010;Soler et al. 2016;Wang et al. 2019, and references therein), whose physical properties are well known.In the diffuse ISM, on the other hand, starlight polarization observations have much lower spatial density (e.g., see Heiles 2000;Panopoulou et al. 2023) and there is little knowledge about the physical properties of diffuse structures (i.e., precise measurements of velocity dispersion and densities).Fortunately, optical polarimetry data from the Interstellar Polarization Survey (IPS) (Magalhães et al. 2005 andMagalhães et al. 2024, in preparation) in the diffuse ISM (IPS-GI, Versteeg et al. 2023) offers a great opportunity to study the properties of the interstellar magnetic fields in diffuse dust close to the Galactic thin disk.
The IPS-GI database (Versteeg et al. 2023) has demonstrated the potential of optical polarization maps, combined with precise Gaia-EDR3 (Gaia Collaboration et al. 2021) distance, optical extinction (Anders et al. 2022), and atomic and molecular line emission to study the magnetized ISM (e.g., see Angarita et al. 2023 andVersteeg et al. 2024).In this paper, we use the IPS-GI optical polarization measurements at intermediatelatitude (|b| > 7. • 5) regions to test the suitability of the classic DCF, the ADF-DCF, and other DCF-like techniques to study the magnetic field properties in the diffuse medium and at small angular scales (< 0. • 3).We report the plane-of-sky magnetic field and ISM properties of the nearby (d < 1, 000 pc) diffuse polarizing clouds deduced from observations.In Section 2, we present the polarization data and ancillary data needed to estimate the properties of the interstellar structures.We explain the DCF, the ADF-DCF, and other methods to estimate the magnetic field strength in Section 3. We describe the polarimetric properties of the selected regions and calculate their physical properties using atomic hydrogen emission and three-dimensional (3D) dust maps in Section 4. We present our results on the turbulent dispersion and magnetic field strength in Section 5. We discuss the turbulent scales observed in our sample and the magnetic field strength obtained from different methods, compare our results with other estimations, and outline the characteristics of the nearby magnetic field in Section 6.Finally, we summarize our work in Section 7.

OBSERVATIONS AND DATA
In this section, we present the optical linear polarization observations and ancillary data used in our research.

Optical polarization data
This research uses data from the Interstellar Polarization Survey (IPS) (Magalhães et al. 2005 andMagalhães et al. 2024, in preparation) in the general ISM, which was reduced and presented in detail by Versteeg et al. (2023).The stellar catalog (IPS-GI) consists of 38 fields of view, each measuring 0. • 3 × 0. • 3, sparsely distributed near and within the Galactic thin disk in the Southern sky.The photometric and polarimetric parameters were obtained with the reduction pipeline SOLVEPOL (Ramírez et al. 2017).We use the quality filters presented in Versteeg et al. (2023) and Angarita et al. (2023), including a strong constraint in the signal-to-noise of P/σ P ≥ 5, which ensures the measurements are not significantly affected by the Ricean bias, as discussed by Ramírez et al. (2017) and Versteeg et al. (2023).The low instrumental polarization, 0.07%, below the median fractional polarization error measured, does not need to be corrected.The final catalog consists of ∼10, 500 stars with high-quality optical polarization measurements, distance, optical extinction, and among other parameters from Anders et al. (2022) and Gaia Collaboration et al. (2021).

H I map
The Galactic All-Sky Survey (GASS, McClure-Griffiths et al. 2009) has atomic hydrogen (H I) spectral line emission data for the entire Southern sky.The survey has a resolution of 16 ′ , comparable with the size of the IPS-GI field-of-view.Therefore, we used positionposition-velocity (PPV) data cubes 1 of neutral atomic hydrogen within an area of 0. • 7 × 0. • 7, centered on the IPS-GI fields' coordinates (Section 4).The velocity range covers from −400 km s −1 to 500 km s −1 with 1 km s −1 of spectral resolution.The average rms noise in the H I map is 0.057 K.

3D Dust maps of the local Galaxy
High-resolution 3D dust maps give us accurate insights into the dust morphology and properties along our lines of sight in the local Galactic neighborhood.We considered three 3D maps: the Vergely et al. (2022) dust map, with 10 pc resolution covering roughly 1.5 kpc and sampling 5 pc; the Leike et al. (2020) map, with radii ranging from about 370 pc to 590 pc and resolution of 1 pc; and the Edenhofer et al. (2023) map, spanning 1.25 kpc with resolutions between 0.4 − 7 pc.We compare the three dust maps with various resolutions to look for consistencies between the morphology of line-of-sight features in each field studied.We primarily use the Vergely et al. (2022) 3D dust map, retrieved through the G-Tomo 2 platform, to estimate the total column density along the line of sight.We queried Leike et al. (2020) and Edenhofer et al. (2023) high-resolution 3D dust maps through the Python Dustmaps module 3 (Green 2018).We use Leike et al. (2020) for the volume density and spatial scale estimations in nearby clouds.Additionally, we use Edenhofer et al. (2023) for comparison with the other maps.We obtained the extinction density, the cumulative extinction profiles with distance, and their respective errors using the central coordinates of the selected IPS-GI fields (see Section 4).

MAGNETIC FIELD STRENGTH CALCULATION
The total magnetic field is assumed to have a regular but not constant component, B pos , and a turbulent component, B t , in the plane of the sky.Polarization measurements combined with interstellar gas properties allow us to estimate the regular and turbulent plane-ofsky magnetic field strength through the DCF method and its modifications under certain considerations.We described the methods and their assumptions in the following sections.

The DCF method
The DCF method estimates the average plane-ofsky magnetic field strength from polarization angle deviations around a mean magnetic field due to small, non-thermal turbulent gas motions (Davis 1951;Chandrasekhar & Fermi 1953).The regular plane-of-sky magnetic field strength is estimated by the DCF method as where σ v and ρ are the velocity dispersion and mass density of the magnetized medium, respectively, f is the fudge factor that corrects the magnetic field strength estimation, and σ θ is the dispersion of the polarization angle.To obtain Equation (1), the method assumes 2 https://explore-platform.eu/sdas 3 https://dustmaps.readthedocs.io/en/latest/index.html equipartition between the gas turbulent kinetic energy and magnetic energy, isotropic turbulence (i.e., σ vpos = √ 2σ v los ), and a small random magnetic field component (i.e., B t ≪ B pos ).Furthermore, the DCF method assumes that the field dispersion is represented entirely by the deviations in polarization angles at the observed scales (i.e., B t /B pos ≈ σ θ ), and this dispersion is solely due to gas turbulent motion composed of Alfvén modes, i.e., incompressible transverse magnetohydrodynamic (MHD) waves (e.g., see Crutcher 2012;Beck & Wielebinski 2013;Jones 2015;Andersson et al. 2015;Liu et al. 2021, and references therein for more details).
It is well-known that interstellar turbulence is anisotropic.Indeed, the theory by Goldreich & Sridhar (1995) predicts larger anisotropies for smaller scales, with turbulent eddies aligning with the local magnetic field.However, as the size of turbulent cells decreases, more cells contribute to the line-of-sight integration, each with its own local magnetic field orientation.Therefore, we assume that the anisotropies average out on these scales, allowing us to use the DCF method for isotropic turbulence.Studies on the DCF method based on synthetic maps created from numerical simulations, which include the natural anisotropies in the MHD turbulence below the injection scale, have demonstrated that under these assumptions, the DCF method for isotropic turbulence provides reasonable results (Cho & Yoo 2016).

The Skalidis & Tassis (2021) method
The Skalidis & Tassis (2021) method (hereinafter ST-DCF) is a modification of the classic DCF method that calculates the plane-of-sky magnetic field strength as Contrary to the classic DCF method, the ST-DCF assumes that the fluctuating magnetic energy is dominated by the cross-term (i.e., B t • B pos ̸ = 0), which represents the compressible modes.Although this assumption may not be valid in all cases, e.g., in dense regions where selfgravity is important (Liu et al. 2021(Liu et al. , 2022)), in the case of the diffuse ISM, ST-DCF may provide reliable estimates of the regular plane-of-sky magnetic field strength.The method demonstrates a mean relative deviation of 17% from the true value in 3D compressible MHD simulations, with its accuracy independent of the turbulence properties.

Angular dispersion function
The ADF is defined as the squared root of the secondorder SF of the polarization angle (Hildebrand et al. 2009;Poidevin et al. 2010), where θ(x) − θ(x + ℓ) is the difference in polarization angle between pairs of stars with angular separation ℓ, and N is the number of pairs per angular separation bin, which has N (N − 1)/2 distinct combinations (see Falceta-Gonçalves et al. 2008).The ∆θ magnitude is constrained in the range [0, 90] • , but in practice, ∆θ is confined to small angles to ensure that the small-angle approximation in the DCF method is valid and to avoid the nπ ambiguity.
As described by Hildebrand et al. (2009), if the observed scales ℓ are much smaller than the maximum scale of variation in B pos , this variation increases approximately linearly with scale length and slope m.If these scales, in addition, are larger than the maximum scale of turbulence δ, the turbulence contributes a constant value denoted by b.Since the two contributions are statistically independent, we can add them in quadrature with the contribution from the measurement uncertainties, σ M,k in the k−th scale-length bin, to find the total observed SF of polarization angle for small ∆θ(ℓ): where Houde et al. (2009) presented a similar approach as in Equation ( 4), departing from the same basic assumptions as Hildebrand et al. (2009) and assuming a Gaussian form for the autocorrelation function for the turbulence.They expanded the study of the dispersion by including the contributions from the signal integration through the thickness of the clouds and over the telescope beam.The method applies to starlight polarization by setting the angular resolution term as W = 0.This approach is valuable for determining the turbulent correlation length δ when the SF or the ADF drops at scales nearing ℓ = 0, possibly suggesting that δ is being resolved (Section 3.3.1),e.g., as in Franco et al. (2010) and Wang et al. (2019).Nevertheless, as reported in Section 5.1, this is not the case for our IPS-GI regions in the diffuse ISM.Therefore, we use the ADF (Hildebrand et al. 2009;Poidevin et al. 2010) corrected for the measurement errors σ M , as for straightforward comparison with the polarization angle dispersion and subsequently fit Equation (4) to the scales at which the ADF increases approximately linearly.

Interpretation of the ADF
A qualitative analysis of the ADF profile as a function of the angular scale length can give us an idea of the properties of the magnetic field and the ISM in the observed region: 1. Flat ADF profiles mean that the fluctuations in polarization angle are scale-independent over the scales probed.This can indicate random noise or that fluctuations in polarization angle exist only on scales smaller than the smallest scale probed by the observations ℓ min (Sun & Han 2004).
2. The ADFs may exhibit multiple bumps or wiggles at different scales due to the contributions of more than one discrete structure that differ in polarization properties (Sun & Han 2004) or due to irregular sampling of sources (Soler et al. 2016).
In such a case, care must be taken in interpreting each slope or bump.
3. Most commonly, ADFs increase with scale ℓ, indicating that the amount of fluctuations grows with scale.This is generally interpreted in either of two ways: (a) The increasing ADF reflects purely turbulence, which shows a power-law behavior in log-log space (e.g., Falceta-Gonçalves et al. 2008).If the ADF flattens out, the scale of this turn-over can be interpreted as the turbulence correlation length δ (Haverkorn et al. 2008).(b) If the measured scales are larger than the turbulence correlation length, i.e., δ < ℓ, the rising ADF is caused by gradients and other large-scale variations in the regular magnetic field.In such a case, the slope represents variations in the field due to a gradual change in polarization angle across the field of view.
A change in the slope of an ADF as a function of scale may indicate detection of the turbulent outer scale (Houde et al. 2009), may reflect properties of the turbulence (Falceta-Gonçalves et al. 2008), or suggest various sources of structure (Haverkorn et al. 2004).
As discussed in Section 4.1.1,the case 3.b may apply to our data.In this case, characteristics of the ADF give a better estimation of polarization angle fluctuations than its dispersion, σ θ (e.g., see Liu et al. 2021).Therefore, the ADFs with DCF analysis combination gives an estimate of the magnetic field in the region (Hildebrand et al. 2009;Houde et al. 2009).We call this modified DCF method as ADF-DCF hereinafter.We describe the generalities of the technique and its difference with DCF in the following section.

Implementing ADFs in the DCF method: ADF-DCF
We can now apply the DCF method, which traditionally is used by equating the turbulent to regular magnetic field ratio ⟨B 2 t ⟩ 1/2 /B pos to the dispersion in polarization angle σ θ .This assumes that the polarization angle dispersion is caused solely by turbulence, which is not completely valid in the case 3.b discussed above.Considering the turbulent contribution in Equation (4), Hildebrand et al. (2009) showed that, in this case, the turbulent to regular magnetic field ratio can be described in terms of the turbulent parameter b as Equation ( 6) can then be used to define the magnetic field strength in terms of the gas velocity dispersion, the density, and the turbulent parameter as If the observed scale is smaller than the turbulent scale (i.e., δ > ℓ max ), we can no longer use Equation ( 7) to calculate the magnetic field strength.Instead, we can implement the classic DCF method (Equation (1), Section 3.1) or the Skalidis & Tassis (2021) method (Equation (2), Section 3.2), assuming that the field of view is larger or equal to the outer scale of fluctuations so that the angle dispersion is a good measure of the magnetic field fluctuations in the DCF approximation.

IPS-GI field selection criteria
The ADF probes the plane-of-sky magnetic field fluctuations by comparing the polarization angles between different lines of sight as a function of the physical distance between them.As we observe angular distances, the 2D ADF defined in Equation (3) can only be applied when we can map angular to physical distances.This is only the case if the dust distribution is similar for all stars, i.e., we probe the same dust in front of all stars in the sample.
In that context, we focus on the IPS-GI fields where starlight passes through a single, clearly defined, diffuse polarizing dust structure.The IPS-GI intermediate Galactic latitude fields, i.e., C2, C16, and perhaps C37, at |b| ≳ 7. • 5, are good candidates for this study.Angarita et al. (2023) found that these fields have, on average, low optical dust extinction (A V < 1 mag).Additionally, Angarita et al. (2023) demonstrated that optical starlight polarization in C2, C16, and C37 originates in the nearby ISM, in or just behind the Local Bubble (LB) wall.This is evident, especially in C2 and C16, as most stars behind the dust structures exhibit a constant degree of polarization regardless of distance.
Figure 1 shows the extinction density profile with distance in the three fields.We show the models from Vergely et al. (2022), Leike et al. (2020), andEdenhofer et al. (2023) (Section 2.3) to assess the consistency of the dust structures at different resolutions.We see dust structures at d < 400 pc in C2 and C16, and d < 1, 000 pc in C37.Optical polarization measurements (the gray dots) are mostly behind these dust structures, confirming the origin of the polarization mea-sured in the foreground.Therefore, we can assume that all lines of sight cross the polarizing dust at comparable distances, and the angular scales between stars correspond to actual physical scales within the cloud (Figure 1).The assumption holds less for C37.Unfortunately, there are only a handful of data points at d < 1, 000 pc in C37 (Figure 1), so we can not disentangle the polarization properties of the individual structures.However, the low density of the dust, the constant optical polarization degree beyond 1, 000 pc, and the highly regular magnetic field in C37 (see Angarita et al. 2023, for more details) make it feasible to assume uniform polarizing properties in all dust structures along a pathlength of 1, 000 pc.

Polarization angle distribution in intermediate-latitude IPS-GI fields
As explained in Section 3.3.1, the interpretation of the ADF depends on the relative observed scales concerning the turbulence correlation length.Interpretation of ADFs is different whether observed fluctuations in polarization angle are part of a turbulent cascade or largescale gradients in the cloud.The IPS-GI field-of-view is exceedingly small (∼1 pc at the distance of the dust peak of Figure 1) but may be comparable to the turbulent outer scales in cold clouds (see Section 6.1).The large-scale distribution of polarization angle may give a clue as to the question of whether observed fluctuations in polarization angle are due to turbulence or gradients.
In Figure 2 (left), we see a gradual change in polarization angle across the fields of view.The smooth change in orientation by approximately 10 • suggests that the dispersion is likely caused by large-scale fluctuations from a structure larger than the field of view, but partially inside the field (i.e., δ < ℓ), rather than pure turbulence (i.e., δ > ℓ max ).Although it may be impossible to make the distinction definitively, the gradual change in the median θ is statistically significant according to the standard error in the right column of Figure 2. Interestingly, the direction of the change in the magnetic field orientation coincides with the direction of the H I emission gradient in both the observations and the models, lending some more credence to the idea that the angle gradient across the field is part of a large-scale structure (see Section 4.2 for details in H I observations and modeling, see also the Figures in Appendix A for comparison).

Extreme values and subsampling issues
The ADFs are a powerful tool to measure fluctuations in the ISM.Still, they are very sensitive to extreme values and non-uniform distributions of the measurements (e.g., see Sun & Han 2004;Soler et al. 2016).Firstly, to account for outlier issues, we removed data points with polarization angles exceeding three times the standard deviation, as they showed no significant correlation with structures along the lines of sight (the red polarization vectors and data points in Figure 3).The atypical polarization values are not correlated with either distance or sky position and are likely due to circumstellar matter causing intrinsic polarization in these stars.Additionally, it can also be caused by observational issues, for instance, in crowded regions where the superposition of different sources can produce mismatches between the ordinary and extraordinary images created by the Savart prism (Magalhães et al. 2005;Versteeg et al. 2023).Secondly, we do not trust angular scales above 2/3 (i.e., ∼0.• 2) of the observed field width since subsampling leads to large uncertainties.Thirdly, polarization angle dispersion enhancements or "bumps" may happen at a certain range of scales, leading to more than one slope in the ADF profile.Given that all fields studied present a gradient in the median polarization angle (Figure 2), we attribute the bumps to irregular sampling of the sources (Soler et al. 2016).

ISM properties in IPS-GI fields
To characterize the magnetic field strength, we need density ρ and velocity dispersion σ v measurements of the ISM structure coupled to the magnetic field.The emission and absorption of atomic and molecular lines such as CO, K I, Ca II, OH, CN, and H I, among others, have been extensively used to estimate the prop- We use H I PPV cubes from GASS to find the gas properties that better correlate with the polarizing dust structure in our lines of sight.To this end, we perform a Gaussian decomposition using the Regularized Optimization tool for Hyper-Spectra Analysis, ROHSA (Marchal et al. 2019), which provides comprehensive information on the multiphase interstellar gas properties.This approach should be used with caution since not all neutral H I gas might be associated with the polarizing structures, even though magnetohydrodynamic (MHD) turbulence permeates all different ISM phases (e.g., see Beresnyak & Lazarian 2019).We obtained models for the cold, lukewarm, and warm neutral medium (CNM, LNM, and WNM, respectively), as well as the total atomic hydrogen emission in a region approximately four times larger than the IPS-GI fields.An example of the models centered on each IPS-GI field is presented in Figure 4.The three fields present a total emission peaking principally at low velocities (|V lsr |∼1.6 km s −1 ), denoting that most of the gas is local.Furthermore, the local H I emission peak comprises at least two ISM phases, with the CNM always being present (Figure 4).Details of the Gaussian decomposition performed on each IPS-GI field can be found in Appendix A. The following sections describe our ρ and σ v calculations.All the measured properties are summarized in Table 1.

Hydrogen column density
We estimate the total column density of GASS H I emission in the IPS-GI regions by integrating the brightness temperature over the entire range of velocity channels observed, where dv is the velocity channel width equal to 0.824 km s −1 (see the maps in the Appendix A, Figure 13).Equation ( 8) is only valid under the assumption of optically thin emission (Dickey & Lockman 1990)4 , which is likely to be the case of the diffuse intermediate-latitude IPS-GI fields with volume densities below 30 cm −3 within a cloud with thickness l ef f (see Table 1 and Section 4.2

.2).
The total atomic hydrogen column density observed, N GASS H , is integrated along the entire Galaxy.However, we need to evaluate the column density only in the region where the polarizing dust screen exists (see Figure 1), i.e., very nearby gas (d < 400 pc in fields C2 and C16, and d < 1, 000 pc in C37 ).ROHSA decomposition of H I spectra accommodates for this by distinguishing low-, intermediate-, and high-velocity clouds (LVCs, IVCs, and HVCs, respectively).We expect the gas associated with the polarizing dust structure to be part of an LVC (with |V lsr | < 20 km s −1 , Wakker 2001).Using the LVC model of the H I decomposition, consisting of all the Gaussian components at |V lsr | < 20 km s −1 (see the insets of Figure 4, also see Appendix A), and Equation (8), we calculated the column density map of the nearby gas (N LV C H , Figure 5).Alternatively, we can find an approximation to the total hydrogen column density (i.e., N HI + 2N H 2 ) associated only with the dust, N A V H , by using the column   2022) (Figure 1).The gas number density in the three fields is calculated as where a V (l) is the extinction density (i.e., extinction per distance unit, also known as differential extinction), l c is the depth of the cloud assumed from the extinction density profiles (Figure 1, also see Table 1), and dl = 5 pc is the distance step.
The values calculated in the assumed total depth of the dust, l c .The good agreement suggests a connection between the modeled atomic hydrogen column density and the dust structure within l c .We, therefore, use N LV C H maps presented in Figure 5 to calculate the density of the polarizing cloud.

Volume density
The volume density maps based on the 3D dust map of Leike et al. (2020) were used by Zucker et al. (2021) to construct the average radial density profiles of local molecular clouds (i.e., local structures mostly within and beyond the LB wall, d∼100 − 400 pc).The radial density profiles were fitted with a two-Gaussian function, leading to the characterization of an inner and outer layer with different standard deviations (referred to as "widths") and peak densities for each molecular cloud considered.Zucker et al. (2021) found inner widths between 2.5 − 4.9 pc, with peak densities of 25 − 52 cm −3 , meanwhile, the outer widths are between 8 − 18 pc and densities of 5 − 15 cm −3 .The 3D volume densities reported by Zucker et al. (2021) are likely lower limits due to the systematic uncertainties propagated from the 3D dust map.The uncertainties propagated from the observed spectral rms noise are in the order of 10 17 cm −2 .(d) Turbulent velocity dispersion scaled to the minimum scale probed, ℓips, the scale of the field of view, ℓF oV , and 1 pc scale, ℓ1pc.* An approximated average peak value rather than a maximum value.† Halfway to the total dust pathlength, lc = 1, 000 pc. Assumed as an average between the distances to all the dust peaks observed in Figure 1, bottom.
distance steps of 5 pc.Meanwhile, we used the conversion n H = 880 cm −3 × S x to convert the G-band extinction density of Leike et al. (2020) into volume density; see more details in Zucker et al. (2021).
Figure 6 presents the total volume density profiles of fields C2, C16, and C37.Considering the resolution limitations of the two dust maps compared (i.e., 2 pc 5 for Leike et al. (2020) and 10 pc for Vergely et al. ( 2022)), 5 Despite the claimed map resolution of 1 pc, Leike et al. (2020) consider that only structures at spatial scales greater than 2 pc should be considered reliable; see also discussion from Zucker et al. (2021).
we approximate the total volume density to 27 cm −3 in fields C2 and C16.This is approximately the average peak value between the main dust structures observed in both dust maps compared.On the other hand, C37 presents a series of diffuse structures along a pathlength of 1, 000 pc depth.As it is mentioned in Section 4.1, we have to assume a single structure with average properties, i.e., a n H ∼3 cm −3 (see the black horizontal line in Figure 6, bottom).
Comparing our observations with the 3D models of Zucker et al. (2021), our field C2 lies approximately within the outer envelope modeled for Ophiuchus molec- −1.5 cm −3 .The latter is consistent with the value obtained for C16 from Figure 6.However, it is more uncertain since the Musca molecular cloud is mostly unresolved by Leike et al. (2020) 3D dust map at small radial distances (i.e., < 2 pc, see explanation in Zucker et al. 2021).On the other hand, the profile of C2 in Figure 6 demonstrates that our field of view may be closer to the amplitude of the inner component of Ophiuchus, which is n H = 31.1 +1.6  −1.5 cm −3 .The n H values obtained from 3D dust models include all the phases and molecular hydrogen (i.e., n H = n HI + 2n H 2 ).However, we showed that the atomic hydrogen column density in the nearby ISM of the IPS-GI regions dominates in a large portion compared to that of the H 2 (see the comparison between N A V H and ⟨N LV C H ⟩ values, Section 4.2.1, and Table 1).Therefore, we expect the total volume density to have a contribution primarily from the neutral atomic medium.

Cloud effective thickness
The cloud effective thickness can be obtained through the ratio between the average value of the LVC column density map (Section 4.2.1) and the approximated volume density (Section 4.2.2) as follows, Using Equation ( 10) and the values of Table 1,

Mass density
The local H I emission in our fields comprises almost all the interstellar gas phases (see Figure 4, also see Appendix A), and the dust in the nearby polarizing screen (Figure 1) is consistent with the gas represented by the LVC model (Section 4.2.1).Therefore, we calculate the mass density map of the local clouds using the LVC column density map and the constant cloud thickness as follows, where m H is the mass of the hydrogen atom.Figure 7 presents the mass density maps per IPS-GI field.The average values per IPS-GI field are shown in Table 1.

Turbulent velocity dispersion
At high and intermediate Galactic latitudes, the cold medium is often associated with dust structures.For instance, Lenz et al. (2017) found that, at low column densities (N H < 4×10 20 cm −2 ), H I is in good agreement with dust maps.Hensley et al. (2022) showed that the CNM fraction strongly correlates with the dust fraction in PAHs.Furthermore, the cold medium also traces the magnetic field properties.For instance, the H I smallscale structures, which are preferentially cold gas, are oriented predominantly with the magnetic field (Clark et al. 2015(Clark et al. , 2019;;Clark & Hensley 2019).Moreover, Lei & Clark (2023) found a strong correlation between the CNM mass fraction and thermal dust polarization, implying that the polarization is created in the densest cold regions.A very short pathlength characterizes this cold medium.On the other hand, the WNM seems to contribute very little (about 4%) to the polarization due to depolarization in the warm medium (see Lei & Clark 2023, for more details).Nevertheless, the internal dynamics of the CNM are related to the dynamics in the WNM according to numerical simulations of the local multiphase ISM (Saury et al. 2014).
In this vein, following the approach from Marchal & Miville-Deschênes (2021), we will treat the cold structures (CNM clouds) as test particles embedded within the warm medium (WNM).This ensures that the velocity dispersion among the cold clouds reflects the turbulent velocity field of the entire multiphase medium.Therefore, we calculate the velocity dispersion of the nearby cloud as the second-moment map of the H I emission using only the CNM components of the LVC model obtained in the H I Gaussian decomposition (Ap- pendix A) as follows, where N gauss is the number of Gaussian components, w i = a i σ v,i is the weight, a i , µ i , and σ v,i are the amplitude, central velocity, and velocity dispersion, respectively, of the ith component, and M 1 is the first-moment map of the H I emission defined as The velocity dispersion measured from H I emission comprises a thermal and a turbulent contribution (e.g., see Marchal & Miville-Deschênes 2021) such that, Approximating the turbulent sonic Mach number, M s , as the turbulent to thermal velocity dispersion ratio, we can use Equations ( 14) and ( 15) to calculate the CNM turbulent velocity dispersion map as follows, The turbulent Mach number for diffuse, cold neutral gas is typically around 3, but with wide variations between 1 − 7, according to Heiles & Troland (2003) and Murray et al. (2015).We, therefore, used the median Mach number, M s = 2.86, from Murray et al. (2015) (see their Figure 5), which is also in agreement with the typical value from Heiles & Troland (2003).The final maps of the turbulent dispersion per field are presented in Figure 7, and the average values are shown in Table 1.
The scale at which the CNM turbulent velocity dispersion was measured in the PPV GASS cubes (i.e., along the line of sight) is larger than the largest spatial scale probed by the IPS-GI observations on the plane of the sky.In Sections 3.3.1 and 4.1.1,we argued that the presence of significant large-scale gradients in polarization angle across the IPS fields, corresponding to the regular magnetic field component, provides evidence for a scenario involving small-scale turbulence within the medium.Therefore, we can scale down the gas velocity dispersion σ v,turb to the scales of the IPS-GI fields as where the minimum angular scale probed by the IPS-GI fields, ℓ ips = 0. • 01, is assumed to be the upper limit of the turbulent correlation length (see Section 5.1), corresponding to 0.024−0.030pc at a distance of 136−172 pc in fields C2 and C16, respectively (Table 1).The distances to the clouds result from the median value of the models in Zucker et al. (2021) (see also d peak values in Figure 6).In the same fields, the scale of the cold neutral medium is assumed to be l ef f ∼10 − 15 pc, i.e., the approximated effective thickness of the cloud (Section 4.2.3).Let us recall that the cold medium in diffuse (N H < 10 21 cm −2 ) high-latitude (|b| > 30 • ) regions is expected to have very short pathlengths and produce most of the polarization with high magnetic alignment according to Lei & Clark (2023).Finally, the exponent q in Equation ( 17) depends on the structure in the medium and is in principle unknown.For a Kolmogorov-like spectrum (Kolmogorov 1941;Frisch 1995), and also in the case of anisotropic MHD Alfvénic turbulence in the strong regime (Goldreich & Sridhar 1995), an exponent q = 1/3 is expected, although slightly different values are observed in H I clouds (Larson 1979) or predicted in models for the scales below the inertial range (Boldyrev 2002).We adopt Wolfire et al. (2003) here, who propose that q = 1/3 is the appropriate value for the smaller length scales in the CNM.To probe our parameters in the scaling law of the turbulent velocity field, we calculated the mean velocity dispersion at 1 pc scale, ⟨σ v,turb ⟩ ℓ1pc (Table 1), finding values between ∼0.5 − 1.3 km s −1 consistent with the average values found in cold regions of the ISM, 0.4 − 1.2 km s −1 (Wolfire et al. 2003;Miville-Deschênes et al. 2017).
The scaling factor of the turbulent velocity dispersion in Equation ( 17) is more difficult to determine in C37 due to the complex distribution of dust structures within a long pathlength (l c = 1, 000 pc, Figure 1).However, to obtain an approximated estimation of the ISM properties in C37, we assume a distance to the cloud roughly at 500 pc, halfway to the total pathlength of the dust observed in Figure 6, bottom.This leads to observed scales (ℓ ips ) around 0.087 pc.Nevertheless, we acknowledge the large uncertainties due to the approximations.

RESULTS
In the following sections, we present our results of the ADFs and the magnetic field strength of the intermediate-latitude IPS-GI fields: C2, C16, and C37.

Angular dispersion function
The resulting ADFs are presented in Figure 8, along with the best fit (Equation ( 4) corrected for measurement errors, as described by Equation ( 5)) in the scale bins at which we see a linear increase in dispersion.The scales above 2/3 of the field width are not considered due to subsampling issues (see histograms in the insets in Figure 8, also see Section 4.1.2).Likewise, small bumps of polarization angle dispersion at small scales are not considered significant; see the discussion of Section 4.1.2.
The ADFs of Figure 8 do not present a steep decline or dip approaching ℓ = 0.This suggests that the turbulent correlation length δ is below the smallest scale probed by our observation, as in the case 3.b of Section 3.3.1.Therefore, the ADFs' moderate slope is consistent with fluctuations from large-scale structures partially within the fields of view (see explanation in Section 4.1.1).
All three fields show an almost linear increasing polarization angle dispersion throughout a large range of angular scales measured, i.e., C16 and C37 between ∼0.• 01 and ∼0.• 2, which correspond to scales between ∼0.03 − 0.6 pc (for C16 ) and ∼0.09 − 1.7 pc (for C37 ), and C2 between ∼0.• 01 and ∼0.• 3, corresponding to ∼0.02 − 0.7 pc, assuming the dust peak distance of Table 1.In field C37, the accuracy of the ADF holds less.The dust in C37, which is distributed within multiple diffuse structures along 1, 000 pc, may have actual scales between 0.01 − 2 pc.However, as mentioned in Section 4.2.3, we assumed a single polarizing screen of 120 pc thick, at an average distance of 500 pc, for simplicity.
The values of the turbulent parameter, b, obtained from fitting the ADFs, are presented in Table 1.We used b to calculate the turbulent to large-scale magnetic field ratios, ⟨B 2 t ⟩ 1/2 /B pos , as in Equation ( 6).The results presented in Table 1 are around ∼0.1 for all fields studied.Therefore, highly regular magnetic fields dominate the IPS-GI intermediate-latitude observations, as suggested in Angarita et al. (2023).This is further confirmed by comparing optical and thermal dust polarization with diffuse dust emission at 12 µm in Section 6.2.2.

Magnetic field strength maps
The average plane-of-sky magnetic field strength calculated with the DCF, ST-DCF, and ADF-DCF methods are presented in Table 2 along with the values obtained with Cho & Yoo (2016) correction to the DCF method (CY-DCF), which are discussed in Section 6.2.1.We used the mass density 2D map, ρ LV C , estimated from the LVC components of the neutral medium (Figure 7,left).Furthermore, assuming a scenario in which the observed scale lengths are below the turbulent correlation length, i.e., δ > ℓ max , we used the DCF and ST-DCF methods with the polarization angle dispersion corrected by the measurement errors, σ θ,corr , and the CNM turbulent velocity dispersion 2D map of the LVC components scaled down to the spatial scale of the field of view, i.e., (σ v,turb ) F oV .In the case where the turbulence is at small scales (see Sections 3.3.1 and 4.1.1),we used the ADF-DCF method with the turbulent parameter b from the fit to the ADFs (Figure 8) and the CNM turbulent velocity dispersion 2D map scaled to the smallest scale probed by the ADFs, i.e., (σ v,turb ) ips (Figure 7, right).See the average values of the 2D maps in Table 1.
The DCF and ST-DCF methods give magnetic field strengths of ∼16 µG and ∼7 µG in C2, ∼9 µG and ∼3 µG in C16, and ∼6 µG and ∼3 µG in C37, respectively.The classic DCF method yields higher values than the ST-DCF method by an approximated factor of 2 in C2 and C37, and 3 in C16.However, the most probable situation is one where the dispersion is mainly attributed to large-scale field variations, and the turbulent correlation length is at small scales.
Figure 9 presents the plane-of-sky magnetic field strength maps calculated with the ADF-DCF method in each IPS-GI field considered.The average B pos ob-Table 2.
Average magnetic field strength in the intermediate-latitude IPS-GI fields.tained with the ADF-DCF method is ∼9 µG, ∼5 µG, and ∼3 µG, for field C2, C16, and C37, respectively (Table 2).The results of the ADF-DCF method are closer to those of the ST-DCF.Nevertheless, the direct comparison is meaningless since both methods are valid for different scenarios of the turbulent medium.
The ADF-DCF method was implemented in this research under the assumption that the observed scales are above the turbulent correlation length, i.e., δ < ℓ min .This is supported by the polarization angle gradients observed in Figure 2 (Section 4.1.1),showing a smooth transition in polarization angle due to a large structure partially within the IPS-GI fields, and the profiles of the ADFs (Figure 8), which do not present a steep dip at the smallest scales.In particular, we highlight the similarities between the morphology of the median polarization angle gradients (Figure 2) and the features observed in the B pos map (Figure 9, obtained from the combination of the ρ and σ v,turb maps of Figure 7).However, due to the limited size of the fields of view, we cannot be completely sure that the turbulence scales are below the observable scales.In the case δ > ℓ max , we recommend considering the results of the classic DCF (Section 3.1) or the ST-DCF (Section 3.2) methods.
Field C37 presents the lowest plane-of-sky magnetic field strength among the sample.Nonetheless, this field has the greatest uncertainties in estimations, primarily because of the numerous assumptions made given the intricate nature of very diffuse dust structures spanning a considerable pathlength.On the other hand, the average values obtained for fields C2 and C16 are in good agreement with previous estimations in the LB wall and the nearby diffuse ISM; see the discussion in Section 6.2.1.

DISCUSSION
We discuss below the ISM properties estimated in the intermediate-latitude IPS-GI fields and compare them with previous estimations in the local, diffuse ISM.Finally, we describe the properties of the magnetic field in the local ISM as observed by optical and thermal dust polarization.

Turbulent scales
We proposed that the ADFs of Figure 8 cannot resolve the turbulent correlation length δ in the three regions studied (Section 5.1).However, it is possible to establish an upper limit at the smallest scale we can demonstrate, i.e., 0. • 01.Assuming the distance to the highest dust peak along the pathlengths of C2 and C16, i.e., d peak approximately behind or within the LB wall, the smallest scale corresponds to ∼24 − 30 mpc.Due to the lack of measurements of this kind in diffuse nearby structures (i.e., with n H < 30 cm −3 ), our results can only be compared to studies in nearby molecular clouds.
The turbulent correlation length was measured with the ADFs method at mpc scales in various molecular clouds.For instance, Houde et al. (2009) found a turbulent correlation length of 16 mpc in OMC-1, assuming a distance of 450 pc.Using starlight polarization, Franco et al. (2010) found δ around 2−5 mpc in different regions of the Pipe nebula assuming d = 160 pc.Wang et al. (2019) used submillimeter polarization in the IC 5146 molecular cloud to calculate the ADF and found a turbulent scale of ∼60 mpc adopting the distance of 813 pc.
Although all the above measurements have been made on dense sources, i.e., n H between 10 3 −10 5 cm −3 (Alves et al. 2008;Karoly et al. 2023;Hwang et al. 2023), the upper limits of the turbulent correlation length in the three IPS-GI diffuse fields are on the same order of magnitude, at mpc scales.This demonstrates significant similarities in the turbulent correlation lengths at low and high densities in the local ISM (d < 1, 000 pc).

Galactic magnetic field in the Solar neighborhood
The intermediate-latitude IPS-GI polarization measurements probe the magnetized, diffuse medium in the nearby ISM (d < 400 pc for C2 and C16, and d < 1, 000 pc for C37, Figure 1).The following sections discuss the consistency of the magnetic field strength results with the values measured in the Solar neighborhood and the estimations made with different methods.Furthermore, we discuss the structure of the magnetic field observed in the plane of the sky with optical polarization and thermal dust polarization and the correlation between the large-scale magnetic field, the diffuse dust emission observed at 12 µm, and the thermal dust emission.

Magnetic field strength
Regardless of the method used, our results (Table 2) are consistent with the Galactic magnetic field strength estimations in the local diffuse ISM according to observations of starlight polarization and the Zeeman effect.For instance, Heiles & Troland (2005) found a median total magnetic field of 6 ± 1.8 µG in the CNM.Crutcher et al. (2010) and Crutcher (2012) presented a maximum line-of-sight magnetic field strength of ∼10 µG in low-density regions of molecular clouds, with typically n H ∼300 cm −3 , using a Bayesian analysis of molecular Zeeman measurements of H I, OH, and CN.Similarly, near-IR starlight polarization observations of the Taurus dark-cloud complex (Chapman et al. 2011) showed plane-of-sky field strengths of ∼10 µG in low-density off-cloud regions (n H ∼100 cm −3 ).Later on, Medan & Andersson (2019) found an average plane-of-sky field strength in the LB wall of ∼8 µG, assuming a distance to the wall of d ≃ 100−315 pc and using optical starlight polarization of OB associations at b > 30 • .
One fundamental assumption of DCF-like methods is the sub-Alfvénic turbulence, i.e., M A < 1, where M A is the Alfvén Mach number.We calculated M A = ⟨σ v,turb ⟩ F oV /V A , where V A = B pos / √ 4πρ is the Alfvén velocity calculated from the plane-of-sky magnetic field strength.In regions C2, C16, and C37, we obtained M A values of 0.18, 0.13, and 0.24, respectively, using parameters from the classic DCF method.Additionally, we calculated M A = σ θ /f , where f is the classical DCF correction factor, yielding M A < 1 for f values within the range [0.1,1].Furthermore, M A can also be estimated using the turbulent-to-total magnetic field ratio obtained from the ADF fit parameters, i.e., M A = ⟨B 2 t ⟩ 1/2 /B pos , which is ∼0.1 in all cases, reaffirming the sub-Alfvénic condition in all three regions.
The classical DCF method does not account for the contributions of multiple turbulent cells along the line of sight.If the number of turbulent cells is large, B pos is likely overestimated.Cho & Yoo (2016) proposed a correction to this issue by using the relation , where δV c is the standard deviation of the line emission centroid velocities observed towards different lines of sight, and N c is the number of turbulent cells.We computed the average V c in each pixel of our CNM model maps using the parameters of the Gaussian LVC components as in Equation ( 13).Subsequently, we estimated δV c using Equation (18) of Cho & Yoo 2016 for each IPS-GI region.Then, we recalculated B DCF pos , hereinafter B CY −DCF pos , using the correction factor ξ = 0.5 for direct comparison with the classical DCF.We found average B CY −DCF pos values of ∼3 µG, ∼7 µG, and ∼1 µG in C2, C16, and C37, respectively.This suggests that DCF results are significantly overestimated in the regions of C2 and C37, whereas for C16, the correction is smaller.
The magnetic field strength calculation depends considerably on knowing the ISM properties very well.In particular, the volume density assumed from 3D dust maps models (Section 4.2.2) is a key parameter in our cloud thickness calculation, which in turn affects the determination of important parameters such as ρ (Section 4.2.4) and σ v,turb (Section 4.2.5).Finding volume densities close to 30 cm −3 in the IPS-GI fields implies that the density, calculated using Equation ( 11), is dom-inated by the cold medium (Draine 2011).However, the H I line emission also shows contributions from the WNM and LNM in the three fields.As we cannot ignore the contributions from the warm medium, the turbulent velocity dispersion, as defined in Section 4.2.5, characterizes the dispersion of the entire multiphase interstellar gas.Consequently, the results are a weighted mixture of ISM phases that adds uncertainty.
Figure 10 presents the magnetic field strength as a function of the cloud thickness, colored by the total volume density, for each IPS-GI field.To construct the functions of each panel in Figure 10, we used the parameters of Table 1 and varied the n H values between 1 − 200 cm −3 (Wolfire et al. 2003;Murray et al. 2015).Let us recall that we are assuming that n H is mainly made up of neutral atomic hydrogen.This is supported by the good agreement between N A V H and ⟨N LV C H ⟩. Furthermore, according to Murray et al. (2015), the n H of the CNM only is expected to be within 5 − 120 cm −3 .
The gray shaded area in Figure 10 shows low volume densities (< 5 cm −3 ) that are unlikely to represent our C2 and C16 fields.On the other hand, the volume density and, thus, the magnetic field in C37 are likely underestimated since we assumed an average peak value rather than the maximum n H .However, the density profiles of Figure 6 (bottom panel) shows n H typically below ∼5 − 7 cm −3 .Furthermore, the red shaded area represents the values obtained for volume densities above 30 cm −3 , which is a good upper limit for the diffuse cold medium (Draine 2011).Zucker et al. (2021) argued that any estimation of n H from Leike et al. (2020) 3D dust map -especially in dense regions -should be considered a lower limit due to the systematic errors of the technique and the resolution limitations of the 3D map.Nevertheless, since our diffuse IPS-GI fields are far off from the center of the molecular clouds, n H may not be larger than 30 cm −3 , i.e., this is also the approximated maximum value of Vergely et al. (2022) dust profiles in C2 and C16 (Figure 6).Therefore, we would not expect our fields to be within the red area either.
The mentioned boundaries enable us to establish approximate values for the lower and upper limits of the magnetic field strength within our fields.The maximum B pos are ∼10 µG, ∼5 µG, and ∼7 µG for C2, C16, and C37, respectively.The minimum B pos are ∼2 µG and ∼1 µG for C2 and C16, respectively.All the calculations were made with the ADF-DCF method.6.2.2.Plane-of-sky magnetic field structure in the local ISM IPS-GI starlight polarization in intermediate-latitude fields probes the magnetic field coupled to nearby dust structures behind and within the LB wall (d < 400 − 1, 000 pc, depending on the field, Figure 1).Figure 11 shows the orientation of the plane-of-sky magnetic field as derived from the Planck Collaboration XII et al. (2020) thermal dust polarization (the black vectors) and the IPS-GI optical polarization (the white vectors) in the surroundings of the studied regions.The thermal dust polarization shows a very regular field towards our diffuse regions with an orientation in excellent agreement with starlight polarization.The consistency with optical starlight polarization also indicates that thermal dust emission in these three regions may come from the nearby ISM.
Furthermore, the high-resolution diffuse dust at 12 µm from WISE (Meisner & Finkbeiner 2014) shows fine structures in the sky 1.5 • around the IPS-GI fields (see the background of Figure 11).Although Planck's thermal dust emission, N τ353 H (the contours), is in good agreement with the diffuse dust, the extremely regular large-scale magnetic field seems to not care about small, diffuse structures and, instead, is highly ordered (Figure 11).Hence, optical starlight and thermal dust polarization probe the large-scale magnetic field, which is highly regular in the surroundings and within the intermediate-latitude IPS-GI fields, in the few hundreds of parsecs within and beyond the LB wall.
This result is consistent with our interpretation of the turbulence in the three IPS-GI fields, i.e., the case 3.b in Section 3.3.1.The ADFs (Figure 8) are tracing the dispersion due to large-scale field fluctuations rather than pure turbulence.Consequently, the turbulent correlation lengths are below the smallest spatial scales observed in the three fields., and C37, using the IPS-GI nomenclature.These fields sample diffuse ISM properties.To this end, we implemented the DCF (Davis 1951;Chandrasekhar & Fermi 1953), the ST-DCF (Skalidis et al. 2018), and the ADF-DCF (Hildebrand et al. 2009) methods to estimate the plane-of-sky magnetic field strength under different regimes of the turbulence.A correction to the DCF method that accounts for the number of turbulent cells along the line of sight (i.e., the CY-DCF method) is also discussed.
The polarization angle dispersion of the three fields studied exhibits a steady increase with angular scale, which may be caused by large-scale field fluctuations with angular scales not smaller than ∼0.• 2 (i.e., 0.6 pc at d∼172 pc).Moreover, the three regions show a smooth change in the polarization orientation, probably due to structures larger than the observed region but partially within it.The above evidence suggests that the turbulent correlation length, δ, is below scales of < 0. • 01 (i.e., < 0.03 pc).The neutral H I emission observed by GASS in the IPS-GI regions is associated with the polarizing dust structures in the nearby ISM (ℓ c = 400 − 1, 000 pc).We thereby produced 2D maps of the physical properties of the nearby clouds, such as ρ and σ v,turb , using models derived from a Gaussian decomposition of the GASS H I line emission towards the IPS-GI fields.The estimation of n H is critical in our case for the accuracy of DCF-like methods.Therefore, we evaluate how the B pos values vary in each observed region depending on the volume density constrained by 3D dust models.
We estimated the plane-of-sky magnetic field strength between ∼3 µG and ∼9 µG with the ADF-DCF method, between ∼3 µG and ∼7 µG with the ST-DCF method, and between ∼1 µG and ∼7 µG with the CY-DCF correction to the classic DCF method in the studied regions.These values align with previous calculations using the Zeeman effect, near-IR, and optical starlight polarization in low-density regions of the local ISM.Furthermore, the turbulent to large-scale magnetic field ratios, ⟨B 2 t ⟩ 1/2 /B pos , is ∼0.1 in all cases, indicating the dominance of the regular component of the magnetic field in the nearby ISM.
Finally, thermal dust polarization at 353 GHz from Planck is in good agreement with the optical starlight polarization observed by the IPS-GI in the three regions studied, which originates in the nearby ISM (d < 400 − 1, 000 pc, depending on the field).Therefore, Planck's polarization should be tracing the local ISM in these three regions as well.The thermal dust and optical polarization show a highly regular plane-of-sky magnetic field orientation despite the fine structures observed in the diffuse dust emission at 12 µm.The latter confirms our conclusion about the low turbulent to large-scale magnetic field ratios.

ACKNOWLEDGMENTS
The authors are grateful to the anonymous referee, whose insightful review contributed to the improvement of this paper.This research has used data, tools, and materials developed as part of the EXPLORE project that has received funding from the European Union's Horizon 2020 research and innovation programme under grant agreement No 101004214.Finally, this work has made use of data from the European Space Agency (ESA) mission Gaia (https://www.cosmos.esa.int/gaia),processed by the Gaia Data Processing and Analysis Consortium (DPAC, https://www.cosmos.esa.int/web/gaia/dpac/consortium).Funding for the DPAC has been provided by national institutions, in particular, the institutions participating in the Gaia Multilateral Agreement.

APPENDIX
A. H I DECOMPOSITION ROHSA (Marchal et al. 2019) is an algorithm that accounts for the spatial coherence of the emission and the multiphase nature of the gas.It performs a non-linear Gaussian regression in which λ a , λ µ , and λ σ are the constant hyper-parameters for the amplitude, a n > 0, the position, µ n , and the standard deviation, σ n , positionposition (2D) maps of the n-th Gaussian component.The hyper-parameters tune the regulation term, R(θ) (with θ n being the set of fit parameters) that penalizes small fluctuations of the parameters in the cost function, Q(θ), i.e., the reduced χ 2 2D map.A fourth term related to the width of the Gaussian components, and thus to the gas thermodynamics, is included in the cost function and parametrized by λ ′ σ .This term would give a higher weight to any solution that could describe any of the Cold, Lukewarm, and Warm Neutral Medium phases (CNM, LNM, and WNM; see Marchal et al. 2019, for more details).
Using the H I data cubes from GASS (Section 2.2) and the set of parameters displayed in Table 3, we performed the Gaussian decomposition on each IPS-GI field studied (Section 4).The sets of parameters were carefully chosen (through a trial-and-error method) in each case to ensure coherent structures in the dispersion-velocity (σ v −v) space weighted by the fraction of the total emission of each Gaussian, as in Figure 16 of Marchal et al. (2019); see Figure 12, right column.The dispersionvelocity heat map is an important metric that allows us to assess the consistency of the Gaussian components and the classification of the same into the different thermal states using the velocity dispersion thresholds (≲ 3 km s −1 for the CNM, 3 − 7 km s −1 for the LNM, and ≳ 7 km s −1 for the WNM; Takakubo 1967;Marchal et al. 2019).So, we can distinguish between the components of the different ISM phases.
To determine the best number of Gaussian components and evaluate the quality of the modeling, we fixed the four hyper-parameters (λ a , λ µ , λ σ , and λ ′ σ ; Table 3) and varied the number of Gaussian components, N Gauss , from 1 to 15 (21 in the case of C37 ).The average of the reduced χ 2 2D map (e.g., Figure 12, middle column) as a function of the number of parameters (Figure 12, left column) helped us to find the N Gauss at which ⟨χ 2 r ⟩ converges towards a value close to one.The final N Gauss per IPS-GI field is shown in Table 3 and the respective ⟨χ 2 r ⟩ is displayed in Figure 12, left column.The final models of each of the ISM phases (CNM, LNM, and WNM) separately are presented in Figure 4, together with the H I emission line and the total model.We also show in the insets in each panel the models made up of the LVC components, with |V lsr | < 20 km s −1 , which are expected to be those associated with the (2) the number of Gaussian components considered for modeling; (3) ROHSA hyper-parameters for an, µn, and σn 2D maps; (4) ROHSA hyper-parameter included in the cost function that accounts for the gas thermodynamics; (5) average reduced χ 2 in the IPS-GI regions (Figure 12).
nearby structure under consideration.The H I emission line in C2 is the simplest.All Gaussian components are centered around V lsr ≈ 1.6 km s −1 , indicative of an LVC.Meanwhile, the H I spectrum of C16 shows intermediate-and high-velocity clouds (IVC and HVC, respectively) that are unlikely related to the polarizing dust structure in the nearby ISM (Figure 1).Similarly, the spectrum of C37 presents IVC and HVC components that are not considered in the study of the local structure.
Figure 13 shows the column density of the observations compared to the modeled column density found with ROHSA and calculated with Equation (8).The residual maps show a good agreement between the observations and the total models with a difference below ∼1%.The CNM, LNM, and WNM models are also presented in Figure 13 with column densities in the order of 10 20 cm −2 .Besides the complexity of C37 's H I spectrum (as well as the A V density profile with distance, Figure 1) and the long pathlength integrated (l c = 1, 000 pc), the column density is comparable with the other fields and the local H I structures are likely cold.

Figure 1 .
Figure 1.Extinction density as a function of the distance retrieved from Vergely et al. (2022) at 10 pc resolution in units of mag pc −1 (the blue solid line), Leike et al. (2020) at 1 pc resolution in units of e-foldings pc −1 (the orange solid line), and Edenhofer et al. (2023) with resolution between 0.4 − 7 pc in the units of the parameter E from Zhang et al. (2023) (the gray solid line) in direction of C2 (top), C16 (middle), and C37 (bottom).The gray dots show the IPS-GI optical polarization as a function of the distance.The vertical black dashed lines show the distance to the dust peaks taken from Zucker et al. (2021).The vertical black solid lines are the total dust depth assumed, lc.

Figure 2 .
Figure 2. Distribution of the median polarization angle (left column) and the standard error (right column) in the field-of-view of C2 (top row), C16 (middle row), and C37 (bottom row).

Figure 3 .
Figure 3. Optical polarization vectors map (left column) and polarization angle as a function of distance (right column) of fields C2 (top row), C16 (middle row), and C37 (bottom row).The red vectors represent unreliable polarization angle measurements above 3σ from the mean distribution and are excluded from the calculation of the ADF.

Figure 4 .
Figure 4. H I emission line of GASS (red solid line) compared to ROHSA's total model (black solid line) and the CNM (blue solid line), LNM (green solid line), and WNM (orange solid line) models in the center of the map.The insets in each panel present the model built with the components of the low-velocity clouds (LVCs, |V lsr | < 20 km s −1 ) from the ROHSA decomposition.

Figure 5 .
Figure 5. Gas column density of the LVC model obtained from ROHSA decomposition in field C2 (left), C16 (middle), and C37 (right).The white pseudo-vectors represent the IPS-GI optical polarization measurements for comparison.density to reddening ratio, N H /E(B −V ), fromBohlin et al. (1978), the total to selective extinction ratio of A V /E(B −V ) = 3.1(Savage & Mathis 1979;Fitzpatrick 2004), and the extinction density profiles fromVergely et al. (2022) (Figure1).The gas number density in the three fields is calculated as Following the basic idea from Zucker et al. (2021), we used the extinction density profiles from Vergely et al. (2022) and Leike et al. (2020) (see Figure 1) to estimate the peak volume density of the local dust structure found along our lines of sight.We followed the same conventions as in Section 4.2.1 to convert the extinction density of Vergely et al. (2022) into volume density considering Note-(a)The uncertainties propagated from the turbulent parameter errors are of the order of 10 −3 or below.(b) Integrated up to the cloud depth lc, the uncertainty is propagated from the error ofVergely et al. 2022 models.(c)

Figure 6 .
Figure 6.Volume density as a function of the distance obtained from Leike et al. (2020) (the gray and orange dots), and Vergely et al. (2022) (the blue dots) 3D dust maps in C2 (top), C16 (middle), and C37 (bottom).The orange solid line is the mean smoothed curve, with kernel 4. The gray bars are the standard deviation.The blue shaded curve is the error from Vergely et al. (2022).The horizontal black line marks the approximated peak volume density.The vertical black dashed line is the distance to the cloud in C2 and C16, which is assumed from Zucker et al. (2021).In the case of C37, the vertical red dashed line is the mean distance of the entire dust distribution.
we found an effective thickness of approximately 15 pc, 10 pc, and 127 pc for C2, C16, and C37, respectively.The first two values are consistent with the low widths modeled by Zucker et al. (2021) in the local clouds (i.e., Ophiuchus and Musca molecular cloud).Since the IPS-GI pathlengths observed are actual segments that cross the 3D structures modeled by Zucker et al. (2021), we expect the effective thickness to be around two times the widths modeled.The values obtained with Equation (10) are within the widths of the highest dust peaks found along the line of sight in the extinction density profiles of Vergely et al. (2022), Leike et al. (2020), and Edenhofer et al. (2023) (Figure 1 and Figure 6).Higher resolution 3D dust maps (e.g., Leike et al. 2020 and Edenhofer et al. 2023) are expected to show smaller spatial scales in the ISM structures than lower resolution maps (e.g., Vergely et al. 2022).Therefore, the Leike et al. (2020) 3D dust map yields lower l ef f values through Zucker et al. (2021) method than those that could be estimated with Vergely et al. (2022) or Edenhofer et al. (2023) 3D dust maps.

Figure 8 .
Figure 8. Dispersion of the polarization angle corrected by the measurement errors as a function of the scale length in field C2 (left), C16 (middle), and C37 (right).The red solid line is the best fit to the scales at which we see a linear increase in dispersion.The insets show a log-scale histogram depicting the number of sampled star pairs per angular scale bin.

Figure 10 .
Figure 10.Magnetic field strength calculated with the DCF (the circles), ST-DCF (the triangles), and ADF-DCF (the squares) methods as a function of the effective thickness, colored by the total volume density.The red and orange horizontal dashed lines represent our Bpos results calculated with ADF-DCF and ST-DCF, respectively.The vertical black dot-dashed line shows the effective thickness of the cloud obtained with the values ofTable 1.
7. SUMMARY We used optical starlight polarization measurements from the IPS-GI catalog (Versteeg et al. 2023), 3D dust maps from Vergely et al. (2022) and Leike et al. (2020), and PPV cubes of H I emission from GASS survey (McClure-Griffiths et al. 2009) to characterize the planeof-sky magnetic field strength and the properties of the local ISM towards three intermediate Galactic latitude fields: C2, C16

Figure 11 .
Figure 11.Left: Diffuse dust emission at 12 µm from WISE (Background) compared with the thermal dust column density (the contours) and polarization (the black pseudo-vectors) at 353 GHz from Planck in the vicinity of fields C2 (top row), C16 (middle row), and C37 (bottom row) marked by the white rectangle in the center.Right: Comparison between the optical polarization (the white pseudo-vectors) and the thermal dust polarization at 353 GHz (the black pseudo-vectors).The orientation of the thermal dust polarization was rotated by 90 • so that both Planck polarization and optical polarization pseudo-vectors indicate the orientation of the plane-of-sky magnetic field.

Figure 12 .
Figure 12.Left column: average reduced χ 2 as a function of the number of Gaussian parameters (Three per Gaussian component: an, µn, and σn).The mean χ 2 r as a function of the number of parameters in the IPS-GI field is shown by the gray line.The red dot represents the mean χ 2 r at the NGauss selected.Middle column: reduced χ 2 map of the total model built with Table3parameters, with the position of IPS-GI measurements overlaid in white.Right column: dispersion-velocity probability distribution function weighted by the fraction of total emission of each Gaussian component (as inMarchal et al.  2019, Figure 16).

Figure 13 .
Figure 13.H I column density maps per intermediate-latitude IPS-GI field in units of 10 21 cm −2 .From left to right: N GASS H and ⟨N LV C

Table 1 .
Properties of the intermediate-latitude IPS-GI fields.
The authors acknowledge the Interstellar Institute's program "II6" and the Paris-Saclay University's Institut Pascal for hosting discussions that nourished the development of the ideas behind this work.Y.A. acknowledges Dr. S. E. Clark, M. Lei (M.Sc.), and Dr. C. Zucker for their contributions to insightful discussions about the properties of the nearby ISM.Over the years, IPS data have been gathered by a number of dedicated observers, to whom the authors are very grateful: Flaviane C. F. Benedito, Alex Carciofi, Cassia Fernandez, Tibério Ferrari, Livia S. C. A. Ferreira, Viviana S. Gabriel, Aiara Lobo-Gomes, Luciana de Matos, Rocio Melgarejo, Antonio Pereyra, Nadili Ribeiro, Marcelo Rubinho, Daiane B. Seriacopi, Fernando Silva, Rodolfo Valentim, and Aline Vidotto.Y.A., M.H., and M.J.F.V. acknowledge funding from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation programme (grant agreement No 772663).C.V.R. acknowledges support from Conselho Nacional de Desenvolvimento Científico e Tecnológico -CNPq (Brazil) (Proc.310930/2021-9).A.M.M.'s work and optical/NIR polarimetry at IAG have been supported over the years by several grants from São Paulo state funding agency FAPESP, especially 01/12589-1 and 10/19694-4.A.M.M. has also been partially supported by the Brazilian agency CNPq (grant 310506/2015-8).A.M.M. graduate students have been provided with grants over the years from the Brazilian agency CAPES.

Table 3 .
Input parameters used in ROHSA per IPS-GI field.