An unambiguous cloudiness index for nonwovens

Cloudiness or formation is a concept routinely used in industry to address deviations from homogeneity in nonwovens and papers. Measuring a cloudiness index based on image data is a common task in industrial quality assurance. The two most popular ways of quantifying cloudiness are based on power spectrum or correlation function on the one hand or the Laplacian pyramid on the other hand. Here, we recall the mathematical basis of the first approach comprehensively, derive a cloudiness index, and demonstrate its practical estimation. We prove that the Laplacian pyramid as well as other quantities characterizing cloudiness like the range of interaction and the intensity of small-angle scattering are very closely related to the power spectrum. Finally, we show that the power spectrum is easy to be measured image analytically and carries more information than the alternatives.


Introduction
Nonwovens, papers, and felts are used in a wide variety of fields, ranging from filtration of gases or liquids over thermal insulation and soundproofing to hygiene and sanitary articles, see e.g. [1]. They consist of fibers of limited length or filaments, that are more or less randomly distributed and forming the solid matter of a macroscopically homogeneous porous material.
The macroscopic properties of nonwovens are significantly determined by the specific weight (i.e. the mean weight per unit area, also called the nominal grammage) and the spatial weight distribution. Moreover, the nonwoven texture -also called 'nonwoven cloudiness' or 'chart cloudiness' -is of high practical relevance. In paper industry, 'formation' and 'flocculation' are used as synonyms for 'chart cloudiness'.
Cloudiness influences filter properties like flow rate, particle retention efficiency, wet strength, porosity, dust holding capability, heat transfer, and sound attenuation [2]. Thus, quantifying cloudiness is important for industrial quality control as well as for the development of new filter materials and manufacturing technologies. Robust estimation of a cloudiness index could help in optimizing processing parameters.
Holding a sheet of nonwoven up against light yields a visual impression of cloudiness as the spatial distribution of darker or brighter regions. Unfortunately, the term 'cloudiness ' has not yet been defined in any industry standard for nonwovens. The industrial standard [3] on paper, board, pulps, and related terms describes 'formation' only very roughly as 'manner in which the fibers are distributed, disposed, and intermixed to constitute the paper' and 'look-through' as 'structural appearance of a sheet of paper observed in diffuse transmitted light'.

arXiv:2201.02011v2 [cs.CV] 11 Apr 2022
The intensive efforts that have been taken to characterize cloudiness become evident by the large variety of dedicated publications, see Section 2 below. Nevertheless, so far, no characteristic or index for quantifying cloudiness could be agreed on. Fixing 'cloudiness' in an industry standard seems to be out of reach. This is due to the wide variety of nonwoven fabrics, the rapid development of physical and optical testing methods, a wide range of image analysis methods that are difficult to survey and, finally, the understandable interest of optical inspection system manufacturers to keep their specific concepts secret.
Here, we follow [4][5][6][7][8][9][10][11][12] in using the power spectrum or -equivalently -the correlation function to measure cloudiness. The power spectrum and related mathematical concepts are carefully introduced in Section 3. Moreover, the power spectrum of the modified Bessel correlation function is found to be a flexible parametric model for the rotation averages of the power spectra of nonwovens. For the sake of better interpretability, we introduce an explicit structure model for nonwovens in the subsequent Section 4. Estimation of the suggested cloudiness index based on images is demonstrated using microscopic images of three nonwoven samples in Section 5.
Finally, in Section 6, we compare our approach to a variety of alternatives: the range of interaction, the Laplacian pyramid, and the intensity of small-angle scattering. The range of interaction promoted as homogeneity measure in [13,14] is just the integral of the correlation function. The Laplacian pyramid decomposes the image of the structure into 'several scales', evaluates the degree of homogeneity on each scale, and finally derives an overall degree of homogeneity [15][16][17]. Small-angle scattering (SAS) has been considered as early as in the 1970s, too, to characterize the cloudiness of paper formation [18]. This is an obvious choice, as the scatter intensity is closely related to the power spectrum. We recall this relation in Section 6.2, as recent devices make SAS attractive for industrial quality control. In general, we contribute to a better overview by showing in Section 6 how closely the quantities listed here are related. Moreover, we describe under which conditions cloudiness can be characterized completely.

