Scattered Light Imaging: Resolving the substructure of nerve fiber crossings in whole brain sections with micrometer resolution

For developing a detailed network model of the brain based on image reconstructions, it is necessary to spatially resolve crossing nerve fibers. The accuracy hereby depends on many factors, including the spatial resolution of the imaging technique. 3D Polarized Light Imaging (3D-PLI) allows the three-dimensional reconstruction of nerve fiber pathways in whole brain sections with micrometer resolution, but leaves uncertainties in pixels containing crossing fibers. Here we introduce Scattered Light Imaging (SLI) to resolve the substructure of nerve fiber crossings. The measurement is performed on the same unstained histological brain sections as in 3D-PLI. By illuminating the brain sections from different angles and measuring the transmitted (scattered) light under normal incidence, SLI provides information about the underlying nerve fiber structure. A fully automated evaluation of the resulting light intensity profiles has been developed, which allows the user to extract the individual fiber orientations in regions with crossing nerve fibers, for each image pixel at once. We validate the reconstructed nerve fiber orientations against results from previous simulation studies, scatterometry measurements, and fiber orientations obtained from 3D-PLI. We demonstrate in different brain samples (human optic tracts, vervet monkey brain, rat brain) that the fiber orientations can be reliably reconstructed for up to three crossing nerve fiber bundles in each image pixel with a resolution of up to 6.5 $\mu$m. We show that SLI also yields reliable fiber directions in brain regions with low 3D-PLI signals coming from regions with a low density of myelinated nerve fibers or out-of-plane fibers. In combination with 3D-PLI, the technique can be used for a full reconstruction of the three-dimensional nerve fiber architecture in the brain.


I. INTRODUCTION
The human brain consists of about 100 billion neurons [1]. To this day, neuroscientists are working on disentangling this gigantic, densely grown network of nerve fibers. So far, diffusion magnetic resonance imaging (dMRI) is the only possibility to map nerve fiber pathways in-vivo. In postmortem brains, dMRI achieves resolutions down to a few hundred micrometers [2,3]. At the same time, many nerve fibers in the white matter have diameters in the order of 1 µm [4,5], which creates significant challenges in reconstructing fiber pathways at the level of single nerve fibers, and solving the problem of crossing nerve fibers within the given spatial resolution. Due to an insufficient knowledge about nerve fiber crossings, fiber tractography algorithms may yield false-positive nerve fiber pathways [6]. The microscopy technique 3D-Polarized Light Imaging (3D-PLI) represents a powerful optical approach to determine the three-dimensional course of nerve fibers in whole, unstained histological brain sections with micrometer resolution [7,8]. The technique has been applied, for example, to reveal previously unknown insights into the connectional architecture of the hippocampal region [9], and it served as cross-validation for fiber tractography algorithms and to improve the interpretation of clinical dMRI data [9,10]. The in-plane resolution reaches a pixel size down to 1.33 µm and allows the mapping of single nerve fiber pathways. Although the course of crossing nerve fibers can be visually tracked in many regions due to a high contrast of the 3D-PLI images, individual nerve fiber orientations cannot be automatically extracted. 3D-PLI yields a single fiber orientation for each image pixel, which creates uncertainties in fiber orientation if the 60 µm thick brain section at this point is comprised of several crossing nerve fibers with different orientations.
Simulation studies of our group [11] have shown that the scattering pattern of light transmitted through artificial nerve fiber constellations reveals the individual nerve fiber orientations in regions with crossing fibers. [12] used coherent Fourier scatterometry with non-focused, normally incident laser light to measure these scattering patterns for real brain tissue samples. These results validate the simulation studies and show that it is indeed possible to determine the individual orientations of crossing nerve fibers from the measured scattering pattern. Although the scatterometry measurement yields a full scattering pattern of the measured brain region (limited only by the numerical aperture of the objective lens), it is only possible to measure a rather small brain region (0.1-1 mm) at a time and the resolution is limited by the laser beam diameter (≥ 100 µm). Along with the simulation results, [11] introduced a scattering measurement where the brain section is illuminated under oblique incidence of light from different angles. This technique can measure only a limited number of scattering angles (along a circle in the full scattering pattern), but it can be used to study the scattering behavior of a whole brain section at once with micrometer resolution. For every measured image pixel, the measurement yields a light intensity profile, which is characteristic for the brain tissue structure at this point. Because of the limited number of measured angles, the line profiles are discretized. Nevertheless, the preliminary measurements showed that the peaks in the line profiles can be used to determine, e. g., the fiber orientations in regions with crossing fibers with an accuracy of 15 • . The measurements were manually evaluated for a small number of selected brain regions (10 × 10 pixels).
Here, we further developed this scattering measurement and introduce Scattered Light Imaging (SLI)-an imaging technique that allows the simultaneous extraction of up to three different (crossing) nerve fiber orientations for whole brain sections with micrometer resolution. We developed the open-source software SLIX (Scattered Light Imaging ToolboX) that allows an automated evaluation of the measurement, and the computation of different parameter maps containing various tissue information. The resulting fiber direction maps can be used, for example, as input for fiber tractography algorithms. The studies were performed in regions of interest of brains of different species (rodent, vervet, and human). We developed a correction procedure that enables the determination of fiber orientations with an accuracy of ±2.4 • for 15 • -discretization of the measured line profiles, keeping the number of required measurements low. By comparing the measured SLI profiles to results obtained from scatterometry measurements and to simulated scattering patterns, we validate our results. We show that the determined fiber orientations are compatible with the fiber orientations obtained from 3D-PLI (in non-crossing regions), and that the SLI measurement yields reliable fiber directions also in other regions with low 3D-PLI signals, i. e. regions with a low density of myelinated nerve fibers and some regions with out-of-plane fibers. Hence, a combined measurement of 3D-PLI and SLI allows a more correct reconstruction of nerve fiber orientations in the brain, especially in regions with crossing nerve fibers.

A. Preparation of brain sections
The measurements were performed on brain sections obtained from a Wistar rat (male, 3 months old), two vervet monkeys (male, 1 and 2.4 years old), and a human (female, 74 years old). The brains were removed from the skull within 24 hours after death in accordance with legal and ethical requirements. 1 They were fixed in a buffered solution of 4 % formaldehyde for several weeks, immersed in 20 % glycerin, deeply frozen, and cut with a cryostat microtome (Polycut CM 3500, Leica Microsystems, Germany). The rat brain and one vervet brain were sectioned coronally, the other vervet brain was sectioned sagittally, into 60 µm thin sections each. The human optic chiasm (including at least 1 cm of the optic nerves and 1 cm of the optic tracts) was removed from the human brain and cut along the fiber tracts of the visual pathway into sections of 60 µm and 30 µm. To study the behavior of crossing fiber bundles with welldefined crossing angles, the sections of the optic chiasm were split into two parts at the median line and the sections of the optic tracts were manually placed on top of each other under different crossing angles (cf. Fig. 7(c)). Each brain sample was mounted onto a glass slide, embedded in a solution of 20 % glycerin, cover-slipped, and measured up to one day (3D-PLI) or several days (SLI) after tissue embedding. To increase tissue scattering, some brain sections were revitalized (freshly embedded in glycerin) after several months. Table I in Appx. C lists the section thickness and the days after tissue embedding/revitalization for each measured brain section.

B. 3D Polarized Light Imaging (3D-PLI)
The 3D-PLI measurements were performed with a polarizing microscope (LMP-1, Taorad GmbH, Germany), see Fig. 1(a). The sample is illuminated by linearly polarized, incoherent light with a wavelength of about 550 nm. The transmitted light passes a circular analyzer (quarter-wave plate and polarizer) and is recorded by a CCD camera (Qimaging Retiga 4000R). During the measurement, the polarizer in front of the sample is rotated by ρ = {0 • , 10 • , . . . , 170 • }. For each rotation angle ρ, the camera records an image of the transmitted light intensity, yielding a series of images and a sinusoidal intensity profile per image pixel (see Fig. 1(b)). The transmittance (average transmitted light intensity over all rotation angles) is a measure of the tissue attenuation; strong sources of scattering and absorption, such as myelinated nerve fibers, appear dark. The phase of the signal indicates the in-plane orientation angle of the nerve fibers (fiber direction ϕ). The amplitude of the signal (retardation | sin δ| ∼ cos 2 α) is a measure of the birefringence strength and decreases with increasing out-of-plane angle of the enclosed nerve fibers (fiber inclination α) [13]. The pixel size in object space is 1.33 × 1.33 µm 2 . The setup is described in more detail in [7,8,11]. The intensity of every image pixel plotted against ρ yields a sinusoidal light intensity profile, which is normalized by the transmittance IT (transmitted light intensity averaged over all rotation angles ρ). The phase of the signal indicates the in-plane fiber direction angle ϕ, the amplitude | sin δ| (retardation) is related to the out-of-plane fiber inclination angle α, see coordinate system at the bottom. (c) Direction and retardation images of a coronal vervet monkey brain section (section 458, see Tab. I).