Concepts of cloudiness of paper and nonwovens
There is a vast variety of publications on characterizing the cloudiness of nonwovens and formation of paper based on transmission of diffuse light, where a transmission light table is applied to ensure homogeneous illumination [6,[19][20][21]. See [22,23] for a comprehensive survey on literature from the late 1980s and [24] for more recent developments.
Robertson's mean flock size [25] is an intuitive characteristic for paper formation. However, there is still no convincing method for segmenting flocks in transmitted light images available. Practically more useful methods are based on estimating the variance of the pixel values or, more general, the co-occurrence matrix of the image data [6,[26][27][28]. In [29], a 'uniformity index' of texture is suggested based on Poisson statistics for the centers of paper flocks (objects). This approach seems a bit inconsistent however. On the one hand, the flock centers cannot be detected robustly. On the other hand, the index is finally computed based on local area fraction variation in a binarized image not using the flock centers at all.
In the carefully written monograph [30], an index of cloudiness is defined as the ratio of two variances: For the nonwoven sample and a reference model, local grammage is averaged over a square of edge length 1 mm. The ratio of the variances of these two characteristics yields the index. This approach dates back to earlier works, see e.g. [31] frequently cited in recent literature like [32], too. Farnood's approach [33,34] is similar, modeling the fluctuation of the local grammage by a Poisson shot noise process (i.e. a dilution process) of spherical flocks. In [35], Farnood's flocculation characteristic is derived from the so-called power spectrum, also known as 'power spectral density' or 'Bartlett spectrum', of the pattern in transmission images.
Woven textiles usually feature a periodic pattern. Hence, applying Fourier methods for quality inspection seems natural, see e.g. [36,37] and references therein. Periodicity of the pattern corresponds to characteristic peaks of the power spectrum, and slight deviations from perfect periodicity can lead to smearing (i.e. broadening) of these peaks [38]. Following this line of thought, transferring these ideas to non-periodic but macroscopically homogeneous patterns and characterizing nonwoven fabrics with the help of spectral analysis seems obvious. In [4-9, 11, 12], the power spectrum (or equivalently the correlation function) is suggested to measure cloudiness, see also Section 2.2 in [10].
Further approaches are based on modeling random structures by Markov random fields and decomposing the image of the structure into 'several scales', evaluating the degree of homogeneity on each scale and computing an overall degree of homogeneity. In [17], this approach is applied originally on the structure of nonwovens. These 'several scales' are also known as the Laplacian or Gaussian-Laplacian or Burt-Adelson pyramid for image data [15][16][17] originally introduced for image compression [39].
Scharcanski [40] suggests a wavelet transform for evaluating sheet formation and cloudiness of nonwovens. Replacing the Fourier transform by a wavelet transform for the analysis of the local grammage in the frequency domain can have computational advantages over the Fourier transform, especially when recursive methods are applied for calculation. Wavelets (in particular Haar wavelets) are very well suited for analyzing piece-wise constant functions like microscopic images of the microstructure of multi-phase materials. The local grammage of nonwoven however is not piece-wise constant. Nevertheless, there are wavelet transforms for a variety of purposes and the choice of one of them influences the estimation of the spectral density. For example, the hybrid method presented in [41] -a combination of a discrete wavelet transform and Gabor filter banks -aims at the segmentation of clouds in a pre-specified size range.
Mathematical modeling of nonwoven structures on a mesoscale also yields a deeper understanding of the phenomenon of cloudiness [6,8,32,[42][43][44]. Model parameters can be interpreted as characteristics of cloudiness (i. e. the uniformity index or the degree of homogeneity of the pattern), given they can easily be estimated from image data. In [12], Gaussian random fields are used for modeling paper structure and a Fourier approach is applied for computing these characteristics, following suggestions from [11,45].

Power spectrum
Obviously, the accuracy of quantities estimated by optical methods always depends on the imaging conditions. Here, we introduce the quantities independent of the imaging parameters. Their limits are discussed later.
Assume the structure of a nonwoven to be captured in a transmitted light image using a point, line, or area detector. In this case, the local grammage of the nonwoven and the detected signal can be expected to be closely related. In fact, often, these two are not really distinguished. Nevertheless, we want to emphasize the difference: Let R 2 be the projection plane orthogonal to the optical axis, and let γ be the linear mass attenuation coefficient of the solid matter. Then the intensity I(x) of the transmitted light at a point x = (x 1 , x 2 ) in the plane R 2 can be related to the local grammage w(x) by Lambert-Beer's law I(x) = I 0 exp(−γw(x)), where I 0 is the initial intensity of the applied light. Thus, the local grammage is See [11,46,47] for the computation of the local grammage from the absorption of visible light. The use of β-and soft X-radiation is suggested in [48][49][50][51], and [52], respectively, and the application of electron beam transmission for estimating local grammage is studied in [53]. Finally, the influence of the choice of radiation on transmittance is investigated for nonwovens in [54] and for paper in [55,56]. Nonwovens feature a random, macroscopically homogeneous structure, i.e. the distribution of the structure is invariant with respect to translations in R 2 . Hence, the local grammage w(x) is a random function on R 2 that is macroscopically homogeneous as well. We denote byw = Ew(x) and σ 2 w = E w(x) −w 2 the mean grammage (also called the 'nominal grammage' or the 'base weight') and the variance, respectively, where E denotes the expectation. Furthermore, let be the normalized local grammage. Then the correlation function k(x) of w(x) is simply given by where the right-hand side of this equation is independent of the position y because of the macroscopic homogeneity of f (x). Finally, the power spectrumk(ξ) of f (x) is defined as the 2-dimensional (2D) Fourier transform of the correlation function k(x),k Vice versa, k(x) is the Fourier co-transform ofk(ξ), The variable ξ = (ξ 1 , ξ 2 ) is known as the circular frequency. Note that Equation (3) is motivated by the projection-slice theorem of the Fourier transform stating that the 2D Fourier transform of the projection of a 3D structure is a slice of the 3D Fourier transform of the structure [57]. Analogously, the 2D power spectrumk(ξ) of f (x) is a slice of the 3D power spectrum of the 3D nonwoven structure. Often, nonwoven structures are not only macroscopically homogeneous, but also isotropic, i.e. the distribution of the structure is invariant with respect to rotations in R 2 . Then the correlation function k(x) and the power spectrumk(ξ) are isotropic as well, which means that they depend only on the radial coordinates r = x = x 2 1 + x 2 2 and = ξ = ξ 2 1 + ξ 2 2 , respectively. So there are functions k 1 (r) and respectively, related to each other by the Bessel transform where J 0 (r) is the Bessel function of the 1st kind and of order 0. See e. g. [58] for a sound introduction to the Fourier transform and related topics. The Bessel transform, also known as the Hankel transform, is inverse to itself, i. e. it is an involution. Consider for example the modified Bessel correlation function with the parameters λ > 0, ν > 0 and the modified Bessel function of the 2nd kind K ν (r) of order ν. Its Bessel transform iŝ see [59].
Graphs of the above power spectra are shown in Figure 2. The power of low frequencies decreases with increasing λ. From Equation (4) one immediately gets i. e. the integral of the power spectrum is constant. From this kind of 'energy conservation law' it follows that decreasing the power of low frequencies (i. e. decreasing cloudiness) leads to increasing the power of higher frequencies. For high frequencies ( λ), the power spectrumk 1 ( ) given by (7) decays as fast as −2(ν+1) . It is often suggested to characterize cloudiness by the variance σ 2 w of the local grammage w(x) or its distribution function. We would like to stress that this is impossible. A simple scaling w(cx) of the local weight of a nonwoven with a constant c > 1 results in the nonwoven appearing more homogeneous due to the strongly decreasing power spectrum. In fact, because of the inverse scaling law of the Fourier transform, the power spectrum of w(cx) isk(ξ/c)/c 2 , andk(ξ/c)/c 2 <k(ξ) for small frequencies, whereas variance and distribution of the local grammage are not affected by the scaling.

Modeling cloudiness
The function defined by (7) is certainly sufficiently flexible to be fit to measured power spectra. However, the parameters λ and ν are hard to interpret. In order to clarify the relation of nonwoven structure and power spectrum, we model the nonwoven structure and derive a formula for the power spectrumk 1 ( ).
We first model the local grammage w(x) of a nonwoven. Let ϕ(s), 0 ≤ s ≤ 1 be a parametric function in the projection plane R 2 . Then the set γ = {ϕ(s), 0 ≤ s ≤ 1} is a curve in R 2 . We assume γ to be of finite length. Its dilation C = γ ⊕B R with the ball B R of radius R yields a filled tube of circular cross-section, see Figure 3a. Now C can serve as a model of a curved fiber of constant thickness 2R. We assume C to be morphologically regular. This requires in particular the curvature of the core γ to be less than 1/(2R). The orthogonal projection p(x) of C on R 2 is the chord length of C at x ∈ R 2 , i. e. the length of the intersection of C with the straight line orthogonal to the projection plane R 2 and hitting the point x, see Figure 3b. Figure 3 Modeling the local grammage of a straight fiber: a) the core γ and the contour of the fiber model γ ⊕ B R , b) the projection of the fiber model on the plane R 2 , c) the core convolved with the filter mask κ(x), d) the difference between projection and filtering. The filter response in x is the length of the core γ intersected with the circle divided by the area of the circle. Now we consider a macroscopically homogeneous and isotropic random system of fiber cores. That is, parametric functions ϕ j (s), 0 ≤ s ≤ 1, j = 1, 2, . . ., where the points ϕ j (0), j = 1, 2, . . . form a macroscopically homogeneous Poisson point field. Finally, the local grammage w(x) is modeled by the macroscopically homogeneous and isotropic random field w m (x) defined as the sum of the orthogonal projections p j (x) of the fiber models C j with the cores γ j = {ϕ i (s), 0 ≤ s ≤ 1} times the specific mass density ρ of the fiber material, See Figure 4 for examples of random fields w m (x) generated by segments γ j with directions uniformly distributed in [0, π) and exponentially distributed lengths. Letw m and σ m be the mean and the standard deviation, respectively, of the random field w m (x). Then the normalized local grammage converges to a Gaussian random field of mean 0 and standard deviation 1 as the mean number of curves per unit area N A goes to infinity due to the central limit theorem [60]. That is, for sufficiently large N A , the random field w m (x) is characterized exclusively by the meanw m , the standard deviation σ m , and the power spectrumk 1,m ( ). As a consequence, we can assert that the power spectrumk 1 ( ) of w(x) characterizes cloudiness uniquely if w m (x) is a suitable model for the local grammage w(x) of the nonwoven.
Aiming at a closed formula for the power spectrumk 1,m ( ) of the model w m (x), we now replace the projections p j (x) of the sphero-cylinders C j by the cores γ j convolved with a suitable filter kernel κ(x). Then p j (x) is approximated by the filter response in the point x ∈ R 2 , see Figure 3c. Practically, this is achieved by filtering the Bresenham line of the segment γ i with the mask of an isotropic mean value filter: The Bessel transform of the radial function κ 1 (r) of κ(x) is given bŷ with the Bessel function J 1 (r) of the first kind and of order 1 [58]. The pair correlation function of the straight line system is see [61]. Its Bessel transform is with the Dirac function δ(x), see [12]. Hence, the Fourier convolution theorem yields that the power spectrumk 1,m ( ) of w m (x) generated by a system of straight fibers of radius R can be approximated aŝ where ψ(λ, R) > 0 given by is a correction factor ensuring thatk 1,m ( ) is actually the Bessel transform of a correlation function k 1,m (r), in particular k 1,m (0) = 1. For high frequencies, the right-hand side of (9) resembles qualitatively the square of sine or cosine waves whose amplitude decays as −4 , and whose first factor is as (7) for ν = − 1 2 . Unfortunately, until now, there is no closed form known for ψ(λ, R). For numerical values see Figure 5. The parameter λ of the exponential distribution is the inverse mean segment length, i. e. 1/λ + 2R is the mean fiber length. Figure 6 shows the power spectrumk 1,m ( ) given by Equation (9). The power of low frequencies (assigning small chart cloudiness) decreases with increasing mean fiber length. The limiting case of infinitely long fibers (λ = 0) yields the smallest cloudiness. The right side of the approximation in (9) converges to the left side as λ → 0 (infinite fiber length, filaments), since the projection and filter response differ only at the fiber ends, see Figure 3d.
Regime changes for the rate of the decay of the correlation function may be observed for fixed fiber length. However, this does not apply here, as fiber lengths are randomly distributed.

Estimation of the power spectrum
The accuracy of quantities estimated by optical methods always depends on the stimulating light (intensity, wavelength, coherence), the size of the investigated  sample (i.e. the window size), and the type of detector (pixel size, modulation transfer function). These metrological parameters have to be adjusted to allow for reliable measurements in any case. Here, we assume sufficient quality of the images throughout. We discuss in this section however details on how to estimate the power spectrum from image data unambiguously.
The pixel value g(x) of a gray-tone image at the pixel position x does never match exactly the intensity I(x) of the transmitted light. The intensity I(x) can be deduced from g(x) only if the imaging conditions are known, including the photo-detector spectral response, its dark current density, and the sensor conversion gain. Thus, estimating mean grammagew and standard deviation σ w objectively requires tedious calibration. However, in order to quantify cloudiness, only the power spectrumk(ξ) of the relative local grammage f (x) (2) has to be estimated.
In the following we assume g(x) > 0 for all pixel positions x of the image. In practice this means that the illumination at the microscope (intensity, exposure time) must be chosen such that the gray-tone is positive almost everywhere. Whenever g(x) and I(x) are approximately proportional, then f (x) can be estimated from cf. Equation (1), where µ and σ are the mean and the standard deviation of the function ln g(x) for all x, respectively.

Estimation based on a finite field of view
In this section, we deal with the practical implications of the fact, that the imaged field of view is always finite. Let W ⊂ R 2 be an image frame or field of view of area A. Let 1 W (x) denote the associated indicator function, 1 W (x) = 1 for x ∈ W , and 1 W (x) = 0 otherwise.
The windowed function f W (x) = f (x)1 W (x), i. e. the image f W (x) of f (x) in W , is square integrable. Thus the correlation of f W (x) with itself is well defined and can be written as the convolution of f W (x) with the reflection f * for all x in the interior W • of W . Here * denotes convolution. Then the correlation function k(x) of f (x) can be rewritten as is the so-called window function. If x is small compared to the window size, then c W (x) is approximately the frame area, c W (x) ≈ c W (0) = A. The use of a Fast Fourier Transform (FFT) induces additionally a periodic extension of f W (x), which together with the Fourier convolution theorem leads to the approximation for estimating the power spectrumk(ξ) by image analysis. See [12,62,63] for exact estimation ofk(ξ).