C. Scattered Light Imaging (SLI)
The setup used for the SLI measurements has been introduced by [11]. Figure 2 shows a schematic drawing of the setup. The light source consists of a green LED panel that emits mostly unpolarized and incoherent light with a wavelength of 525 nm. A diffuser plate is placed on top of the LED panel to ensure uniform illumination. The sample is placed on a specimen stage at distance h above the light source. A mask with a hole is placed on top of the diffuser plate so that the center of the sample is illuminated under a large polar angle θ (see Fig. 2(a)). The transmitted light intensity is recorded by a CCD camera (AVT, Oscar F-810C, Germany). The distance between sample and camera objective (l) can be varied to obtain a certain field of view. To ensure that only light that is scattered by the sample and falls mostly vertically onto the camera is recorded, the distance was chosen to be large enough (around l ≈ 30 cm). To avoid the detection of ambient light, the light path between the sample and the camera objective was sheltered by a dark tube (see Fig. 2(a)). Two different types of masks were used for most measurements (see Fig. 2(b)): one mask with circular holes with distance d = 11.5 cm to the center, diameters of 3.5 cm, and steps of 22.5 • (left)-and another mask with rectangular holes with distance d = 12 cm to the center, sizes of 2.4 × 4.0 cm 2 , and steps of 15 • (right). If not otherwise stated, the measurements were performed with the mask with rectangular holes and 15 • -steps. Some regions of interest ( Fig. 19 in Appx. D) were measured in 5 • -steps using rectangular holes with sizes of 0.4 × 4.0 cm 2 . Table I in Appx. C lists the settings of the SLI measurements (mask, sample/camera distances, exposure time, pixel size) for all measured samples. To uniquely identify each measurement (also for measurements of the same sample at different times), a unique identifier (e. g. #123) has been assigned to each measurement. The illumination angle θ is defined by the distance between the center of the mask and the outer holes d and the distance between mask and sample h via: θ = arctan(d/h), cf. Fig. 2(a). The measurements were performed for illumination angles between θ ≈ 46 • and 48 • , and for a pixel size in object space between 6.5 µm and 15.0 µm. The circular hole in the middle of the masks (cf. Fig. 2(b)) was used as a reference in order to align mask, probe, and camera, and to adjust the field of view of the camera in the bright field. The region of interest should be illuminated by the middle hole to guarantee a constant polar angle of illumination θ during the scattering measurement.
During the measurement, all holes but one were covered, and an image was taken for each of the outer holes, starting at the twelve o'clock position (φ = 0 • ) and rotating clock-wise. The resulting images form a series in which each image pixel contains a light intensity profile (SLI profile, I(φ)) with respect to the direction of illumination (azimuthal angle, φ), cf. Fig. 3. To reduce noise, four images were taken for each direction of illumination and averaged. In the current prototype setup, the position of the mask needs to be adjusted manually for each illumination angle, so that one SLI measurement (with ∆φ = 15 • ) takes about half an hour.

D. Image registration
In order to compare SLI and 3D-PLI measurements of the same tissue sample (Sec. III C), the transmittance images (average transmitted light intensity in a 3D-PLI measurement) were registered onto the corresponding SLI images (average transmitted light intensity in an SLI measurement) using in-house software tools based on the software packages [14], ELASTIX [15], and ANTs [16][17][18], which perform linear and non-linear transformations. The same transformation was applied to all other 3D-PLI modalities (direction and retardation images).
In order to compare SLI measurements performed at different times after tissue embedding (Sec. III B 2), the corresponding SLI images of the brain section (average transmitted light intensity I(φ) ) were registered onto each other as described above, and the same transformations were applied to the whole image stacks ( Figure 3 demonstrates how the line profiles from the SLI measurements (SLI profiles) were evaluated, for selected regions in a coronal section of a vervet brain. For each rotation angle of the rectangular mask (φ = {0 • , 15 • , . . . , 345 • }), an image was taken (see Fig. 3(a)). Figure 3(b) shows the resulting SLI profiles (average intensity of 10 × 10 pixels plotted against the direction of illumination φ) for four selected regions: (i) a region with mostly parallel in-plane fibers, (ii) a region with slightly inclined fibers, (iii) a region with steep out-of-plane fibers with large inclination angle α, and (iv) two in-plane crossing fiber bundles.

E. Evaluation of SLI profiles
As described by [11] (see also Fig. 5), parallel in-plane nerve fibers yield a line profile with two distinct peaks that lie approximately 180 • apart, cf. Fig. 3(b)(i). With increasing fiber inclination α, the two peaks move closer together until they merge in one broad peak, cf. Fig. 3(b)(iii). The middle position between the two peaks indicates the in-plane orientation (direction angle ϕ) of the nerve fibers. If two in-plane fiber bundles cross each other, the resulting line (i) Figure 3. Evaluation of SLI profiles. (a) Image series obtained from an SLI measurement of a coronal vervet brain section (upper left corner, see Tab. I#90). The green arrowheads at the image borders indicate the direction of illumination φ. (b) On the left: Schematic drawing of the nerve fiber geometries associated with the evaluated brain regions (in-plane, inclined, steep, and crossing nerve fibers). The angle ϕ denotes the orientation of the fibers in the xy-plane, the angle α the out-of-plane inclination angle. On the right: Corresponding SLI profiles (average intensity plotted against the direction of illumination φ) obtained from a region of 10 × 10 pixels, indicated by the red arrows in (a). The dashed vertical lines indicate the corrected positions of the determined peaks, the solid vertical lines the computed fiber direction angles ϕ.
profile shows four distinct peaks. Each pair of peaks (lying approximately 180 • apart) belongs to one fiber bundle, respectively, cf. Fig. 3(b)(iv). Hence, by investigating the peaks of the measured SLI profiles, it is possible to obtain information about the out-of-plane and-more importantly-the in-plane crossing angles of the nerve fibers in an investigated brain section.
To enable a better comparison between measurements, the intensity values for each illumination angle φ were divided by the average intensity over all illumination angles for each image pixel, yielding a normalized SLI profile The position of the peaks was computed with Python using the SciPy package (version 1.4.1) and the function scipy.signal.find peaks [19], which compares neighboring values to identify the local maxima in a one-dimensional array: To start with, each point in the SLI profile at which the slope changes from positive to negative was identified as peak. To include peaks at the outer boundaries of the SLI profile, the line profiles were extended taking the 360 • periodicity of the signal into account. To avoid that small irregularities in the SLI profiles are detected as peaks, only peaks with a certain prominence were considered as peaks. The prominence (scipy.signal.peak prominences) indicates how strongly a peak stands out from the background signal: Between the peak and the next higher points on the left and right side of the peak, the minimum values (min1 and min2) are determined. The prominence is the vertical distance between the top of the peak and the higher of the two minima (see Fig. 3(b)(ii) in red). If not otherwise stated, only peaks with a prominence ≥ 8 % of the total signal amplitude (I max − I min ) were considered, see Appx. A for derivation. Hence, the fifth peak in Fig. 3(b)(iv), for example, is assumed to be a noise artifact. To account for inaccuracies introduced by the 15 • -discretization of the SLI profiles, the determined peak positions were corrected by calculating the geometric center of the peak tips (with a height corresponding to 6 % of the total signal amplitude, see Appx. B). Hereby, we assumed that the scattering peaks generated by nerve fibers are symmetric and do not influence each other, so that the geometric center of the peak gives a better estimate of the mean fiber direction than the peaks in the discretized SLI profile. With this correction method, the peak positions can be determined with a Gaussian discretization error of ±2.4 • (see Fig. 15), which is much less than without correction (approx. uniform distribution between ±15 • /2).
The in-plane fiber orientations (direction angles ϕ, see solid vertical lines in blue/magenta/green in Fig. 3(b)) were computed from the corrected peak positions (dashed vertical lines in corresponding colors): In regions with one peak (iii), the fiber orientation was defined as the peak position itself (1). In regions with two peaks ((i) and (ii)), the fiber orientation was computed by the average position between the two peaks (1 + 2). In regions with four peaks (iv), the direction angle of the first bundle ϕ 1 was computed by the average position between the first and third peak (1 + 3), the direction angle of the second bundle ϕ 2 by the average position between the second and fourth peak (2 + 4), provided that the peaks lie 180 • ± 35 • apart. For six peaks with a pair-wise distance of 180 • ± 35 • , the direction angles for the three fiber bundles were computed equivalently.
To enable a fully automated evaluation of the SLI profiles for whole brain tissue samples, we developed the opensource software SLIX (Scattered Light Imaging ToolboX) which is available on GitHub (https://github.com/3d-pli/ SLIX). The software evaluates the peaks of all measured SLI profiles and generates different parameter mapscontaining information about the number of (prominent) peaks, the average peak prominence, width, and distance, as well as the derived in-plane nerve fiber orientations (cf. Fig. 8).

F. Comparison of SLI and scatterometry profiles
To validate the results from Scattered Light Imaging, the SLI profiles were compared to scatterometry profiles obtained from coherent Fourier scatterometry with non-focused, normally incident laser light [12]. In the SLI measurement, a whole brain section (field of view of several centimeters) can be studied at once with micrometer resolution, but only for a limited number of scattering angles, yielding a discretized SLI profile (with steps ≥ 15 • ). In the scatterometry measurement, on the other hand, it is only possible to image small spots of the brain section (with 0.1-1 mm diameter), but the measurement yields the full scattering pattern of the investigated region (only limited by the numerical aperture of the objective lens), see Fig. 4. The measured scattering patterns are therefore a good validation for the measured SLI profiles.
In the scatterometry measurement, a collimated, non-focused laser beam with 633 nm wavelength is transmitted through the sample (brain section) under normal incidence. The scattered light behind the sample is collected by an objective lens and measured by a CCD camera which records the Fourier transform of the image plane (scattering pattern), i. e. the intensity of scattered light on a hemisphere behind the sample projected onto the xy-plane (cf. Fig. 4(c)). For the comparison, we used measurement data from [12] (data set: [20]). The data set contains measured scattering patterns of 232 spots in different brain samples (circles in Fig. 16), measured with a laser beam diameter of 1.12 mm, a numerical aperture of 0.4, and an exposure time of 30 ms.
The illumination angle in the SLI measurement (θ ≈ 46 • -48 • ) is much larger than the maximum scattering angle recorded in the scatterometry measurement (θ ≤ arcsin(NA) ≈ 23 • ). To enable a comparison between SLI profiles and scatterometry measurements while taking the information of the whole scattering pattern into account, the intensity values were integrated from the center to the outer border of the pattern for each azimuthal angle φ. As the scattering pattern represents a distribution on a hemisphere, the integrals were computed as spherical (polar) integrals and plotted against φ. To not depend on the very details of single, individual nerve fibers, the resulting curves were smoothed out using a Savitzky-Golay filter with 45 sampling points and a second order polynomial [21], see [12] for more details. These smoothed polar integrals will be referred to as scatterometry profiles in the following. Figure 4(b) shows the SLI and scatterometry profile of the same tissue spot for two crossing optic tracts.
From the scatterometry profiles, we computed the position and prominence of the peaks (as described in Sec. II E) and compared the determined peak positions to the peak positions of the SLI profiles for the same tissue spots. To enable a direct comparison, it is crucial that both SLI and scatterometry profiles are generated for exactly the same spot. In the scatterometry measurement, marker dots were placed on the brain section to determine the initial position of the laser spot during the measurement, and micrometer screws were used to move the sample in the x/y-direction in steps of 0.5 mm or 1 mm. By imaging the marker dots on the sample with a digital microscope, the position of the measured spots in the brain sections could be retraced (see [12] for more details). To evaluate the same tissue regions with SLI, the image of the average transmitted light intensity I(φ) in the SLI measurement was aligned   Fig. 3(iv)). The SLI profile was computed by averaging the transmitted light intensity in the considered brain region for each angle φ. The scatterometry profile was obtained by computing the polar integral of the measured scattering pattern (φ = 0 • is defined along the y-axis, in clock-wise direction), and smoothing the resulting curve. (c) In the scatterometry measurement, a collimated laser beam with a diameter of 1.12 mm and a wavelength of 633 nm is transmitted through the sample under normal incidence, and the scattered light on a hemisphere behind the sample is measured (lower image). The resulting scattering pattern shows the scattered light intensity projected onto the xy-plane (upper image). The measurement was performed with a numerical aperture NA = 0.4 so that the measured scattering pattern contains only scattering angles ≤ arcsin(0.4) = 23.6 • (see red circle), which is much less than the illumination angle θ ≈ 48 • used in the SLI measurement (in orange).
with the corresponding image of the digital microscope, and the determined positions of the spots were transferred to the SLI image. First, the alignment was done manually using anatomical structures as reference, and the SLI profiles were computed by averaging the transmitted light intensity in each 1.12 mm tissue spot and for each angle φ in the SLI measurement. To improve the alignment of the images, the determined peak positions of the SLI profiles were compared to the peak positions of the corresponding scatterometry profiles in order to determine systematic errors. To ensure that only correct peak pairs are compared to each other, the peaks in the scatterometry profiles (with prominence ≥ 8 %, cf. Sec. II E) were only compared to the peaks in the SLI profiles if both line profiles have the same number of peaks. When the difference between the peak positions showed a significant systematic shift for all line profiles in one brain section, i. e. the mean is not compatible with zero (|mean| > 3σ/ √ N, with standard deviation σ and number of peak differences N ), the SLI images of the brain section were rotated by this value. The positions of the scatterometry spots in the SLI images were corrected accordingly. Figure 16 in Appx. C shows all evaluated spots in the brain sections after the correction, and Tab. I lists the rotation angles ρ by which the SLI images were rotated with respect to the scatterometry measurements. To enable a better comparison between SLI and scatterometry profiles, the minimum intensity values were subtracted from all profiles and the resulting intensity values were divided by the average signal intensity, respectively, yielding normalized line profiles I N (φ).

G. Comparison to simulated scattering patterns
In some brain regions (Sec. III C 2), we compared the SLI profiles to simulated scattering patterns [11]. These scattering patterns were computed by modeling the propagation of light through artificial nerve fiber constellations with a massively parallel 3D Maxwell Solver (software TMDE3D TM [22,23]) based on a conditionally stable finitedifference time-domain algorithm [24,25]: Space and time are discretized and the electromagnetic field components of the light wave are computed by approximating Maxwell's equations by second-order central differences (see [26] for more details). The simulations were performed for a plane wave with coherent, circularly polarized light (λ = 550 nm) and normal incidence. The nerve fibers were modeled by an inner axon and a surrounding myelin sheath with different refractive indices and an average fiber diameter of 1 µm, in a volume of 30 × 30 × 30 µm 3 (see Fig. 5, top left). The scattering patterns (see Fig. 5, top right) were computed from the intensity distribution of the scattered light wave behind the sample on a hemisphere projected onto the xy-plane (cf. Fig. 4(c)).
In contrast to the scattering patterns obtained from the scatterometry measurements, the simulated scattering patterns contain all possible scattering angles for the transmitted light (θ = [0 • , 90 • ]) and are not limited by the numerical aperture. Therefore, the simulated scattering patterns can be directly sampled along the circle corresponding to the polar angle of illumination (θ ≈ 48 • ) used in the SLI measurement.
To take the finite size of the illumination beam angle into account, the simulated scattering patterns were blurred with a Gaussian filter with a diameter corresponding to the diameter of the circular hole in the mask (8 • in angular space). To consider all possible discretizations of the SLI profile, the line profiles were then evaluated without discretization along a circle in the (blurred) scattering pattern corresponding to θ = 48 • , starting at the twelve o'clock position (along the y-axis) and running clock-wise. Figure 5 shows the simulated scattering patterns for models of parallel in-plane fibers, parallel out-of-plane fibers, and in-plane crossing fibers, as well as the resulting line profiles (taken from [11], fig. 9(b)). Parallel in-plane fibers Parallel out-of-plane fibers Crossing in-plane fibers

A. Comparison to scatterometry measurements
To validate our results obtained from Scattered Light Imaging, we compared the (normalized) SLI profiles of some brain tissue spots (with 1.12 mm diameter) in two/three crossing sections of human optic tracts and coronal/sagittal vervet brain sections to the corresponding scatterometry profiles obtained from coherent Fourier scatterometry (see Sec. II F). From the 232 tissue spots marked in Fig. 16, we used 200 spots for the comparison (non-white circles): We excluded tissue spots that do not lie completely inside the tissue or that belong to anatomically known regions with out-of-plane nerve fibers (like fornix or cingulum in the coronal vervet brain sections). The reason to exclude these tissue spots is that the scattering patterns are not expected to be radially symmetric (cf. middle image in Fig. 5 and [11] fig. 7(a)). As the measured scattering patterns contain only small scattering angles (θ ≤ 23 • ), the scatterometry profiles (polar integrals) of out-of-plane fibers are expected to be very different from the corresponding SLI profiles which are obtained from a much larger scattering angle (θ ≈ 48 • ). For this reason, the SLI profiles for out-of-plane fibers were not compared to the corresponding scatterometry profiles, but only to simulated line profiles (see later, Sec. III C 2).
For all 200 tissue spots, the determined peak positions of the SLI profiles were compared to the peak positions of the corresponding scatterometry profiles. As the peaks of the scatterometry profiles were shown to be reliable [12], also smaller peaks (with prominence ≥ 3 %) were considered for the comparison. For the SLI profiles, only peaks with prominence ≥ 8 % were considered. To increase statistics, all SLI and scatterometry profiles were compared to each other, even if the number of peaks differs. To ensure that only peaks belonging to the same structures were compared to each other, we matched them by maximizing the sum of inverse squared differences between SLI and scatterometry profiles. That is, we optimized the matching of peaks by putting more weight on peak pairs with a small difference in peak position. Left-over peaks (without partner) were not considered (cf. dashed vertical lines in Fig. 6(b)). In total, 523 peaks could be compared to each other. Apart from the difference between peak number and position, we also determined the sum of absolute distances between the normalized SLI and scatterometry profiles φ |I N (φ) SLI − I N (φ) Scatt | (in steps of ∆φ = 1 • , SLI linearly interpolated) to study the overlap between the two curves. When comparing the SLI and scatterometry measurements, it should be taken into account that the two types of measurement show differences: Up to several weeks elapsed between the measurements (cf. Tab. I). Furthermore, the measurement setups differ substantially (SLI: single large scattering angle with high spatial resolution and unpolarized incoherent light (monochromatic wavelength 525 ± 25 nm); scatterometry: full scattering pattern with smaller scattering angles, low spatial resolution, and polarized coherent laser light (monofrequent wavelength 633 nm)). Finally, the localization of the tissue spots measured with scatterometry is only accurate to a certain extent. Nevertheless, the results show that the SLI and scatterometry profiles are compatible with each other, hence validating the SLI approach.
When considering the sum of distances between the curves ( Fig. 6(a)-(c)), it is noticeable that tissue spots with in-plane parallel nerve fibers (cc, in blue) show a better correspondence of the line profiles than tissue spots with crossing fibers. Especially in very inhomogeneous regions like in the corona radiata of the vervet brain, we observe larger differences (cr, in magenta), see also Fig. 16.

B. Optic tracts as model system for crossing nerve fiber bundles
To test the SLI measurements and demonstrate that they can be used to extract additional information like the individual fiber orientations in regions with crossing nerve fibers, we considered model systems of two and three crossing fiber bundles with well-defined crossing angles. For this purpose, up to three sections of optic tracts obtained from a human optic chiasm were manually placed on top of each other under different crossing angles (see Fig. 7(a),(c)). The same model systems have also been used in [11,12]. First, we show in different studies that the fiber orientations and other tissue parameters can be reliably and automatically reconstructed, not only for selected brain regions but for all image pixels of the sample (Sec. III B 1). Then, we study how the time after tissue embedding influences the scattering of the brain tissue and the contrast obtained in SLI measurements (Sec. III B 2).

Reconstruction of nerve fiber orientations
2x Human Optic Tracts s0036 (#162) Vervet493_corona2_x-4mm_y+0mm (orange) Vervet493_corona1_x+1mm_y-1mm (magenta) between the solid and dashed curves demonstrates the conformity of scattering in pixels from such areas allowing for the selection of small representative regions. In regions with parallel nerve fibers (i), the SLI profiles show two distinct peaks that lie 180 • apart; the arithmetic mean value of the peak positions corresponds to the orientations of the fibers in the section of the optic tract (indicated by green/magenta/yellow lines in (c)). The images in (b) show that one tract lights up when being illuminated from the broadside, thus revealing its underlying structure. In regions with crossing fibers (iii), the SLI profiles show several distinct peaks, where each 180 • -peak pair corresponds to one fiber tract direction. The crossing angles of all fiber tracts can be correctly determined in the crossing region, not only for two crossing sections of optic tracts with 90 • and 70 • crossing angles but also for three crossing sections with 55 • crossing angle. A comparison between the graphs in (i) and (ii) shows that the signal in the crossing region is a superposition of the signals in the corresponding single tissue layer: The peaks in (i) and (ii) occur at similar positions, see vertical dashed lines. Figure 19 in Appx. D shows the SLI profiles for the three crossing sections of optic tracts obtained from an SLI measurement with 15 • and 5 • -steps. Although the SLI profiles show more details for 5 • -steps (red curves), the determined peak positions correspond very well to the peak positions obtained from the SLI profiles with 15 • -steps (black curves). To reduce measurement time, all following SLI measurements were therefore performed with 15 • -steps.    Figure 8 shows the different parameter maps that were generated from an SLI measurement of two crossing sections of optic tracts, using the SLIX software (see Sec. II E): For each image pixel, the corresponding SLI profile was normalized, the peaks were determined, and the peak prominences and peak positions were computed. In addition, the peak width was determined as the full width of the peak at a height corresponding to the peak height minus half of the peak prominence. For regions with two peaks with prominence ≥ 8 %, the distance between the peaks was computed. For regions with two and four prominent peaks with pair-wise distances of 180 • ±35 • , the fiber orientations (ϕ 1 and ϕ 2 ) were computed from the arithmetic mean values of the peak pair positions, as described in Sec. II E.
The section of the optic chiasm (section no. 36, see also Fig. 17 in Appx. D) consists not only of the optic tracts, but also of two other anatomical regions, labeled in Fig. 8(a): While the optic tracts (ot) contain mostly parallel, densely packed myelinated nerve fibers of the visual pathway, the supraoptic decussation (sod) is located dorso-medial to the optic tract and contains less myelinated nerve fibers. The supraoptic nucleus (SON) is part of the hypothalamus located on the superior and the medial side of the optic tract and is occupied by a large portion of cell bodies. Figure 18 in Appx. D shows polar plots of the SLI profiles for the two crossing sections of optic tracts (sampled at each 24 th image pixel). The morphology of the resulting shapes clearly shows the transition between parallel and crossing nerve fibers.
The average intensity of the SLI profiles in Fig. 8(a) is a measure for the scattering of the sample. While the optic tracts are highly scattering, the supraoptic decussation, which contains less myelinated nerve fibers, is less scattering. In regions belonging to the supraoptic nucleus, which contains mostly cell bodies, the SLI profiles show a large number of randomly distributed peaks with a small average peak prominence (see Fig. 8(b),(c) on the left) and cannot be reliably evaluated. When taking all peaks into account (left image in Fig. 8(b)), many regions are not correctly assigned (regions with parallel fibers show three peaks or more; regions with crossing fibers show five peaks ore more). When considering only peaks with prominences ≥ 8 % (right image in (b)), regions with parallel and crossing fibers are more correctly assigned (parallel fibers with two peaks are shown in magenta, two crossing fiber bundles with four peaks are shown in green), as also found in Appx. A. The average prominence of the peaks is larger in regions with parallel fibers and lower in regions with crossing fibers, the same holds also for the average peak width (c). The distance between the peaks, evaluated for regions with two prominent peaks (d), lies mostly between 170 • and 180 • , showing that most fibers are orientated in-plane, as expected. At the image borders, some values lie between 160 • and 170 • , suggesting that the fibers are more inclined there. The direction maps in (e) show the two fiber orientations (ϕ 1 and ϕ 2 ) computed from the different peak pairs in regions with two and four prominent peaks. The individual nerve fiber orientations in the two sections of the optic tracts are clearly visible in the crossing region, and transition smoothly to the regions with parallel fibers (single tissue layer). The image on the right shows the in-plane fiber orientations as blue and green lines computed from the direction angles (for every 24 th image pixel). The lines nicely show the course of the fibers in the optic tracts, not only for regions with parallel fibers (single tissue layer), but also within the crossing region (double tissue layer). Figure 9(a-b) shows the reconstructed (in-plane) nerve fiber orientations from 3D-PLI and SLI measurements for the two crossing sections of optic tracts. 2 The high-resolution parameter maps can be found in the data repository [27]. The images on the left (i) show the fiber directions for each image pixel (in different colors according to the color legend in (ii)). For image pixels in which the SLI measurement yields different fiber direction angles (ϕ 1,2,3 , cf. Fig. 8(e)), only a single fiber direction angle is displayed. To make individual nerve fiber orientations in the whole sample visible, the images in the middle (ii) show the fiber direction (a) for every 150 th image pixel, and (b) for every 30 th image pixel. The SLI fiber directions were replaced by the median 3 fiber direction for every 3 × 3 pixels beforehand. The right images (iii) show a zoom-in of (i) for a small region marked by the arrow, respectively.
In contrast to SLI, the 3D-PLI measurement yields only a single nerve fiber orientation for each measured image pixel. In regions with non-crossing nerve fibers (single tissue layer), 3D-PLI and SLI yield similar (in-plane) nerve fiber orientations. In regions with crossing nerve fiber layers, the fiber directions from 3D-PLI correspond to the fiber directions in one or the other fiber layer (magenta/dark blue or yellow/green), and in some regions to an intermediate fiber direction of both (cyan), see Fig. 9(a)(ii) and (iii). This observation is consistent with simulations by [28] who showed that in regions of two crossing fiber bundles with nearly 90 • crossing angle, the fiber direction obtained from 3D-PLI corresponds to the direction of the fiber bundle that contains > 50 % of the fibers. For a slightly different crossing angle, the resulting fiber direction is an intermediate direction of the two fiber bundles.
While the 3D-PLI measurement (a) yields only a single nerve fiber orientation, also in regions with crossing fibers, the SLI measurement (b) clearly shows the individual nerve fiber orientations of the two crossing optic tracts. The SLI measurement reliably reconstructs the individual nerve fiber orientations also in regions with three crossing fiber bundles/layers (c). The course of the three individual fiber bundles is clearly visible both in regions with two and three crossing sections of optic tracts (see zoomed-in region on the right).

Time-dependent changes in scattering
To increase the birefringence and scattering contrast, all brain tissue samples were embedded in a glycerin solution before being cover-slipped (see Sec. II A). To study the effect of embedding time on the scattering properties of the brain tissue, we measured the same brain tissue samples at different times after tissue embedding (from 3 days up to 8 weeks after embedding). Figure 10(b) shows the (non-)normalized SLI profiles for two crossing sections of optic tracts, evaluated for a region with parallel fibers (i) and two regions with crossing fibers ((ii),(iii)) obtained from SLI measurements at different times after tissue embedding.
With increasing time elapsed after tissue embedding, the scattering decreases notably, especially in the first few weeks (see (iii)). Within 8 weeks, the average light intensity of the SLI signal decreases by more than 70 %, after 5 months it is reduced by more than 85 %. Figure 10(d) shows the average intensity evaluated along a geometric line profile (yellow line in (a)). The vertical dashed lines indicate the borders of the crossing region. While the crossing region (double tissue layers) shows stronger scattering than the region with parallel fibers (single tissue layer) right after the tissue embedding, it is the other way around after 2 weeks. Especially interesting is the comparison of the normalized SLI profiles: While the profiles in regions with parallel fibers remain almost the same (see lower image in (i)), the peaks in regions with crossings fibers become more prominent with increasing time after tissue embedding (see lower image in (ii)), i. e. the amplitude of the intensity profile does not scale down proportionally to the lowering of its average, but less. Figure 10(c) shows the brain areas for which four peaks are detected (interpreted as two crossing fiber bundles) in different colors, depending on the time elapsed after tissue embedding. Three days after tissue embedding, only a small part within the crossing region is correctly classified as two crossing fiber bundles (in red). This region grows as time goes by (yellow and blue). After 8 weeks, four peaks are detected in the majority of the pixels in the crossing region (cyan), but also in the gray matter (SON). As the prominence of the peaks in the SLI profiles is a crucial criterion for a correct peak detection, a longer lapse of time between embedding and measurement leads to a more reliable fiber crossing classification even at a lower level of scattering (see Fig. 10(b)(ii)). One possible explanation why the fiber orientations in regions with crossing fibers are more correctly determined for brain sections with longer embedding time and reduced scattering might be that strongly scattering fiber bundles influence each other (higher-order scattering), reducing the prominence/contrast of the resulting peaks in the SLI profile.

C. Nerve fiber architectures in whole brain sections
So far, we have considered model systems of two and three crossing nerve fiber bundles. In the following, we show that SLI measurements can be applied to whole brain tissue samples and reliably reconstruct the (crossing) nerve fiber orientations in the brain. First, we compare SLI and 3D-PLI measurements of the same brain section and show that SLI yields similar in-plane fiber orientations as 3D-PLI in regions with non-crossing nerve fibers (Sec. III C 1). Then, we study the SLI profiles for regions with inclined nerve fibers and show that they also agree with predictions of simulated scattering patterns (Sec. III C 2). Finally, we show that crossing nerve fibers in complex brain regions, like the corona radiata of the vervet brain, can be reliably reconstructed (Sec. III C 3).

Comparison between SLI and 3D-PLI
To enable a direct comparison between the nerve fiber orientations obtained from SLI and 3D-PLI, the 3D-PLI measurement of the same brain section was performed directly after the SLI measurement, and the images (average transmitted light intensities) were registered onto each other (as described in Sec. II D). Figure 11 shows the comparison between SLI and 3D-PLI for a coronal rat brain section. For the comparison, we considered only white matter regions  with non-crossing nerve fibers, i. e. regions with one or two prominent peaks in the SLI profiles (prominence ≥ 8 %, see Sec. II E). Subfigures (c) and (e) show the differences between the in-plane fiber direction angles ϕ obtained from SLI and 3D-PLI. Note that ϕ = 0 • is here defined along the positive x-axis, running counter-clockwise (while the direction of illumination φ = 0 • is defined along the positive y-axis, running clockwise). The 2D histogram in (e) shows that the differences are approximately normally distributed with a standard deviation of about 23 • (for about 766 300 evaluated image pixels). The distribution of differences strongly depends on the out-of-plane angle α of the nerve fibers and, thus, on the retardation | sin δ| (δ ∝ cos 2 α): In regions with steep out-of-plane nerve fibers (i. e. regions with large α and low retardation, respectively), SLI and 3D-PLI yield very different direction angles (red and blue dots in (c)), e. g. in the cingulum (cg), dorsal fornix (df ), or stria medullaris (sm). The standard deviation is more than 20 • for regions with a retardation below 0.3 and up to 40 • for regions with a retardation below 0.1 (see lower graph in (e)). In regions with in-plane nerve fibers (i. e. regions with small α and large retardation values), SLI and 3D-PLI yield very similar direction angles (gray in (c)), e. g. in the corpus callosum (cc), fimbria (fi), or optic tract (ot). The standard deviation is about 5 • or less for regions with a retardation of 0.7 or higher. Taking into account that the error induced by the 15 • -discretization of the SLI profiles is already ±2.4 • (see Appx. B), the fiber orientations reconstructed from the SLI measurement correspond very well to the in-plane fiber orientations determined by 3D-PLI and can be considered as reliable.
Subfigure (f) shows the direction angles of the SLI and 3D-PLI measurements in the stria medullaris and the graph at the bottom their distribution in the delineated region. In the stria medullaris, a fiber tract is running steeply with respect to the coronal section plane. In such regions, 3D-PLI is prone to random in-plane direction angle distributions due to the low retardation signal (b). While the 3D-PLI measurement (red curve) yields randomly distributed direction angles, the SLI measurement (blue curve) yields one dominant fiber direction around ϕ SLI ≈ 70 • . This result is in accordance with the course of the fiber bundles in the stria medullaris which are expected to have a vertical but no lateral component in the coronal section plane, slightly declining on their rostral course. Hence, in some regions with out-of-plane nerve fibers, the SLI measurement patches the lack of reliability in 3D-PLI measurements and can be used as corrective. However, the SLI peaks merge for out-of-plane fibers and become smaller (see Fig. 12). This behavior increases the risk of peaks not being detected and fiber directions determined from SLI profiles being shifted to wrong values. Figure 11(d) and (g) show the distance between two peaks in the SLI profiles (evaluated in regions with one or two prominent peaks) in dependence on the retardation. A peak distance of 0 • means that only a single peak in the SLI profiles was detected. While regions with high retardation values (in-plane fibers) show always large peak distances between 160 • and 180 • as expected, regions with low retardation values (out-of-plane fibers or regions with less myelinated fibers) exhibit a broad distribution of peak distances from about 30 • to 180 • . The standard deviation ranges from 7.7 • (at retardations above 0.8) to 59.4 • (at retardations below 0.1).
Subfigure (h) shows a zoomed-in region of the retardation image (optic tract), where regions with large peak distances (160 • -180 • ) and small retardation values (≤ 0.5) are highlighted in red. Although the optic tract (in white) has much higher retardation values than the neighboring fiber bundle (supraoptic decussation, in red), the SLI profiles (here shown for a region with retardation 0.4) show similar peak distances (170 • -180 • ) for both fiber bundles, suggesting that the two fiber bundles contain only moderately inclining nerve fibers. In fact, the supraoptic decussation follows mostly the course of the optic tract in this region [29], but it contains less myelinated nerve fibers than the optic tract, which explains why we obtain lower retardation values. Without any correction, this would lead to an over-estimation of the out-of-plane nerve fiber inclinations in 3D-PLI. Therefore, the inclination angles in 3D-PLI are corrected by taking the transmittance 4 into account. When using a tiltable specimen stage in the 3D-PLI measurement, the fiber inclination can be determined even more reliably, however at lower resolution (see later Sec. IV D.) The fact that the retardation is not always directly related to the out-of-plane inclination of the nerve fibers also explains why not all values in the 2D histogram in (g) follow the black curve in the histogram, which shows the simulated peak distances for fiber bundles with different inclination angles (see later, Fig. 12(d)), assuming that the retardation | sin δ| is directly related to the inclination (δ ∼ cos 2 α) without applying any correction. Taking into account that the fiber inclination might be over-estimated by less myelinated fibers (low retardation values), the prediction of the simulation is setting a lower limit, which is also observed in the 2D histogram. Hence, SLI yields not only a good reference for the fiber directions in regions with out-of-plane nerve fibers, but also for the fiber inclinations in regions with lower myelination. Figure 12 shows the SLI profiles for nerve fibers with different out-of-plane inclination angles at the transition zone between in-plane and out-of-plane fibers in the corpus callosum (cc) and cingulum (cg) of a coronal and sagittal vervet brain section. In the coronal (sagittal) section, the fibers in the corpus callosum are mostly oriented in-plane (out-of-plane). For fibers in the cingulum, it is the other way around. To study the transition between in-plane and out-of-plane fibers, the SLI profiles were evaluated for a chain of neighboring regions (averaged over 3 × 3 pixels) along a straight line (see rainbow-colored lines in the zoomed-in areas in (a)). The blue/violet curves belong to regions with in-plane fibers, the yellow/red curves to regions with out-of-plane fibers, the green curves to the transition zone. In-plane fibers (coronal: cc, sagittal: cg) show two distinct peaks lying 180 • apart; out-of-plane fibers (coronal: cg, sagittal: cc) show one broad peak. With increasing inclination, the distance between the peaks and the peak prominences decrease, while the peak widths increase. Figure   prominence in comparison to the retardation for the coronal rat brain section introduced in Fig. 11. As expected, the average peak width increases and the average peak prominence decreases with increasing inclination (decreasing retardation).

Inclined nerve fibers
The prediction of the simulation studies by [11] is shown in Fig. 12(c) and (d): For different inclination angles (α = {0 • , 10 • , . . . , 90 • }), the authors generated simulated scattering patterns ( fig. 7(a) and S3(c)), which were evaluated as described in Sec. II G (average intensity values evaluated along the white dashed circle in (c)). Subfigure (d) shows the simulated line profiles in one plot. As observed for the SLI profiles (b), the peak distances in the simulated line profiles decrease, the peak widths increase, and the peak prominences decrease with increasing fiber inclination, until the two peaks merge into one broad peak for α ≥ 80 • . In Fig. 11(g), the peak distances of the simulated line profiles are plotted against the measured retardation of the rat brain section (| sin δ| = | sin(arcsin(ret max ) cos 2 α)|, where ret max is the maximum retardation value of the brain section, connecting the tissue retardation to the fiber inclination. The inclination dependence of the measured SLI profiles corresponds very well to the description of the scattering profiles as a function of fiber inclination in the simulation studies, validating again the benefit of the light scattering approach for the elucidation of the nerve fiber architecture in the brain. Figure 13 shows the nerve fiber orientations in two coronal vervet brain sections reconstructed from an SLI measurement (c),(e). The fiber orientations are encoded in different colors, according to the color bubble in (f). The colored images show the fiber orientations for each image pixel; the vector maps in (c) and (e) show the fiber orientation for every 40 th and 50 th pixel in the SLI measurements. The SLI fiber directions were replaced by the median fiber direction for every 3 × 3 pixels beforehand. The yellow rectangles mark the region of the corona radiata shown in the zoomed-in areas. Subfigure (f) shows a sketch of the crossing nerve fiber pathways in the corona radiata known from vervet brain atlases [30]. Figure 13(b) shows the nerve fiber orientations reconstructed from a 3D-PLI measurement of the same brain section shown in (c).

Crossing nerve fibers
The pathways of the crossing fibers in the corona radiata (which are already visible in the structure of the retardation image) become clearly visible in the SLI images (white ellipse). For better reference, the retardation image is shown in the background of the vector maps. As expected, regions with small retardation values correspond to regions with crossing nerve fibers, where SLI yields additional information. The fiber pathways reconstructed from SLI correspond to the course of the fiber pathways shown in (f): In addition to the course of the corpus callosum (5)-(4), the vector map in (c) also shows the fiber bundle starting at (5) and fanning out to (2) and (3) in red/magenta, as well as the fiber pathways between (2)-(3) and (1)-(4) in green/yellow. The same pattern of crossing nerve fibers is also visible in the SLI vector map of the neighboring vervet brain section (e).

IV. DISCUSSION
For building a detailed network model of the brain, a correct reconstruction of crossing nerve fiber pathways is crucial. Here, we introduced Scattered Light Imaging (SLI), a neuroimaging technique that allows to reconstruct the individual orientations of (crossing) nerve fiber pathways with micrometer resolution in whole brain tissue samples.

A. Evaluation of SLI profiles
Previous work [11,12]) has shown that the distribution of light scattered on nerve fiber configurations (scattering pattern) reveals the individual fiber orientations in regions with crossing nerve fibers. In [12], scatterometry is used to measure the scattering patterns for different brain regions: a brain section is illuminated by collimated, normally incident laser light with a minimum beam diameter of 100 µm (limited by diffraction), yielding a comparable resolution to the one achieved by postmortem dMRI [2]). In SLI, we use an inverse setup, where the brain section is illuminated under different angles and the scattered light behind the sample is measured under normal incidence. In this way, we can overcome the limitations of the scatterometry measurement and measure whole brain tissue samples (with a field-of-view of several centimeters) and with a spatial resolution that is only limited by the pixel size in object space (currently: ≥ 6.5 µm).
However, in contrast to the scatterometry measurements, SLI takes only some scattering angles into account, yielding discretized SLI profiles (here: 24 data points in 15 • -steps). When analyzing the SLI profiles, the discretization (limited number of data points) leads to two major sources of errors: First, the peak positions (and thus the derived nerve fiber directions) can be determined only with limited accuracy. With the correction procedure presented in Appx. B, we could notably increase the accuracy of the detected peak positions and reduce the discretization error to ±2.4 • , keeping the number of required measurements low. Second, the signal is more susceptible to noise. While the fiber directions in 3D-PLI can be determined from a sinusoidal fit (requiring only a small number of data points, cf. Fig. 1(b)), there exists no simple fit function for the SLI profiles. The interpretation of the measured signals depends on the number of detected peaks (cf. Fig. 3(b)) and can be significantly changed by single outliers. When taking all detected peaks into account (left image in Fig. 8(b)), many regions (especially regions with crossing nerve fibers) are not correctly assigned. When considering only peaks with a prominence ≥ 8 % of the total signal amplitude (right image in Fig. 8(b)), we could notably increase the number of correctly assigned regions (96.5 % for in-plane parallel nerve fibers, 65.2 % for two crossing fiber bundles, see Appx. A). In regions with three crossing fiber bundles, only 36.3 % of the regions are correctly assigned, but in 96.7 % of the cases there are at least two different fiber orientations detected. Note that when some brain regions are not correctly assigned, this means that some existing fiber orientations are not determined, but the assigned fiber orientations can safely be assumed as correct. Therefore, missing fiber orientations can be estimated from the fiber orientations in neighboring pixels. As a consequence, the effective resolution is reduced, but the achievable resolution (e. g. 13 µm for 2 × 2 pixels) is still much higher than for dMRI.

B. Validation of results
To demonstrate that SLI yields reliable results, we compared the fiber directions obtained from SLI to fiber directions obtained from scatterometry (Sec. III A) and 3D-PLI (Sec. III C 1), and showed that the differences in fiber directions are normally distributed with a standard deviation of ±8.8 • (scatterometry) and ±5 • (3D-PLI, for regions with mostly non-crossing in-plane nerve fibers with retardations above 0.7). Taking into account the different measurement setups and a discretization error of ±2.4 • , this is a good correspondence. For regions with out-of-plane nerve fibers, we could show that the SLI profiles are compatible with simulated scattering patterns (Sec. III C 2). As the distance between the peaks decreases with increasing nerve fiber inclination until the peaks finally merge (see Fig. 12), the interpretation of the SLI profiles is very susceptible to noise so that an exact computation of the out-of-plane fiber inclination remains challenging.
Another validation of our measurement technique are the model systems of crossing sections of optic tracts with well-defined crossing angles (Sec. III B). We could show that SLI yields reliable (in-plane) nerve fiber orientations for each image pixel, not only in regions with parallel nerve fibers, but also in regions with (up to three) crossing nerve fiber bundles/layers, yielding valuable additional information to the (single) nerve fiber orientations obtained from 3D-PLI (cf. Fig. 9(a)).

C. Reconstruction of crossing nerve fibers
In Sec. III C, we demonstrate that SLI measurements can be applied to whole brain sections and enable the reconstruction of individual nerve fiber orientations in brain regions with highly complex, crossing nerve fibers, like the corona radiata of the vervet brain (Fig. 13). The reconstructed fiber directions correspond to anatomically known fiber pathways and to visible structures in the high-resolution 3D-PLI images.
In principle, the structural information of 3D-PLI images can be used as an estimate for the underlying nerve fiber organization in the measured brain sections. In many brain regions, the course of crossing nerve fibers can be visually tracked by means of visible structures in the high-resolution 3D-PLI images. However, as 3D-PLI yields only a single fiber orientation vector for each measured tissue voxel, individual nerve fiber orientations in regions with crossing fibers cannot be automatically extracted. By taking directional information of neighboring image pixels into account, e. g. by creating orientation distribution functions from high-resolution 3D-PLI vector data [31], the orientations of crossing nerve fibers can be estimated at the expense of resolution. A "super-pixel" formed from 5 × 5 pixels in 3D-PLI (px = 1.33 µm) has a similar size as one image pixel in SLI (px = 6.5 µm). In regions with crossing nerve fibers, 3D-PLI often yields erroneous (intermediate) nerve fiber orientations, which correspond to neither fiber orientation in the crossing region (cf. Fig. 9(a)(ii) and (iii) in cyan). In this case, the super-pixel also contains erroneous fiber orientations. Moreover, the minimum size of the super-pixels is fundamentally limited (they must contain several pixels and the pixel size in 3D-PLI cannot be much below 1.33 µm [13]). With SLI, the individual fiber orientations in regions with (up to three) crossing nerve fibers can be reliably reconstructed within one measured image pixel, and the resolution can be easily increased to a comparable resolution as 3D-PLI by using a different objective lens.
In a long-term study (Sec. III B 2), we could show that even for weakly scattering tissue, the fiber orientations are still reliably determined by SLI, especially in regions with crossing nerve fibers: The scattering decreases with increasing time after tissue embedding, thus increasing the signal-to-noise ratio in regions with crossing nerve fibers.
A possible explanation is that the glycerin solution in which the brain sections are embedded (see Sec. II A) slowly diffuses over time, which reduces the refractive index contrast of the tissue and, thus, the scattering. When "old" brain sections are freshly embedded in glycerin solution ("revitalization"), the amount of scattering is comparable to what was observed right after the initial tissue embedding.

D. Combination of SLI and 3D-PLI
The SLI measurements do not only provide additional information in regions with crossing nerve fibers, but also yield a good reference and corrective for fiber directions in other brain regions with small 3D-PLI signals. In some regions with steep out-of-plane nerve fibers, the in-plane fiber orientations derived from SLI were shown to be more reliable than those obtained from 3D-PLI (see Fig. 11(f)). In some regions with less myelinated nerve fibers (sod in Fig. 11(h)), we could show that the SLI measurement yields also a good estimate for the out-of-plane angle, which can be used as validation for the fiber inclinations determined by 3D-PLI, which need to be corrected e. g. by transmittance weighting to compensate for low myelination.
While previous studies of our group [11] have shown that the transmittance images obtained from a 3D-PLI measurement can be used to classify regions with low birefringence signals into regions with in-plane crossing, outof-plane, and less myelinated fibers, our new technique reliably reconstructs the in-plane fiber orientations in these regions and gives an estimate of the out-of-plane fiber inclinations. When using a tiltable specimen stage [32,33], the out-of-plane nerve fiber inclinations can be reliably determined with 3D-PLI, also in regions with low fiber densities. However, tilting is not (yet) available for microscopic resolutions ¡ 10 µm.
Because of the similar setup (cf. Fig. 1(a) and Fig. 2(a)), the SLI measurement can be easily combined with a 3D-PLI measurement, which significantly improves the high-resolution 3D-reconstruction of nerve fibers: The directions of in-plane crossings fibers (determined from SLI) can be used, for example, as a-priori information for fiber tractography algorithms. As SLI and 3D-PLI have different sources of error, SLI can serve as cross-validation for 3D-PLI, also when using a tiltable specimen stage. 3D-PLI yields always a (single) fiber orientation for each image pixel. SLI can be used to indicate how reliable the determined orientations are and identify problematic regions. The orientations of parallel in-and out-of-plane fibers can be determined accurately with 3D-PLI, especially in regions with a high density of myelinated fibers. In regions with less myelinated fibers or crossing regions, the fiber directions from SLI can be used as additional information. Also in regions with out-of-plane nerve fibers, SLI can be used as cross-validation to determine the correct in-plane direction angles of the nerve fibers and get an estimate of the out-of-plane angles.

E. Perspectives for further improvements
The SLI measurements yield valuable additional information, improving the state-of-the-art of nerve fiber reconstruction in regions with crossing nerve fibers, steep nerve fibers, or regions with a low density of myelinated fibers. However, there is potential for further development. One major source of error is the correct detection of peaks in the SLI profiles. Future studies should focus on how to combine different parameter maps (cf. Fig. 8) to improve the reconstruction of nerve fiber orientations also in regions with less prominent peaks, e. g. by including the peak prominence as a measure for the signal-to-noise ratio.
For larger brain sections like the vervet brain, only the inner part of the measured brain region can be used for evaluation: As the image borders are not illuminated under the same polar angles while rotating the hole of illumination, these regions get more light for some rotation angles than for others, causing a bias in the resulting SLI profiles and an erroneous reconstruction of nerve fiber orientations (cf. black arrow in Fig. 13(c)). In the future, it should be studied how neighboring images of the same brain region can be combined to avoid this effect. When using a different type of illumination that generates an approximately plane wave (e. g. using apertures or lenses), it would be possible to illuminate the whole sample under the same polar angle so that these problems do not occur.
To enable more precise measurements, the SLI measurement should be fully automated so that a larger number of scattering angles can be measured. It should also be studied if the measurements of other scattering angles (using different types of masks) yield additional information about the substructure of the tissue.
The pixel size in the SLI measurements (here ≥ 6.5 µm) can still be improved by using a different objective lens. To avoid artifacts caused by image registration and different resolutions, the SLI and 3D-PLI measurements should be integrated in the same setup. The results obtained from the SLI measurements can be used to validate and correct the nerve fiber orientations reconstructed by 3D-PLI and notably improve the nerve fiber tractography.

V. CONCLUSION
The detailed reconstruction of crossing nerve fibers in the brain poses a major challenge for many neuroimaging techniques. Here, we introduced Scattered Light Imaging (SLI), a neuroimaging technique that allows to reconstruct the individual orientations of (crossing) nerve fibers with micrometer resolution over a field of view of several centimeters. In measurements of different species (human, vervet, rat), we verified that SLI yields reliable fiber orientations, in particular in regions with crossing nerve fibers. SLI can be easily combined with 3D-Polarized Light Imaging (3D-PLI). This allows to benefit both from the high-precision three-dimensional nerve fiber orientations from 3D-PLI and the full structural information from SLI, greatly enhancing the reconstruction of the highly complex nerve fiber architecture in the brain. In conclusion, SLI is a major step towards a full reconstruction of the brain's nerve fiber network with micrometer resolution, especially in regions with crossing nerve fibers.

DATA AND CODE AVAILABILITY STATEMENT
All images used in this study (original image stacks and generated parameter maps from SLI, as well as the transmittance, direction, and retardation images from 3D-PLI) are available in original resolution on the EBRAINS data repository [27]. The code used in this study has been published as open-source software SLIX (Scattered Light Imaging ToolboX) and is available on GitHub (https://github.com/3d-pli/SLIX). optic tracts with different crossing angles, and out-of-plane fibers in the cingulum/fornix of a coronal vervet monkey brain section. The SLI measurements were performed as described in Sec. II C with a mask with rectangular holes and 15 • -steps. For each region marked in Fig. 14(a), the line profiles were computed for each image pixel and the peaks in the line profiles were determined as described in Sec. II E for different minimum prominence values: {0 %, 1 %, . . . , 20 %} of the total signal amplitude. The computed number of peaks was then compared to the expected number of peaks for each image pixel in the different types of regions. Figure 14(b) shows the fraction of correctly determined regions (i. e. the number of pixels where the computed peak number corresponds to the expected peak number, divided by the total number of pixels in the considered type of region) for different minimum prominence values. In regions with parallel in-/out-of-plane fibers (green/blue curves) where one or two large peaks are expected, the fraction of correctly determined regions becomes largest when only peaks with larger prominence values are considered: when only peaks with a prominence values ≥ 20 % are considered, 99.7 % (95.6 %) of the regions with in-plane (out-of-plane) fibers are correctly determined. In regions with three crossing fiber bundles/layers (orange curve) where the distance between the resulting peaks is small, the fraction of correctly determined regions becomes largest (55.7 %) when all peaks (i. e. with prominences ≥ 0 %) are considered.

CREDIT AUTHORSHIP CONTRIBUTION STATEMENT
In regions with two crossing fiber bundles/layers (magenta curve) where the distance between the resulting peaks is larger, the fraction of correctly determined regions becomes largest (66.7 %) for a minimum prominence value of 5 %: For smaller values, too many (i. e. also erroneous) peaks are considered. For larger values, not enough (less than four) peaks are considered.
To find a good compromise between regions with different fiber constellations, the curves for in-plane (crossing) and out-of-plane fibers were averaged with equal weight (see black curve in Fig. 14(b)). For a minimum prominence value of 8 %, the fraction of correctly determined regions becomes largest (69.4 %). For this value, 96.5 % of in-plane parallel fibers, 79.8 % of out-of-plane parallel fibers, 65.2 % of two crossing fiber bundles, and 36.3 % of three crossing fiber bundles are correctly determined. The low value for the three crossing fiber bundles is due to the 15 • -discretization of the SLI profiles and the relatively small amplitude of the resulting peaks (cf. Fig. 7 on the right). Most of these regions (96.7 %) are still determined as crossing fibers (4, 5, or 6 peaks, see orange dotted curve). Therefore, when evaluating line profiles obtained from SLI measurements, we used a minimum prominence value of 8 % for the peak detection, as described in Sec. II E.

Appendix B: Correction of determined peak positions
During an SLI measurement, the mask is rotated in steps of ∆φ = 15 • or 22.5 • (see Sec. II C), yielding a discretized line profile I(φ) for each image pixel. The peak positions in the discretized line profiles usually do not correspond to the underlying nerve fiber orientations in the measured tissue sample (cf. Fig. 15(a)). Assuming that most peaks are symmetric, we can use the geometric centers (centroids) of the peaks as an estimate for the "real" peak positions of the ideal, non-discretized line profile:φ . (B1) In this notation,φ uncorr is the peak position in the discretized line profile,φ corr is the estimated peak position of the non-discretized line profile, φ L and φ R are the left and right limits used for the centroid calculation, and I(φ) is the measured light intensity for the illumination angle φ.
To find the optimum values for the left and right limits, we compared different methods (see Fig. 15(b)): (i) No correction: The peak positions were determined as described in Sec. II E, without any correction (blue circles).
(iii) Centroid between two minima: The centroids of the determined peaks (pink arrow) were computed with eq. (B1), using the nearest left and right minima (with prominence ≥ 8 %) as left and right limits. When there was no such minimum on one side, a distance of two steps (30 • ) was chosen as limit.
(iv) Centroid of peak tip: All 15 • -steps in the discretized line profile were divided into 100 equidistant points and linearly interpolated between the points, yielding steps of ∆φ = 0.15 • . From each determined peak position, the points in the line profiles were followed to the left and the right until the distance to the peak becomes larger than two steps (30 • ), a minimum (with prominence ≥ 8 %) is reached, or the intensity value at one point reaches a value below {0 %, 1 %, . . . , 20 %} of the total amplitude (I max − I min ) of the line profile. The centroid of the peak (green arrow) was computed with eq. (B1), using these limits (here: for a peak tip height of 5 % of the total signal amplitude).
As we only want to correct for discretization artifacts, the peak positions were shifted by maximally one step size: To identify the optimum correction method, we considered 232 line profiles obtained from scatterometry measurements of 9 different brain tissue samples (see all circles in Fig. 16). The scatterometry profiles were obtained by computing the (smoothed) polar integrals of the respective scattering patterns (see Sec. II F). The resulting scatterometry profiles can be considered as non-discretized and be used as model systems to test the different correction methods. To obtain a discretized line profile, the scatterometry profiles were sampled in steps of 15 • (see Fig. 15(a)). In real brain tissue, the fiber orientations are not known with < 15 • precision, i. e. each fiber direction (peak position) is equally probable within 15 • . Although the scatterometry profiles belong to different fiber configurations, the underlying fiber orientations are not uniformly distributed within 15 • . To generate an unbiased, uniformly distributed data set reflecting all possible fiber orientations (samplings) that might occur in an SLI measurement, 15 different starting points {0 • , 1 • , . . . , 14 • } were used for the sampling, yielding 15 discretized line profiles for each scatterometry profile. For all 15 × 232 discretized line profiles, the corrected peak positionsφ corr were computed, using the four correction methods described above (colored arrows in Fig. 15(b)), and compared to the original peak positionsφ of the scatterometry profile (orange circles in Fig. 15(a)). To ensure that only positions belonging to the same peak are compared to each other, we considered only line profiles with the same number of peaks. In total, 8759 peaks in 3301 discretized line profiles could be used for comparison. For 4.7 % (0.4 %) of the line profiles, one (two) peaks were missing after sampling.
The histograms in Fig. 15(c) show the difference between original and corrected peak positions for the different correction methods. All curves are centered around 0 • (the means lie between −0.01 • and −0.65 • ), but they differ in terms of absolute mean and standard deviation: (i) Without any correction (blue curve), the differences are almost uniformly distributed between ±5 • (absolute mean: 3.9 • , standard deviation: 4.8 • ). As every peak is sampled in 15 •steps with 15 different offsets, the resulting peak positions of the sampled line profiles are expected to be uniformly distributed around the original peak position within ±(15 • /2) in case of symmetric peaks. In case of asymmetric peaks, the differences can be larger. (ii) When correcting the peak positions by computing the centroid between neighboring minima (pink curve), the kurtosis increases, but the corrected peak positions can still differ a lot from the original peak positions (absolute mean: 4.9 • , standard deviation: 6.5 • ). (iii) When computing the centroid over ±2 steps (orange curve), the distribution becomes much narrower and resembles a normal distribution (absolute mean: 2.6 • , standard deviation: 3.7 • ). (iv) When computing the centroid of the peak tip with a height of 5 % of the total signal amplitude (green curve, cf. Fig. 15(b) in green), the distribution has an even smaller standard deviation (absolute mean: 1.6 • , standard deviation: 2.5 • ). Figure 15(d) shows the absolute mean and standard deviation for method (iv) for different peak tip heights ({0 %, 1 %, . . . , 20 %} of the total amplitude). If the height becomes too small, there are not enough points of the sampled line profile included in the centroid calculation (for 0 %, the peak position is not corrected at all). If the height becomes too large, the peak tip becomes less symmetric (in case of closely neighboring peaks) and the centroid calculation does not yield the correct peak position. For a peak tip height of 6 %, the standard deviation becomes smallest (2.4 • ) with an absolute mean of 1.5 • .
When evaluating the SLI line profiles, the determined peak positions were therefore corrected with correction method (iv), using a peak tip height of 6 % of the total signal amplitude.
Appendix C: Brain tissue samples Figure 16 shows the brain tissue samples (average transmitted light intensity of SLI measurements) that were used to compare the SLI profiles to scatterometry profiles of the same measured tissue spots. The scatterometry measurements were performed as described in Sec. II F (with beam diameter ∅ = 1.12 mm, numerical aperture NA = 0.4, and exposure time t = 30 ms, for more informations see [12]). All tissue spots that were used for the comparison are marked by a circle (diameter: 1.12 mm). The colored circles (in green, blue, black, orange, magenta) indicate the tissue spots for which the differences between SLI and scatterometry profiles have been computed (see Sec. III A). In spots surrounded by a green circle, the differences are small, in spots surrounded by a magenta circle, the differences are large. Table I lists the sample properties and measurement parameters for all investigated brain tissue samples.  Figure 16. Brain tissue samples used for comparing SLI and scatterometry profiles (coronal/sagittal vervet brain sections, two/three crossing sections of human optic tracts). The images show the average transmitted light intensity of the SLI measurements, labeled by the sample, section number, and measurement ID #. Table I contains more details about the measurement, like the rotation angle by which the SLI images were rotated with respect to the scatterometry measurement.
The circles indicate the positions of the tissue spots (with 1.12 mm diameter) that were measured by scatterometry. All 232 spots were used to determine the optimum correction method to reduce discretization artifacts in SLI profiles (see Appx. B). The colored circles (in green, blue, black, orange, magenta) indicate how good the SLI profiles correspond to the scatterometry profiles in the respective tissue spot (green: small differences, magenta: large differences, see Fig. 6). In the vervet brain sections, anatomical regions are labeled (cr = corona radiata; cg = cingulum; cc = corpus callosum; f = fornix).   Table I. List of sample properties and measurement parameters for all investigated brain tissue samples: sample, section number, section thickness (T ), measurement ID (#), dates of SLI and scatterometry measurements (in days after tissue embedding), rotation angle (ρ) by which the SLI images were rotated with respect to the scatterometry measurement (Fig. 16); settings used for the SLI measurements (see Fig. 2): mask (circular/rectangular), focal length of camera objective (f ), distance between light source and sample (h), distance between sample and camera objective (l), exposure time (t), and pixel size in object space (px). The asterisk indicates that the given number of days are the days after the revitalization of the tissue (not after the initial embedding). The coronal rat brain section was obtained from a Wistar rat (male, 3 months old), the vervet brain sections from two African green monkeys (Chlorocebus pygerythrus, male, 1 and 2.4 years old), and the sections of the optic chiasm from a human (female, 74 years old). Figure 17 shows the (normalized) transmittance and retardation images obtained from a 3D-PLI measurement of section no. 36 of the human optic chiasm (left: whole section; right: sections of optic tracts crossing each other). Figure 18 shows the SLI profiles for the two crossing sections of optic tracts drawn along a circle (for every 24 th image pixel). The SLI profiles were drawn along circles (for every 24 th image pixel). The profiles were shifted by 90 • beforehand so that the position of the peak-pairs corresponds to the in-plane orientations of the fibers. The zoomed-in region shows the transition between parallel and crossing fibers, which is visible in the morphology of the shapes. The image in the background shows the average intensity of the SLI profiles for each image pixel. Figure 19 shows the SLI profiles for three crossing sections of optic tracts obtained from an SLI measurement with 15 • and 5 • -steps. Figure 20 shows the average peak width and peak prominence of the coronal rat brain section (Fig. 11) in comparison to the retardation. . Average peak width and prominence of the coronal rat brain section (Tab. I#78). For peaks with a prominence ≥ 8 % of the total signal amplitude, the peak width and prominence were computed from the normalized SLI profiles for each image pixel (as described in Secs. II E and III B 1) and averaged for each SLI profile. The 2D histograms show the average peak width and prominence plotted against the retardation for each image pixel. The graphs below show the standard deviation σ plotted against the retardation.