Averaging over several images
The limited field of view implied by the imaging method can be compensated by taking images of several disjoint fields of view. The linearity of the Fourier transform and the macroscopic homogeneity of nonwovens result in the following important fact: For several images f W,j (x), j = 1, . . . , m taken at disjoint positions of W , the power spectrum of the mean f W (x) = 1 m m j=1 f W,j (x) coincides with the power spectra of the individual images. That means, we can first pixelwise average the f W,j (x) and then estimate the power spectrum of the mean.
As a consequence, instead of m time-consuming FFTs, only one FFT is needed for estimating the power spectrum of f W (x). Thus, the number m of images can be increased arbitrarily while the computing time remains almost constant.

Application to microscopic images of nonwovens
We consider three nonwoven samples. Our first example is a polypropylene nonwoven used as one layer in face masks. Specimen FL1 is produced by Sandler AG, Schwarzenbach/Saale, Germany. It has a nominal grammage of 100 gm −2 , thickness of 0.98 mm at 0.05 Ncm −1 , initial pressure drop of 163 Pa, and air permeability of 98 Lm −2 s −1 at 200 Pa according to the industrial standard [64]. About 35 % of the fibers have a thickness d smaller than 1.25 µm, 33 % with 1.25 µm≤ d <2.5 µm, 29 % with 2.5 µm≤ d <5 µm, and 3 % with d ≥ 5 µm.
The structure is imaged in transmission mode using the Leica stereo light microscope MZ16. The detailed image conditions are: planapochromatic 1×-objective, numerical aperture 0.14, working distance 55 mm; light source Phlox-LedRGB-BL-100 × 100-S-Q-IR-24V, white, exposure time 2 ms; detector Basler acA4112-8gm with 4 096 × 3 000 pixels; sensor IMX304 CMOS, GigE, mono. Figure 1 shows m = 4 gray-tone images of the nonwoven taken at the same lateral resolution (2 048 × 2 048 out of totally 4 096 × 3 000 pixels of size 3.367 µm). Figure 7b shows the pixelwise mean f W (x) of the f W,j (x) obtained from the four images shown in Figure 1.   Figure 7a. The isotropy of FL1 is passed on tok(ξ), which obviously is invariant with respect to rotations about the origin. Thus the radial functionk 1 ( ) contains the same information ask(ξ). It is therefore obvious to determinek 1 ( ) as the rotation average of the estimate ofk(ξ), which also reduces the statistical error. Figure 9 shows the power spectrumk 1 ( ) of the nonwoven FL1 estimated based on images at varying lateral resolution (Figures 7a to 7d). Additionally, the model function (7) is fit to the estimated power spectrum.
Total area and lateral resolution have to be chosen such thatk 1 ( ) can be estimated in a given spectral band with a sufficiently small total error. The statistical error increases for decreasing , and it diverges as ↓ 0. For higher frequencies the systematic error increases as fine details of the nonwoven structure are not resolved anymore. Thus, the lateral resolution and the size of one field of view determine the   Table 1 for the parameters λ and ν of the adapted model function from Equation (7).

Cloudiness index
We now derive the cloudiness index CLI as the total power in a medium frequency band [ 0 , 1 ] for 0 < 0 < 1 , Equation (8) yields that 0 ≤ CLI ≤ 1. Thus, the CLI can be given in percent.
The thus defined cloudiness index depends solely on the chosen frequency band which is determined to a high extent by the imaging conditions as explained in the previous Section 5. 3. The lower bound 0 is bounded by the inverse length of the diagonal of the imaged field of view. The upper bound 1 is restricted by the lateral resolution, that is the doubled inverse pixel size. Decreasing 0 or increasing 1 requires enlarging the field of view or increasing the resolution, respectively. Moreover, an expansion of the range does not necessarily lead to a better assessment of cloudiness due to the power conservation law. For industrial quality control, this range has to be specified according to the application case.
For application example FL1, we use 0 = 0.02 µm −1 and 1 = 0.10 µm −1 , corresponding to wavelengths between 62.8 µm and 314.2 µm, and obtain a cloudiness of CLI = 25.7 %. See Table 1 for the corresponding results for examples FL2 and FL3. Obviously, computing the cloudiness index CLI from the rotation meank 1 ( ) Table 1 Estimates of the parameters λ and ν for the fitted function (7), the cloudiness index CLI, and the range of interaction RI. The specified error bounds are derived as the standard deviation of the mean of the estimates based on the individual images. Total sample area is A. removes the information on anisotropy of a nonwoven. A cloudiness index including directional information could be derived by considering the power spectrumk(ξ) as a function of the polar coordinates r = ξ and ϕ = arctan ξ2 ξ1 for ξ 1 = 0 and ϕ = π sgn ξ 2 , otherwise, where sgn ξ 2 is the sign of ξ 2 . Restricting the angle ϕ to a sub-range can then reveal directional dependency of the structure. For example, for sample FL2 from Figure 10, the cloudiness index is 17.1 % for the sector |ϕ| ≤ π 4 (i. e. for directions ϕ close to the processing direction of the nonwoven) and 45.5 % for |ϕ| > π 4 . Alternatives are e. g. the orthogonal projection of the power spectrumk(ξ) onto the x-and y-axis, respectively, or fitting a 2D function generalizing (7) tok(ξ). In this case, the parameters λ and ν depend on the angle ϕ. Characterization of anisotropy is however subject of further research.
In inspection systems for industrial quality control, line cameras are often used to scan the material continuously in the production line, see e. g. [65,66]. The focus is on defect detection. However, cloudiness of nonwovens and formation of paper can be measured with such systems, too. The 1D correlation function of the camera signal is a linear section through the 2D correlation function k(x) through the origin and in line direction. In the isotropic case, this 1D correlation function coincides with the rotation mean k 1 (r). For the power spectrum, however, this is not true. The 1D power spectrum of the camera signal obtained with the 1D Fourier transform is, in the isotropic case, the cosine transform of k 1 (r) but not the Bessel transform.

Relation to other quantities characterizing homogeneity
Many measures of structural homogeneity or deviation from it have been suggested, see Section 2. Here, we shed light on how they are interrelated with particular focus on our CLI based on the power spectrum and on the Laplacian pyramid as promoted in [15][16][17].
First, we recall that the correlation function k(x) and the power spectrumk(ξ) can be converted into each other by the Fourier transform. That is, these two functions carry the same information about the local grammage w(x). Therefore it seems a matter of taste, whether the cloudiness of a nonwoven or a paper is characterized based on k(x) ork(ξ).
However, there are tangible advantages of characterizing cloudiness based onk(ξ). On the one hand, cloudiness is the total power of a medium frequency band. On the other hand, there are various approaches to estimatek(ξ) directly by scattering. Formally, to estimate the correlation function k(x) effectively, one has to go via frequency space anyway. This involves two time-consuming Fourier transforms, while only one suffices for estimatingk(ξ).

The range of interaction
In [13] and [14] the range of interaction RI, i.e. integral of the correlation function, is suggested as a cloudiness index for nonwovens, see also [67]. Analogously to Equation (8), it follows from Equation (3) that RI is up to a constant the power spectrum at 0, For example, the range of interaction for the modified Bessel correlation (6) is From the parameters of the fitted function in Figure 9, we easily obtain the estimate RI = 0.046 mm 2 for the nonwoven FL1, see also Table 1 for the other examples. Estimating the range of the correlation is apparently very simple as demonstrated in this example. However, the errors should be taken into account. Unavoidable large errors in estimating the correlation function k 1 (r) for r near the maximum width d max of the window W induce large errors in RI when using Equation (14) as the basis for estimating RI. Moreover, for r ≥ d max , k 1 (r) cannot be estimated at all. This effect can also be seen in Equation (15). Direct estimation ofk(0) is impossible. Determiningk(0) by extrapolating a fitted function instead bears the risk of inducing another error due to inadequate fitting. Finally, it is not clear, whether RI really captures the cloudiness of a nonwoven: For very small , the power spectrumk 1 ( ) represents nothing but a gradient in the nonwoven structure or a shading of the image e. g. due to uneven illumination of the field of view.

Small-angle scattering
Small-angle scattering has a long tradition to measure cloudiness, in particular for paper. In the past, its practical use was hampered by the high operating expenses. Currently, this obstacle weakens however with cheaper and more robust hardware rendering small-angle scattering more attractive. In this section, we outline briefly how to estimate the power spectrum by small-angle scattering of partially coherent light. In theses experiments, the scattering intensity I s is collected depending on the so-called scattering angle 2θ, where θ is the Bragg angle [18].
Let the light have wave length λ. The scattering intensity is a function of the magnitude q = (4π sin θ)/λ of the momentum transfer scattering vector, I s = I s (q). For spatially isotropic structures, the scattering intensity I s (q) is, up to a multiplicative constant c, the sine transformation of the spatial correlation seen as a function of the radial coordinate [68,69]. For small-angle X-ray scattering (SAXS), the factor c can be interpreted physically, see e. g. [18].
For flat structures like thin nonwovens, the orthogonal projection of the spatial correlation function onto the xy-plane is approximately the same as the correlation function of the projection w(x) of the structure. As a consequence, q ≈ and c · I s (q) is equivalent to the power spectrumk 1 ( ). That means, c · I s (q) is the Bessel transform of k 1 (r).
SAXS covers a sample area too small to characterize the cloudiness of nonwovens. This disadvantage is easily overcome by modern scanning methods like scanning small-angle X-ray scattering (sSAXS, [70]) enabling fast scanning of larger sample areas. The high exerimental effort of SAXS and sSAXS compared to optical methods limits their use however to nonwovens made of completely opaque fibers.
Attempts to use small-angle scattering of visible light to determine paper formation, see e. g. [71], did not succeed because of the strong absorption.

Laplacian pyramid
Characterizing cloudiness of nonwovens by the Laplacian pyramid is based on filtering images with differences x ∈ R 2 , σ 1 < σ 2 , of Gaussian functions given by In electrical engineering, 1D versions of these filters are called band pass filters. The function h(x) is known as the mask of a difference of Gaussian (DoG) filter, see Figure 15 for an example. In [15,16] a sequence h j (x) = g σj (x) − g σj+1 (x) of DoG filters with σ j = 2 j−1 2 is applied to the image f W (x). The filter response f W * h j (x) is evaluated by the mean square of the pixel values MSP j . In our notation, this is the L 2 norm of the function f W * h j (x) over the window area, for integer j. In [17], a weighted sum of the MSP j is suggested as cloudiness index for nonwovens. The weights are chosen empirically, driven by customer's demands. We use Plancherel's equation to shed light on the close relationship of the Laplace pyramid and the power spectrum. Plancherel's equation says that for all functions of fast decay, the norm of the function is the same as the norm of its Fourier transform [58]. Together with the Fourier convolution theorem this yields where the Fourier transform of the maskĥ j (x) is known as the transfer function of the DoG filter. See Figure 16 for an example ofĥ j (ξ). Define the radial functionĥ 2 1,j ( ) ofĥ 2 j (ξ) byĥ 2 1,j ( ξ ) = h 2 j (ξ) for ξ ∈ R 2 , see Figure 17. The maskĥ j (ξ) is positive for all ξ. Consequently, This means, the square of the transfer function of the filter cuts a frequency band from the power spectrum of f (x). Subsequently, the total power is calculated in this band. The transfer functionĥ 1,j ( ) attains its maximum at max = 1.177/σ j . The square of the norm ofĥ 1,j ( ) is For different j, the functionsĥ 2 1,j ( ) differ considerably in the full width at half maximum and they overlap strongly. Therefore, the MSP j do not contain independent information about the spectrum of the nonwoven. Changing the variances σ 2 j and σ 2 j+1 of consecutive Gaussian functions in h j (x) such that their difference is smaller, does not remedy this. Table 2 yields the estimated MSP j for indices j with position max of the maximum transfer functionĥ 1,j ( ) in the frequency band [ 0 , 1 ] chosen for estimating the cloudiness index CLI. The quantities MSP j and k 1 ( ) are completely different by definition. However, they are actually comparable in the following sense: The ratio 2πMSP j / ĥ j 2 is in the same order of magnitude aŝ k 1 ( max ). For example, for the nonwoven FL1 we have 2πMSP 9 / ĥ 9 2 = 26.3 µm 2 andk 1 ( max ) = 25.7 µm 2 . Finally, note that a spectral band could be separated sharply by the filter mask The function (x) is known as the mask of the Bessel band pass filter, see Figure 18 for an example. Its transfer function is the indicator function of the circular ring with radii 0 and 1 , i. e.ˆ (ξ) = 1 B 1 (ξ) − 1 B 0 (ξ) for ξ ∈ R 2 , where 0 and 1 act like sharp frequency bounds. The result would be equivalent to evaluating cloudiness based on our CLI from (13). Table 2 Estimates of MSP j for j = 9, . . . , 12. The error is the standard deviation of the mean of the estimates from the individual images. The position max of the maximum and the norm ĥ j of the transfer functionĥ j (ξ) are included to facilitate interpretation.

Discussion and conclusion
We define the cloudiness of nonwovens as the total power of the normalized local grammage f (x) in a frequency band. Our cloudiness index CLI is therefore independent of the nominal grammagew and the standard deviation σ w of the local grammage w(x). Estimatingw and σ w requires careful optical densitometry measurements, whereas estimating CLI does not. CLI thus simplifies estimating cloudiness based on light-optical transmission images considerably. Moreover, the associated technical effort is low. As a consequence, CLI is very attractive for quality control in industry. The power spectrum determined by image analysis is closely related to other cloudiness measures, as discussed in Section 6: The range of interaction is -up to the factor 2π -nothing but the power spectrum at the origin (Section 6.1). Physical methods (SAXS, sSAXS) yield the scattering intensity as a function of the magnitude of the momentum transfer scattering vector. This corresponds -up to a multiplicative constant -to the power spectrum (Section 6.2). Finally, the normalized mean square of the pixel values of a transmission image processed with a DoG filter (Laplacian pyramid) is the power spectrum at a fixed frequency (Section 6.3). Hence, all these methods intuitively amount to the same and obtaining comparable cloudiness values just requires careful adjustment of the respective parameters.
The definition of the power spectrum is tied to the definition of the continuous Fourier transform. Commonly used definitions of the Fourier transform differ in the use of the normalization, the sign of the exponents, and using circular frequency instead of frequency. Another source of deviation is the step from continuous to discrete Fourier transform. Here, we understand the discrete FFT as a numerical version of the continuous Fourier transform, see e. g. [72, page 607ff]. Care is needed to estimate the power spectrum based on Equation (12) reproducibly, in particular with respect to normalization as various versions of the FFT are used in the literature and in software packages.

Declarations
List of abbreviations FFT fast Fourier transform SAS small-angle scattering SAXS small-angle X-ray scattering sSAXS scanning small-angle X-ray scattering DoG difference of Gaussian (filter)

Availability of data and materials
The images used in this paper along with a ToolIP [73] graph estimating the CLI will be provided on https://www.itwm.fraunhofer.de/cloudiness.

Funding
This research was supported by the project "Product and quality optimization for protective clothing against infections made from nonwovens" (ProQuIV) within the Fraunhofer Society's anti-Corona programme.       Figures 7a, 11a and 11b, respectively, and the function from (7) fitted to the estimated power spectra. Linear scale.