Error Estimation of Linear Polarization Data from Coronagraphs – Application to STEREO-A/SECCHI-COR1 Observations

We attempt to quantitatively study the uncertainties of the polarized brightness and for the polarization angle which are to be expected for measurements from noisy image detectors in classical coronagraphs. We derive the probability density functions (PDF) which apply to polarization observations with polarization filters at 0∘, 60∘ and 120∘. The noise in the directly observed image intensities is assumed to be normally distributed. We find that for low and medium signal-to-noise ratios the polarized brightness obeys a distribution with a strongly biased mean which, if not taken account of, leads to an overestimation of the polarized brightness and degree of polarization. The PDFs are compared with data from the SECCHI-COR1 coronagraph onboard STEREO-A in order to detect systematic or random perturbations of the polarized brightness and the polarization angle beyond the unavoidable photon and detector hardware noise. This noise is estimated from two successive filter sequences taken in-flight during calm coronal conditions on 18 May 2008 and is expressed in the form of an intensity–variance relation. Two small deviations between the measured distributions and the predicted PDF for the polarization angle were found. The standard deviation of the polarization angle error decreases with increasing signal-to-noise ratio of the polarized brightness. For ratios larger than about 8 this decrease was found not as steep anymore as predicted which could hint to a small additional noise source. Next, we found a systematic constant deviation of the polarization angle by −1∘\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$-1^{\circ }$\end{document} for all signal-to-noise ratios of the polarized brightness. Besides these small discrepancies, our theoretically derived PDFs agree quite well with the distributions of measured brightnesses in test regions of the images. The PDFs we present here can equally be applied to similarly measured data from other coronagraphs and may help to quantify uncertainty limits of the derived polarization. They can be used for in-flight health checks of an instrument, are useful when separating unpolarized stray light from the polarized K-corona and when comparing the observed polarization data with results from model simulations.


Introduction
The white light observed by coronagraphs is composed of solar light scattered at coronal electrons (K-corona), interplanetary dust (F-corona) and emissions from coronal ions like Extended author information available on the last page of the article Hα, He D 3 , Fe X, Fe XIV, Ca XV, etc., in optical wavelengths (E-corona). These different light sources are difficult to separate unless the spectrum variation of the light is measured in parallel. Instead, many coronagraphs measure the linear polarization of the incoming signal. The polarized component should come only from the K-and E-corona, while the F-corona inside a field of view of about 5 solar radii can be assumed unpolarized (e.g., Koutchmy and Lamy, 1985;Levasseur-Regourd et al., 2001;Morgan and Habbal, 2007). Since the E-corona is restricted to few emission lines, the major part of the polarized broad-band white-light coronal emission is K-corona for most coronagraphs.
The polarization of the K-corona is caused by the directional anisotropy of the sunlight incident on coronal electrons (e.g., Minnaert, 1930;van de Hulst, 1950;Billings, 1966;Inhester, 2015). The degree of polarization depends on the distance from the solar surface and on the scattering angle to the observer. The resulting orientation of the observed polarization is parallel to the closest solar limb. Deviations from this orientation are expected only in rare cases, e.g., in the presence of relativistic electron beams in the corona which may occur in the vicinity of strong flares (Molodensky, 1973;Inhester, 2015) or if the Thomson-scattered light is contaminated by some E-coronal lines scattered by magnetically sensitive ions (e.g., House, Querfeld, and Rees, 1982;Heinzel et al., 2020).
Because the observed polarized brightness is much less sensitive to stray light, it has widely been used to estimate the electron density of the average corona and of the total mass of transient coronal mass ejections (CME) (e.g., Quémerais and Lamy, 2002;Feng et al., 2015). The degree of polarization observed from CMEs allows one to estimate their propagation angle off the plane-of-the-sky (e.g., Moran and Davila, 2004;Dere, Wang, and Howard, 2005;Vourlidas and Howard, 2006;Moran, Davila, and Thompson, 2010;DeForest, de Koning, and Elliot, 2017). It may also give bounds to the azimuthal widths of CMEs along the line-of-sight if observations from two vantage points are available (Bemporad and Pagano, 2015). The estimations of the CME mass, size and propagation direction have important implications for space weather forecasts. Another field where precise polarization data is required is to detect and interpret E-coronal line contributions in the white-light signal. Examples are the suppression of the polarization by resonant scattering of the Hα line (Mierla et al., 2011), or, as recently suggested by Heinzel et al. (2020), modifications of the degree of polarization and the polarization angle by a contamination from the helium D 3 -line at 587.7 nm. Both lines may be emitted from prominences and patches of cool plasma in the core of CMEs and allow inferences about the temperature and possibly of the magnetic field at the scattering site. Moreover, the decomposition into polarized and unpolarized contributions in the observed white light helps to estimate the internal, presumably unpolarized, stray light of the instrument. This stray light is an unwanted but unavoidable constituent of raw coronagraph images due the huge contrast between the coronal signal and direct sunlight of about 9 orders of magnitude.
In spite of this wide-spread use of the polarization data from coronagraphs, there are only relatively few investigations which attempt to quantify the error bounds of their results. One exception is Moran et al. (2006), who employ the observed deviation of the orientation of the polarization from the closest limb direction to test the performance of the instrument's polarization capability and derive corrections to compensate instrumental effects. In this paper, we attempt to investigate these bounds more closely. Besides possible systematic instrumental effects, the major source of the uncertainty is the image noise due to the quantization of the detected photons (shot noise) and from the detector hardware. We intend to quantitatively determine the influence of this noise on the uncertainty of the polarized brightness and the polarization angle from observations of coronagraphs such as SECCHI-COR1 onboard STEREO-A spacecraft.
The results we obtain can be useful in several ways. For one, they may help to detect temporal variations of the instrument performance during the mission, anomalies or systematic errors which are not caused by the expected image noise. The better we can quantitatively characterize the latter, the better are the chances to discover the former. In fact, this study was motivated by exploring ways to check in-flight the quality of the polarization data from the new coronagraph PROBA3/ASPIICS to be launched in 2023 (Lamy et al., 2010;Renotte et al., 2015;Galano et al., 2018). Another application of our results we see in the determination of error bounds in the various types of model calculations which use the observed degree of polarization or polarization angle which were mentioned above. E.g., before one of these quantities is interpreted in terms of a physical model parameter, e.g., CME propagation distance from the plane-of-the-sky and its line-of-sight extent from the measured degree of polarization, or the D 3 line contribution from the polarization angle, a reliable estimate of the uncertainties of the observed polarization data is highly desirable. Since Thompson scattering scales linearly with the coronal electron density, coronagraph observations are well suited for applying tomography to a series of coronagraph images from one or better more vantage points. An essential precondition for tomography is a good brightness scaling of the images (Frazin et al., 2012). Again, a quantitative error estimate is a prerequisite. If the observed polarization data and model simulations are compared by means of an inverse problem, the noise estimates could be decisive for the amount of regularization to be applied.
The estimation of measurement errors and the uncertainties of derived data products is a fundamental issue and has in a much wider context been discussed in Thomas et al. (2015). The attempt to integrate data and data errors in a common FITS-concept (Flexible Image Transport System) reveals how closely they are connected and how much a quantitative assessment of data errors is needed. This paper was inspired by this discussion and hopefully contributes to it within the realm of coronagraphy.
The plan of the paper is as follows. In Section 2, we derive the probability density functions (PDF) for the polarized brightness and the polarization angle as derived from typical coronagraph measurements such as the COR1 instruments of STEREO/SECCHI. We assume that the noise of the directly observed images is normally distributed and that this noise is the only source for the uncertainties. We here define this noise as the photon shot noise plus the hardware noise produced by the image sensor which includes dark current and read noise. It does not include non-random perturbations like stray light, ghosts or internal instrumental reflections. However, these enhance the image mean signal and thereby also add to its noise variance. While at least in principle we have a chance to subtract the mean of these artifacts, provided we can measure or model them perfectly, their contribution to the noise variance cannot be eliminated. The PDFs depend strongly on the level of this noise and they only approach the often invoked normal distribution when the signal-to-noise ratio is sufficiently small. In Section 3, the PDFs are compared to the noise distribution found in SECCHI-COR1 data observed on 18 May 2008. We first empirically estimate the noise variance for every pixel and express its dependence on the observed pixel intensity by an analytic fit which is used as a basis for the analysis that follows. In the remaining part of Section 3, we compare the statistical behavior of the observed polarized brightness and the polarization angle with the PDFs derived in Section 2. A major emphasis here is to find how well the observed mean values and standard deviations agree with the predictions from the PDFs. In a final section we discuss the results and limitations of our study.

PDF of the Polarized Brightness and the Polarization Angle
In the COR1 coronagraph, like in many other conventional coronagraphs, the primary observation consists of three successive intensity measurements through polarization filters with maximum transmission at 0 • , 60 • and 120 • . The primary measurement is therefore a vector I = (I 0 , I 60 , I 120 ) for each image pixel. The orientation of a linear polarization angle is unique only modulo 180 • , so that the 60 • polarization angle is also referred to by 240 • . For what follows it is, however, convenient to restrict the polarization angles to a range on which they are uniquely defined. We also omit the pixel indices here because they do not matter for our study except that we concentrate sometimes on limited image regions. A basic assumption is that each pixel is considered statistically independent.
To derive the polarization properties, the measurement vector is converted into a vector of Stokes' components S = (S I , S Q , S U ) by the demodulation matrix D, The first component of the Stokes vector is the total brightness S I . From the other two components, the polarized brightness S P and the polarization angle α are obtained by respectively. For consistency of our notation, we do not use the conventional nomenclature B and pB for the total and polarized brightness here. We do not expect any circular polarization in the light scattered from the coronal electrons. Circular polarization may, however, also be produced inside the optical path of a telescope (Sánchez Almeida and Martinez Pillet, 1992). These authors showed that an axisymmetric instrument design and small view angles off the optical axis prevent a noticeable cross-talk between polarizations. Most coronagraphs satisfy these conditions and we therefore omit the circular Stokes component here. Also, the above expression for the demodulation matrix D only holds in the ideal case that the three linear polarizers have unit maximum transmission and contrast and that their relative orientations are perfectly separated by 60 • . Deviations from this ideal can be fully compensated by adapting the elements of D. However, not knowing any better we assume that Equation 1 holds for COR1.
Due to the linearity of the relation between S and I , the properties of the noise in I are straightforwardly inherited by the noise in S. Especially, if the measurements in I are normally distributed random variables, the elements of the Stokes vector are normally distributed random variables, too. For a given pixel, we call I = (I 0 , I 60 , I 120 ) the true values we would have measured without noise. In this sense, we denote by the overbar the deterministic noiseless value also for other variables below. Since we assume that the intensities are normally distributed, I also represents the mean values of I but this may not necessarily hold for other variables. The respective variances of the measured (I 0 , I 60 , I 120 ) are σ 2 I = (σ 2 0 , σ 2 60 , σ 2 120 ). The noiseless components S are then related to I by the same relation as Equation 1 and they are again identical to the respective means of S. The variances of S are given by the sum of the variances of the terms which contribute to the respective Stokes component.
Hence the normal distribution of S is characterized by a mean and a covariance of The superscript denotes transposition and the brackets statistical averaging. We made use of the obvious fact that the noise of different components in I − I are independent. The elements of the covariance matrix W are explicitly The covariance matrix is symmetric and becomes diagonal only if the variances σ 2 i of I i are equal. If photon noise dominates, however, the respective σ 2 i tend to increase with signal intensity which yields different σ 2 i for different I i . As a consequence, we will in general also have different noise levels for S I , S Q and S U .
The probability distribution for the three-dimensional Stokes vector S given the mean I and its variances σ 2 I then is the trivariate Gaussian However, the quantities of Equation 2 we are interested in only depend on the two components S P = (S Q , S U ) of S. The probability density for these two components alone is obtained straightforwardly by marginalization, i.e., integration of Equation 6 over S I . For multivariate normal distributions the result is obtained systematically Figure 1 Contour lines of the bivariate PDF of Equation 7 in Stokes space (light blue). This distribution needs to be integrated along the dashed red circle to obtain the PDF P (S P ) of the polarized brightness S P or along the dashed dark blue radial line to obtain the PDF P (α) of the polarization angle α. We denote by α the true polarization angle of S P and u ± are the eigenvectors of V and .
where = V −1 and V is the two-dimensional covariance given by the Q, U -submatrix of W , The normalization in Equation 7 is To be consistent with Equation 2, we define a random state vector in 2D Stokes space and similarly, the deterministic noiseless Stokes vector, by S P =(S Q , S U ) = S P (cos 2α, sin 2α), S P =(S Q , S U ) = S P (cos 2α, sin 2α), so that (S P , 2α) and (S P , 2α) represent the conventional cylindrical coordinates of their respective Stokes vector. The cylindrical azimuth angle φ = atan(S U /S Q ) of S P with the Qaxis is therefore just half the polarization angle, analogously for S P . Note that α describes a unique polarization state only for values in [0, π]. Therefore every angle φ = 2α, as usual in [0, 2π], can be attributed to an independent polarization state. From Equations 7, 8 and 9 we can obtain the desired PDFs for S P and α by suitable integration in Stokes space as illustrated in Figure 1. We will see below that this integration is easily performed by making use of the cylindrical coordinates introduced in Equation 10. An example for a P 2D of Equation 7 is sketched in Figure 1 by means of contour lines.
Obviously, the results of the integrations for the desired PDFs depend strongly on the contours induced by . In order to scale the results in a meaningful way, we introduce a representative scalar variance parameter σ 2 P such that σ −2 P is the average of e e for unit vectors e of all directions in Stokes space. By this choice, the isotropic approximation → 1σ −2 P leaves the area inside a contour of a given contour level invariant. Then σ −2 P is also the mean of the eigenvalues λ ± of , i.e., where v ± = λ −1 ± are the respective eigenvalues of V . Equivalently, σ 2 P is the harmonic mean of v ± . Due to the symmetry of V , QU so that we can express our noise measure in terms of the covariance of S P − S P as Most often, the true intensity I and the variance σ 2 I of the noise added to it in a measurement are not independent. From our understanding of CCD and CMOS detectors, the detector output signal with an exposure time τ after subtraction of the detector bias and dark current is related to the incident irradiance E by I = gq Eτ + noise. Here, g is the electron conversion gain (units: DN/e − , note in the astronomical literature often the inverse is considered as "gain", see e.g., Astier et al., 2019) and q the quantum efficiency of the detector. The noise ideally has zero mean and a variance σ 2 = g 2 q Eτ + c = gI + c. We may apply this relation to the intensity of each polarization filter to Assuming that this relation holds, P 2D Equation 7 depends only on the true intensity vector I we would measure if its detection was noiseless. Applying this PTC relation systematically in Equations 3 and 4, the inverse two-dimensional covariance becomes Inserting S P = (S Q , S U ) = S P (cos 2α, sin 2α) with polarization angle α we easily see that the eigenvectors of become In other words, as the center S P of the contour ellipses in Figure 1 rotates with 2α, the ellipse orientation rotates oppositely with −α and the angle between u + and S P is β = 3α (see Figure 1). However, this behavior strictly holds only if Equation 12 is valid.

Polarized Brightness
As described above, the probability density for S P for a given polarized brightness S P = (S Q , S U ) and covariance −1 is obtained by the integration of the bivariate PDF of Equation 7 along a circle in Stokes space with radius S P . In the cylindrical coordinates Equation 10, the probability of the infinitesimal Stokes space element with 2α replaced by φ is P 2D (S Q , S U ) dS Q dS U = P 2D (S P (cos φ, sin φ)) S P dS P dφ. Integration over φ and omitting the integration measure dS P for P (S P ) yield P S P | I , σ 2 I = S P P 2D S P cos φ sin φ I , σ 2 I dφ (13) Figure 2 The PDF of Equation 13 of the normalized polarized brightness S P /σ P for some values of S P /σ P as given in the legend according to color. The colored ranges represent the numerically derived PDF values for all polarization angles α and degrees of polarization between p = 0 and 0.4. The solid curves give a specific example for p = 0.4 and a polarization angle so that the angle between the -eigenvector u + and S P is β = π/4 (see Figure 1). The dashed curves represent the respective isotropic approximation Equation 14 which is independent of α. For the PDF with the largest S P /σ P = 6, we have also indicated as green dashed curve the asymptotic normal Equation 16. For this value of S P /σ P , it is hardly distinguishable from Equation 14.
where the integration is over the full circle φ = 0 to 2π . This way we maintain the normalization of the probability, i.e., P (S P ) dS P = 1. Unfortunately, this integration in its general form has to be evaluated numerically. We choose a Gauss-Legendre integration rule which gives an error better than 10 −6 for 12 integration points if the integration range is restricted to the range of φ where the integrand is significant. A simple analytic form can be obtained when we approximate the covariance −1 by the isotropic matrix σ 2 P 1. This simplification holds if in the PCT relation Equation 12, either the non-photonic noise in parameter c dominates or if the polarization is low and all I i are about equal. In Appendix A we derive for this case Here, I 0 is the modified Bessel function of first kind (Abramowitz and Stegun, 1964). The ratios S P /σ P and S P /σ P can be considered as a signal-to-noise ratio of the respective polarized brightness.
In Figure 2, we show examples of the PDF Equation 13 and its approximation Equation 14. For the calculation of the former we assumed Equation 12 to hold with c = 0 and g and I tuned so that the desired ratio S P /σ P and degree of polarization p were obtained. The colored ranges in Figure 2 cover the values of the PDF for a degree of polarization p = 0 to 0.4 and all polarization angles α so that the angle β of the contour ellipse axis in Figure 1 rotates by 2π . The boundaries of the colored ranges are reached for angles β = 0, π/2 and π . The solid curves show a specific example with orientation angle β = π/4 for which the PDF is very close to the isotropic approximation Equation 14, represented by the respective dashed curves inside the colored ranges.
Since the PDF of Equation 14 with isotropic covariance is a decent approximation, we can use it to discuss the two limiting cases S P /σ P → 0 and → ∞. In the limit of vanishing Left: The polarized brightness argument S P max /σ P of the PDF maximum (black and gray) and mean value S P /σ P of the PDF of Equation 13 (red and rose), both as function of true polarized brightness S P /σ P . The colored ranges were obtained from Equation 13 with the same parameter ranges for p and β as in Figure 2. The black and red dashed curves again represent the respective quantities for the isotropic approximation Equation 14. The green dashed curve gives the maximum likelihood polarized brightness S P /σ P on the abcissa for a measured value S P max /σ P on the ordinate. The solid diagonal is drawn to show the small but finite offset of the above curves for larger arguments. Right: Standard deviation (black and gray) and rms distance to the true S P /σ P (red and rose) for the PDF of Equation 13 as function of the true polarized brightness S P /σ P . The colored ranges were obtained again from Equation 13 with the same parameter ranges for p and β as in Figure 2 and in the left panel. The dashed curves again represent the respective quantities for the isotropic approximation Equation 14.
polarized brightness, the PDF becomes which complies with the χ 2 2 -distribution for two degrees of freedom for the two random variables S Q and S U if each has zero mean.
For large values of S P /σ P , we find a normal distribution For details of the derivation, see again Appendix A. The PDFs of Equations 13 and 14 have a heavily biased mean value and also the argument of the PDF maximum differs considerably from the true polarized brightness. This has consequences for possible observations. For a single measurement, e.g., from a single image pixel, the most probable value of S P we will measure is the value S P max = argmax P (S P ) where the PDF has its maximum. Even if the true S P is zero, we most probably will measure a finite value for S P around σ P . If we have many independent, statistically equivalent measurements of S P , we may use them to reconstruct the PDF itself or average them to obtain approximately the expectation mean value S P = ∞ 0 S P P (S P )d S P . If, e.g., again S P = 0, a mean of even S P = √ π/2 σ P will result. For large S P /σ P , the bias of the PDF maximum argument and of the PDF mean value shrink less rapidly than one might expect from the relatively fast evolution of the PDF shape towards a normal distribution. In Figure 3, we show in the left diagram the most probably measured polarized brightness S P max (black and gray) and expection mean value S P (red and rose) from Equation 13 (colored ranges) and the isotropic approximation Equation 14 (dashed curve). Both values approach the true S P (represented by the exact diagonal) only very slowly. Even for medium signal-to-noise ratios of the polarized brightness, we tend to overestimate the true polarized brightness unless we correct for the bias.
We propose a correction of the bias for the above two estimators of S P /σ P : i) We have a single measurement of the polarized brightness, S P 1 /σ P and we assume that this measurement is the most likely one we could have made, i.e., it maximizes the PDF. This implies that S P 1 /σ P > 1. ii) We have a number of equivalent measurements and determined their mean to S P ∞ . Here we assume that this sample mean estimate agrees with the exact mean S P P (S P , /σ P )dS P . This implies that S P ∞ /σ P > √ π/2.
Replacing either S P 1 /σ P or S P ∞ /σ P by x, an improved value for S P /σ P can be found from Another way to deal with the bias is the maximum likelihood estimate of S P . This estimate seeks the most likely value of S P for a given measurement S P , i.e., argmax S P (I ) P (S P | I , σ 2 I ). Applied to the PDF of Equation 14 with isotropic covariance, this estimate is shown by the green dashed curve in Figure 3. For consistency of the plot axes, the measured S P is taken from the ordinate, the resulting maximum likelihood S P is assigned to the abcissa.
The left diagram of Figure 3 shows some second moments of Equations 13 and 14. The standard deviation, (S P − S P ) 2 1/2 /σ P is drawn in black and gray, the rms distance (S P −S P ) 2 1/2 /σ P from the true S P in red and rose. For large S P /σ P , the standard deviation is about σ P which proves that Equation 11 is a useful measure of the variance. For S P /σ P less than about 2, The standard deviation shrinks because S P is restricted to positive values. This is also the reason why the measured polarized brightnesses deviates from the true value by a considerable fraction of the measurement standard deviation.

Polarization Angle
To obtain the probability density for the polarization angle α for given polarized brightness S P = (S Q , S U ) and covariance −1 we need to integrate the bivariate PDF of Equation 7 along a radial line from the origin of Stokes space in direction e = (cos 2α, sin 2α). We again use the cylindrical coordinates Equation 10 to transform the infinitesimal Stokes space element P 2D (S Q , S U ) dS Q dS U = P 2D (S P (cos 2α, sin 2α)) S P dS P 2dα. Integration over S P and omitting the integration measure dα for P (α) yields Figure 4 The PDF of Equation 19 for the deviation of the polarization angle α from its true value α for five values of the parameter S P /σ P as given in the legend according to color. The parameters for the colored ranges are the same as in Figure 2. The solid curves show a specific example for p = 0.4 and a polarization angle so that the angle between the -eigenvector u + and S P is β = π/4 (see Figure 1). The dashed curves again represent the isotropic approximation Equation 14.
The integral can be evaluated by elementary means. To shorten the notation, we abbreviate the direction of integration in Equation 18 by (cos 2α, sin 2α) = e(α) and similarly (cos 2α, sin 2α) = e for the direction of the true Stokes vector so that S P = S P e and S P = S P e. With these definitions, the result can be written as (see Appendix B) where S 0 = S P e e e e , We note that by construction, Equation 19 is periodic on [α − π/2, α + π/2]. If we extend the range of α to 2π , Equation 19 for S P = 0 becomes the Bingham distribution on a circle (Bingham, 1974). Figure 4 in the same manner as for Figure 2. Note that the PDF is not necessarily symmetric as the special case demonstrates (see solid curves in Figure 4) where the parameters are such that the angle between the -eigenvector u + and S P is β = π/4. The PDF is symmetric to α only for orientation angles β = 0, π/2 and π . The PDF then coincides with one of the boundaries of the colored ranges in Figure 4.

Examples of Equation 19 are shown in
Again, the PDF considerably simplifies for some limiting cases. From Figure 4 clearly Equation 19 approaches a normal distribution when S P /σ P is large. In this case erf→ sign(e e) and we can set K = 2π/e e S 0 where e e > 0. So for large S P /σ P and e e > 0, For negative e e, we can in this limit safely set P (α | I , σ 2 I ) → 0. It can be shown that the exponent in Equation 20 minimizes at e = e, i.e. at α = α, but depending on the orientation of , also Equation 20 may not be completely symmetric about α.
Perfect symmetry is obtained only when we approximate the covariance −1 by the isotropic matrix σ 2 P 1. In this case, e e → σ −2 P , e e → σ −2 P cos 2(α − α) and the PDF of Equation 19 becomes which for large S P /σ P and |α − α| < π/4 further simplifies to The isotropic approximation Equation 21 is again shown as dashed curves inside the colored ranges in Figure 4. Of particular interest here is the variance of the PDF of Equation 19 because it gives us the statistical error we have to expect for the determination of the polarization angle. From the large S P /σ P -limit in Equation 20 and its isotropic approximation, the variance approaches The variance has units of rad 2 and the factor 1/4 reflects the ratio 1/2 between the azimuth angle in Stokes space and the polarization angle. For small and medium values of S P /σ P , the expectation values of mean and variance have to be calculated from suitable moment integrals over the PDF of Equation 19 which take account of the periodicity of the random variable (e.g. Mardia and Jupp, 2000). We here adopt the concept of the circular variance for the first and second moment. These moments are replaced the variance function The expectation mean value α for the polarization angle then minimizes V (ψ) and the corresponding minimum value V ( α ) is the circular variance of the PDF, i.e., In Appendix C, we propose a simple procedure by which Equation 25 can be solved.
In Figure 5 we show these moments in black for the same range of parameters used before. We compare α with the most probably measured polarization angle α max = argmax P (α) in the left diagram and in the right diagram the standard deviation σ α with the upper and lower quartile distances from α = α drawn in red. The deviation of α and in particular of α max from α is a result of the fact that the PDF for most angles β is not symmetric (see specific example in Figure 4). It is largest for angles β near π/4 and 3π/4 (see Figure 1). However, both are much smaller than the standard deviation and are probably only relevant for high precision determinations of the polarization angle. From the right diagram we see that σ 2 α,∞ in Equation 23 approximates the variance quite well for values of S P /σ P as low as ≈ 1.

Observations
The data we use to test the PDFs derived above were obtained from the COR1 instrument on board STEREO-A, one of the two Solar TErrestial RElations Observatory spacecraft (STEREO, see Kaiser et al., 2008). COR1 is the innermost coronagraph of the Sun Earth Connection Coronal and Heliospheric Investigation (SECCHI) instrument suite (Howard et al., 2008). It observes the corona through three linear polarizers at 0 • , 60 • , and 120 • in a 22.5 nm wide white-light waveband centered at the Hα line at 656 nm (Thompson and Reginald, 2008). The telescope has a field of view from 1.4 to 4 R but here, we will use the data only between 1.5 to 3 R . The image data was reduced onboard the spacecraft to 1024×1024 pixels resulting in a pixel resolution equivalent to 7.5".
The exposure time was 1.7 s for each polarized image, one full polarization sequence was taken within 18 s. The polarized image sequences were repeated after 10 min. In this work we investigate data from a polarization scan taken at 18 May 2008 at t 1 = 1405 UT and t 2 = 1415 UT. Within these 10 minutes, the Sun did not show any activity so that the measured intensities in equally polarized images varied very little.
For each individual image, the raw data was first cleaned using the SolarSoft secchi_prep routine. By default, secchi_prep applies the following steps: subtraction of the CCD bias, applying flat-field, vignetting and exposure time correction and an adjustment of the image orientation so that the y-axis points to solar north. Here, a roll angle of 7.5 • had to be applied. The standard processing steps comprise two more actions we did not readily apply. We did not multiply the image intensities with the calibration factor which converts the digital image values to fractions of mean solar brightness (MSB) but we rather prefer to continue in units DN s −1 for the components of the measured and preprocessed intensity vector. Another processing step is the subtraction of a monthly minimum background to clean the images from stray light, some internal reflections and other instrument artifacts . Here, we keep both, the image before background subtraction with intensities J = (J 0 , J 60 , J 120 ) and after subtraction with intensities I = (I 0 , I 60 , I 120 ). The background J − I is the same for the images at 1405 and 1415 UT but it differs for different polarizations.
In our subsequent analysis we only consider the restricted field of view from 1.50 to 2.94 R around the occulter center because for the quiescent time of the observations, the intensities at larger distances were very low. As a result, we have four polarization image sets J (t) and I (t) with and without background, respectively, at two observation times t 1 = 1405 UT and t 2 = 1415 UT. In Figure 6, we show for the limited field-of-view considered here all three components of I and, as an example, one polarization of J for one observation time. Note the strongly enhanced intensity level in J as compared to I .

Intensity-Variance Relation
So far we have treated the parameters of the PDFs, I and σ 2 I , as if they were independent. For most detectors, they are related, e.g., by the ideal relation Equation 12 with constants g and c from a pre-flight calibration.
Here we chose another approach. We have to keep in mind that the available COR1 data J and I have already been processed to some extent. For this reason we will work with an empirical intensity-noise relation between J and its variance. Since the detector cannot distinguish between light from the corona or from stray light, a unique relation between intensity and noise can only be expected to exist for the intensity J which still includes the background made up of stray light, internal instrument reflections and ghosts. The background correction between J and I is determined from data of a whole month around the observation and can be assumed to contain much less statistical noise than the individual image. So J and I can well be assumed to have the same variance. Moreover, since the CCD detector should not be sensitive to different polarizations, we expect that the noise variance does not depend on the polarization of incident intensity.
To find the intensity-variance relation for J , we make use of the fact that the coronal conditions changed very little during the 10 min between the two polarization sequence observed so that the difference images are dominated by the same level of image noise. With this prerequisite in mind, we then use for each pixel the difference between images at t 1 and t 2 of equal polarizer orientation as a noise sample and the respective average as the corresponding intensity estimate. For each pixel and each polarization i = 0, 60, 120 we form The factor 1/ √ 2 in J i accounts for the fact that the intensity difference has twice the noise variance of a single image. An example for one such difference image is shown in Figure 7.
For each pixel, the pair of numbers in Equation 26 only represents a random sample of the desired intensity-variance relation. The relation can be estimated by suitable averaging Figure 6 Images of the intensities I 0 , I 60 and I 120 without background (top row and bottom left) and one of the measured intensities J 120 which still includes the background (bottom right). All images were observed at 1415 UT. Note the different color code for the intensities I and J . The inner ring represents the position of the solar disk, the regions 1 to 5 in image I 120 are marked here for later reference.
where a trade-off has to be sought between a sufficiently large number of pixels and statistical homogeneity for each average. The compromise between these two requirements was found experimentally. We segment the images in small tiles with a size of 3 pixels in radial direction and 2.5 • in azimuth. These 9600 tiles have between 25 and 50 pixels each. Within each tile, we form tile-averages J i = tile-mean( J i ) and In Figure 8, we plot in color-code the relative number of occurrence of pairs (J i , J 2 i ) summed over all polarizations i = 0, 60 and 120 in discretized 2D classes of J and J 2 . The size of each color tile in Figure 8 represents the respective class range.
This fit is overplotted as smooth curve on top of the distribution in Figure 8.
Equation 27 does not have the expected PTC-like form Equation 12. Instead, we had to add a term proportional to J 3 in order to account for the steep increase of the observed variance in Figure 8 for logJ > about 3.5. A quadratic term would not have been sufficient here.
Higher order terms in a PTC-relation are known and some are well understood (see, e.g., Astier et al., 2019). The most common reason is fixed-pattern-noise which introduces an additional J 2 -dependence when the PTC curve is determined from the variance with respect to an intensity average over the entire field-of-view (fixed-pattern-noise, see, e.g., Janesick, 2007). In Equation 27, only intensity differences from the same pixel enter the variance and in addition, the image was flat-fielded before which should have drastically reduced the effect of the fixed-pattern-noise in Equation 27. Some of the tiles we average over, especially in the bright areas near the occulter edge, may be intersected by the bright internal reflections visible in J (see lower right of Figure 6). This may cause some enhancement of the variance, but we checked that this effect is small by reproducing Equation 27 with different tile-sizes. We rather assume that the steep increase in Equation 27 is due to the vignetting component in the flat-field correction (Thompson et al., 2011) which causes an enhancement of the detector counts proportional to the vignetting correction factor ≥ 1. The local variance is consequently enhanced by this factor squared. This correction especially affects the region near the occulter edge where intensities and vignetting corrections are large. A solid explanation of Equation 27 is beyond the scope of this paper. Instead, we just use it as an empirically-determined, effective PTC-like relation.
As a check of Equation 27, we have binned the J i values from Equation 26 of all polarizations into logarithmically equidistant classes and determined the mean values of ( J i ) 2 and their standard deviations for each class. These tile-independent variance estimates are overplotted in Figure 8 as a step-wise curve with vertical bars indicating the uncertainty of the respective class mean estimate. The horizontal width of the steps indicate the class width of log J i . The mean variances obtained follow Equation 27 quite well. A small but significant deviation occurs between log J i = 3.6 to 3.9 where the distribution of ( J i ) 2 within a class is skewed so that mean and maximum deviate by up to a factor 1.4. This makes at most a factor 1.2 in the standard deviation estimate of the image noise.
In the following we will use the background subtracted pixel intensities I to calculate samples of the polarized brightness and the polarization angle using Equations 1 and 2 as in the conventional polarization analysis. The expected variance of the pixel intensities in Equation 4 is approximated by Equation 27 substituting for each pixel and polarization the intensity J i for the argument J . Since the background correction between J and I is determined from data of a whole month, it can be assumed to contain much less statistical noise than the individual image. For this reason, it is safe to use σ 2 J ( J i ) also as noise variance for I i for the same pixel and the same polarization i. The variance vector in Equation 4 for a given pixel is accordingly given by where J i , i = 0, 60, 120 are the three components of Equation 26 at the respective pixel.

Polarized Brightness
We test the polarized brightness PDF in Equations 13 and 14 with the above observations in three small regions 1 to 3 which we select so that they represent different levels of the normalized polarized brightness S P /σ P . In region 1, the observed signal is almost pure noise, with region 2 we have tried to cover a faint plume, better visible in the image of the polarization angle in the next chapter. Region 3 is representative of a large signal strength and is located near the axis of a streamer. The regions are displayed in the upper right panel of Figure 6 and also in Figure 9 which gives an overview of the observed normalized polarized brightness for the entire field-of-view. The values of S P were obtained for every pixel individually from Equations 1 and 2. Similarly, the variance parameter σ 2 P was calculated for every pixel from Equation 11 where Figure 9 Normalized polarized brightness S P /σ P . The color code is logarithmic. Inside the three regions, the statistical distribution of the polarized brightness values is analyzed in more detail.
σ 2 QQ , σ 2 UU and σ 2 QU were derived from Equations 28 and 5. The fringe-like structures in Figure 6 near the occulter edge are due to the instrumental reflections and stray light patterns in the observed images J and hence in σ 2 P . The requirement for a statistically homogeneous sample from the regions had to be compromised with the size of the regions and a sufficient number N of pixels to reduce the statistical uncertainty. The pixel numbers are N = 1909 (region 1), 589 (region 2) and 445 (region 3). We reduce the small inhomogeneity of the data in each region by fitting a 2nddegree 2D polynomial to each component of I and subtracting the variable part of the respective fit relative to the region center. This has practically no effect except for region 3 which sits on the axis of a streamer. It includes a non-negligible radial intensity gradient which is compensated by this detrending procedure.
To compare the observed distribution of the normalized polarized brightness S P /σ P with the polarized brightness PDF from Equation 13, we use as PDF parameters I its regional mean and for σ 2 I the variances from Equation 28 with the regional means of J as argument. Since I and J are normally distributed, their mean values are unbiased and their variances are only 1/N times the variance of a single pixel intensity as in Equation 27. Figure 10 for region 1 represents the case where the K-corona is negligible and the observed intensity is unpolarized. S Q and S U have practically zero mean while the total brightness S I has a finite mean which is probably the result of residual stray light. The fact that the mean values of S Q and S U almost vanish gives us some confidence that the different backgrounds subtracted to obtain I 0 , I 60 , I 120 is consistent at least in region 1. The distributions for S Q , S U and S I comply almost perfectly with normal distributions with the parameters derived form the regional means of I and J as explained above. These parameters are always given in the legend of the respective figure.
The same holds for the observed distribution of S P /σ P in the right diagram which obeys closely the almost identical PDFs in Equations 13 and 14. Note that the best estimate of S P /σ P in region 1 is 0.02. However, the most probably observed pixel value is near ≈ 1. This explains the dominant blue to light blue color in the polar sectors in Figure 9. It corresponds to a value of log 10 S P /σ P near 0 even though the probably true polarized brightness there is much smaller and practically below detectability.
With region 2 we have tried to cover a faint plume above the southern polar region. Therefore the region was placed close to the occulter edge, it is only 5 • or 20 pixels wide and has 30 pixels in radial direction. Since the region 2 is almost exactly south of the Sun  center, we expect for the polarized components in S Q and S U a polarization angle close to 90 • with respect to the vertical. A signal polarized at this angle has a 3D Stokes vector along (1, −1, 0). The negative shift of S Q relative to S U is a clear indication of a finite polarized brightness. The relatively larger observed S I accounts for a low degree of polarization. The best estimate of the polarized brightness is S P = 0.82 σ P = 10.3 DN/s even though the most probably observed pixel value is about twice as much (right diagram of Figure 11). The observed distribution has much more scatter than region 1 which we attribute to the reduced number of pixels to about one third and also to the influence of the stray light and reflection features in J so close to the occulter edge.  Region 3 is an example for a strong K-corona signal. The observed and expected distributions are shown in Figure 12. Since the region has an azimuth from 111 • to 116 • , the polarization angle should be centered at 23.5 • and the polarized component of the signal should have a 3D Stokes vector in direction (1, 0.68, 0.73). The light observed from this region also has a decent unpolarized component so that S I is larger than (S 2 Q + S 2 U ) 1/2 = 153 DN/s. Remarkable is the fact that different from this expected Stokes vector direction we observe S Q > S U which causes a deviation of the observed polarization angle from 23.5 • by about −1 • (this issue is discussed further in the next section). The distribution of the observed polarized brightness S P /σ P again agrees with the theoretical expectation, except that our variance parameter σ 2 P slightly overestimates the uncertainty here. This leads to a little wider theoretical distribution than observed. The best estimate of the polarized brightness in this region is S P = 11.2 σ P = 153 DN/s. For this large value, the PDF in Equation 13 comes close to a normal distribution and averaging S P /σ P over all pixels reproduces in this case the best estimate S P /σ P for this region quite well.
For all three comparisons of the measured S P /σ P distributions with the PDF we have performed χ 2 goodness-of-fit tests (e.g., Fisher, 1934;Barlow, 1989). We used the isotropic approximation Equation 14 of the PDF here (short dashed curves in Figure 12) because it has only a single free parameter, S P /σ P and the differences to Equation 13 are only small. The χ 2 test is widely applied in different science fields but the procedures and definition of indices slightly differ. In appendix D we therefore have briefly summarized how we performed the test and define the comparison indices, effect size and p-value, we use here. These indices are easily calculated but are much less easily interpreted (see e.g., Halsey et al., 2015).
A common rule is that effect sizes above about 0.1 are a hint that the measured distribution has been "affected" so that it differs sufficiently from the theoretical PDF. Larger effect sizes corresponds to larger discrepancies. The p-value is a measure for the probability that the observed sample was drawn from numbers distributed like the comparison PDF. This hypothesis is usually considered acceptable if the p-value is above about 0.05 to 0.1.
The outcomes of our tests are listed in Table 1. Except for region 1, the PDFs do not appear to comply with the observed distributions according to the common interpretation of the test. A difference between our test to their common use is that we did not fit the PDF parameter S P /σ P but, as explained above, S P and σ P were determined form the regional means of J and I . In addition, the inhomogeneities in J especially in region 2 probably devalidate the homoscedasticity requirement (equal variance of measured data) for the χ 2 test. Given these problems with the measured data, we still consider the observed distribution in acceptable agreement with the theoretical PDFs, except that its variance parameter σ 2 P seems to be somewhat underestimated in region 2 and overestimated in region 3.

Polarization Angle
To compare the PDF in Equation 19 of the polarization angle with the observed data, we use regions 4 and 5 (see Figure 6), which more or less cover the two streamer regions above the east and west limb. For every pixel in these regions we calculate the discrepancy α − α from Equations 1 and 2 and also the pixel value of normalized polarized brightness S P /σ P . We then bin the angle errors in classes of different polarized brightness levels and determine their distribution inside each class. The choice of the regions 4 and 5 insures that we have a wide range of polarized brightness levels.
Since the regions cover a considerable azimuthal range, we use the closest local limb orientation as estimate of α for each pixel, Here, φ is the pixel's azimuth angle with respect to the Sun center which is displaced from the occulter center. The azimuth is defined clockwise with φ = 0 towards north and a branch cut in the south. With this definition, φ here varies in the range [−π, π] and the polarization angle here is confined to the range [−π/2, π/2]. The polarization angles α and α are shown in the left panel of Figure 13, the theoretically expected standard deviation σ α for α − α calculated for each pixel using Equation 25 in the right panel. The smallest standard deviation σ α of ≈ 1.3 • is predicted just above the occulter in the bright streamer regions where S P /σ P rises up to 22. Since the artificial reflections in the background part of J enter into σ 2 QQ , σ 2 UU and σ 2 QU their imprint can also be found in σ P in Equation 11 and therefore also in the standard deviation σ α . They are particularly visible near the occulter edge.
To estimate the empirical distribution of the error in the polarization angle, we can compare α pixel-by-pixel with its theoretical value α of Equation 29. This method has successfully been applied by Moran et al. (2006) to data from the SOHO/LASCO-C2 coronagraph to correct systematic errors of the coronagraph polarization data. In the upper panel of Figure 14 we plot the polarization angle for every pixel with S P < 0.71 σ P versus the azimuth angle φ of the pixel position relative to the Sun center. The ideal result should be a twoteeth sawtooth curve represented by the two oblique solid lines. It is obvious from the color code, that pixels with large values of S P /σ P (red dots) come closer to their nominal values than the pixels with low S P /σ P values (blue). In the panel below, we show the same data Figure 13 Left: Polarization angle α for each pixel. The angles for pixels with S P < 0.71σ P were omitted. The color code inside the occulter indicates the angle α of the closest limb orientation for reference. Right: Standard deviation σ P derived from Equation 25 and the measured intensity J for each pixel.

Figure 14
Top: Scatterplot of the pixel-wise polarization angle vs pixel azimuth with respect to the Sun center. The color code represents the magnitude of parameter S P /σ P . Pixels with S P < 0.71σ P were omitted. The solid vertical lines separate the four azimuth quadrants, the oblique lines indicate the theoretically expected α. The dashed vertical lines represent the azimuths of the region boundaries 4 and 5 at 1.5 R . Bottom: The above data with the same color code in spatial cylindrical coordinates centered at the sun center. Due to the transformation from Cartesian image to cylindrical coordinates, most of the spatially isolated pixels with S P < 0.71 σ P were lost due to interpolation. The transformed region boundaries defined with respect to the occulter center are drawn as dashed curves.
with the same color code but with the polarization angle replaced by the distance from the Sun center to help identify from which image location the polarization angles were obtained.

Figure 15
Logarithm of the distribution of the angular deviation α − α of observed polarization angles in Figure 13 in regions 4 and 5 for different classes of the normalized polarized brightness S P /σ P . The color refers to the respective brightness class as given in the legend. The dashed lines are the respective theoretical PDFs from Equation 21 with parameter values S P /σ P one third from the bottom and two thirds from the top boundary of the respective brightness interval. For S P /σ P > 12, we used S P /σ P = 15 as parameter in Equation 21.
In the streamer regions, even the low S P /σ P pixels near the outer edge see some polarized light so that the relation S Q and S U is not entirely random. In Figure 13 we can therefore determine meaningful polarization angles almost to the outer edges of these regions. In Figure 14, the low S P /σ P -valued pixels (blue) in these regions remain still somewhat concentrated near their nominal positions.
In the polar azimuth sections, on the other hand, I is practically unpolarized. Yet there are some few pixels with S P just above the 0.71 σ P threshold, but they are mostly scattered over all polarization angles. We assume that these pixels arise from statistics because for a given pixel the value of S P /σ P has a considerable uncertainty and there is a finite probability that some pixels with S P exceeding 0.71σ P may occasionally be observed even in the polar regions with very little K-corona and practically unpolarized light. In these polar sections, the relation between S Q and S U should be completely random. However, the polarization angle distribution shows an enhancement roughly by a factor 2 near α = 0 • and ±60 • . These pixels can be traced back to intensity vectors I for which one element happens to be much larger than the other two. The enhanced probability of these pixels is again a statistical artifact and we could qualitatively reproduce for unpolarized I , discarding S P /σ P -values below a fixed threshold and σ 2 P rising with I as for photon noise. We have classified the data from every pixel in regions 4 and 5 according to their value S P /σ P and constructed for each class the distribution of the error α − α of the polarization angle. The result is displayed in Figure 15. We also compare the measured distributions with the theoretical PDF in Equation 21 indicated by the dashed curves in Figure 15. The important polarization brightness parameter S P /σ P in Equation 19 was in each class set to s = 2/3 s low + 1/3 s up where s low and s up are the lower and upper bounds, respectively, because the number of pixels decreases with S P /σ P . This weighting roughly represents the gradient of the S P /σ P occurrence distribution at medium values S P /σ P ≈ 3. The S P /σ Prange of the respective classes is given in the legend of Figure 15. We used the PDF with isotropic covariance here, because the pixels in regions 4 and 5 span a wide range of polarization angles. So some of the discrepancy between observed and predicted polarization angle errors might be due to the fact that the polarization angles are not equally averaged over in regions 4 and 5. The predictions match again well except that the observed distributions are slightly shifted in α − α and their widths appear slightly wider than observed. We again applied χ 2 goodness-of-fit tests to compare the observed and predicted distributions. The results are listed in Table 2. All p-values are way below 0.01, but the effect sizes at least for lower S P /σ P -values promise a closer agreement between the observed distribution and the theoretical PDF. Again, the major reason for the discrepancy is apparent from Figure 15: the variance parameter σ P derived from J slightly underestimates the width and the observed distributions are slightly shifted from the center α − α = 0. In view of the uncertainty of the intensity-variance relation in Equation 27 which directly affects the normalization of the observed polarized brightness S P /σ P , we still find that the agreement between the observed angular distributions and the predicted PDF in Equation 19 is acceptable.
For a more quantitative comparison of the mean polarization angle error and its standard deviation, we have binned the brightness in regions 4 and 5 more closely in classes from S P /σ P = 0.71 to 15.21 in steps of 0.5. We determined the empirical distribution of the polarization angle for each class and calculated its circular mean and standard deviation as defined in Equation 25. The result is presented in Figure 16 by the step-like curves. The error bars correspond to the standard deviation of the respective measured estimate assuming that the polarization angle errors are normally distributed.
In the upper diagram we show the measured circular standard deviation of the polarization angle error derived according to Equation 25. It is slightly larger than the theoretical values derived from Equation 21 or from its asymptotic approximation Equation 23 (overplotted as black and red dashed curves, respectively). As to be expected, they all decline with the mean of S P /σ P for the respective class. For S P /σ P > about 8, however, the decrease of the measured σ α does not continue quite as steep as expected which leads to an enhanced separation between measured and expected standard deviation.
In the lower diagram of Figure 16, we present the deviation of the mean polarization angle for each class from the local true angle α. We find a systematic angle error of about −1 • at all levels of the polarized brightness. The error bars here represent the standard deviation of the measured mean angle error and is about 1/ √ n times the standard deviation from the diagram above with n the number of pixels contributing to the mean. The small error bars and the fact that the mean angle error persists over the whole range of polarized brightnesses indicate that the −1 • deviation is clearly significant. Note that the measured standard deviation in the top diagram is the deviation from the measured mean shown in the bottom diagram. It therefore does not include the −1 • shift from α and the increased distance between measured and theoretical σ α cannot be explained by this shift.

Summary and Discussion
We have derived PDFs of the polarized brightness S P and of the deviation of the polarization angle α − α from the expected orientation parallel to the closest limb. The PDFs are applicable to data from a coronagraph which measures the polarization state by three successive intensity measurements I = (I 0 , I 60 , I 120 ) through linear polarization filters mutually rotated by 60 • . The basic assumption is that the uncertainty of the polarized brightness and polarization angle is due to the noise of the directly observed image intensities alone. In general, the PDFs depend on the vector of the true polarized intensities I and their variances σ 2 I . Since only the polarized brightness matters here, an equivalent set of parameters are the true Stokes components S P = (S Q , S U ) and their covariance matrix −1 . If high precision is not required or the polarization is low, an isotropic approximation −1 → σ 2 P 1 still leads to useful representations of the PDFs. For the random variables S P /σ P or α − α they depend only on the signal-to-noise parameter S P /σ P . The variance parameter σ 2 P is defined in Equation 11 and approaches the true variance of S P for large S P /σ P . The asymptotic variance of α − α is σ 2 α,∞ = σ 2 P /4S 2 P (in rad 2 ). The PDFs allow one to estimate the statistical uncertainty of polarization data products from coronagraphs and, by comparison with the observed variance of the polarization data, can also help to detect error sources of the instrument which do not originate in the presumed noise from the image detection.
A basic result which is obvious but often yet unnoticed is that the polarized brightness estimate from the data of a single pixel or from an average over a small pixel area can be Figure 17 Distribution of the degree of polarization p = S P /S I for regions 1, 2 and 3. The vertical line marks the best estimate of the degree of polarization for the respective region obtained from averaging I and calculating p from these averages.
heavily biased towards larger values. With the PDF in Equation 13, this effect can be quantified and we derived corrections Equation 17 which reduce the bias of the PDF mean and the argmax PDF values to a large extent, although not the uncertainty of S P . The bias affects most strongly the polarized brightness when S P < 2σ P . Here we recommend to improve the signal-to-noise ratio by a local or temporal averaging of the primary measurement vector I rather than averaging the S P values, because I are normally distributed and their mean values are unbiased. The averaging should include so many pixels that the effective parameter σ P shrinks below about S P /2.
The problem is even worse for the degree of polarization p = S P /S I which depends on a further random variable, the total brightness S I . Since p is often used to estimate the distance of the coronal scattering site off the plane-of-the-sky, it is a key parameter in 3D coronal density reconstructions. To illustrate the problem, we show in Figure 17 the distributions of pixel values of p in regions 1 to 3 and compare it with the much more reliable value obtained by calculating the degree of polarization from the regional averages of I . For noise-free measurements there is the obvious constraint S 2 P = S 2 Q + S 2 U ≤ S 2 I and an analog constraint for the elements of I since it is connected to S by Equation 1. If noise is added, the elements of I and similarly of S become statistically independent and the above constraint cannot be guaranteed any more. As a consequence, the measured degree of polarization may exceed unity.
In this paper, we have compared the PDFs with distributions of the polarized brightness and the polarization angle determined from data of the STEREO-A/ SECCHI-COR1 instrument. The observations were made during 18 May 2008 when the corona showed quiescent streamers above the eastern and western equator.
For many detectors, the σ 2 I depend on the measured intensities I by an intensity-variance relation characteristic for the detector, often represented by the photon transfer curve (PTC). For the observations used here, this relationship is complicated by the fact that a relatively large background component has to be subtracted from the detected images to obtain the science data I . The intensity J including the more or less unpolarized background had to be used to establish the intensity-variance relation. Since J is only little polarized, the resulting variance vectors σ 2 I are well balanced and yield almost isotropic covariance matrices −1 with the off-diagonal element well below 10% of the smallest diagonal element even where the polarization of I is large.
We find not unexpectedly that the image noise estimated from the intensity-noise relation dominates the uncertainty of the polarization data for STEREO-A/SECCHI-COR1. The observed distributions of the polarized brightness and the polarization angle agree in general with the expected PDFs besides a slight underestimate of the predicted uncertainty in most cases. We consider the agreement as adequate even though the χ 2 goodness-of-fit test failed in all but one comparison. We assume that the data is too inhomogeneous to meet the restrictive requirements of the test and moreover, the PDF parameters were not optimized for the comparison.
We have found two noticeable discrepancies of the polarization angle error from the theoretical PDF. The expected decrease of the polarization angle error with normalized polarized brightness agrees with our calculation up to about S P /σ P = 10. Beyond this value the standard deviation still continues to decrease but less rapidly than expected. At S P /σ P = 15, the measured standard deviation goes down to about ±2.5 • , the value which should be reached in accordance with our noise estimate is ±1.8 • . This small discrepancy could be caused by an additional noise source not included in the intensity-variance relation Equation 27 and becomes only visible where the influence of our noise estimate becomes sufficiently small. Alternatively, this effect may in principal also have been caused by underestimating the noise in the intensity-variance relation Equation 27 for the higher pixel intensities. We do not see evidence for this in the fit curve of Equation 27 in Figure 8.
We also found a systematic deviation in the polarization angle of −1 • from the closest limb orientation. This value exceeds the uncertainty of the mean polarization angle estimate for most of the of S P /σ P classes we investigated. This offset can simply be remedied if we modify the ideal demodulation matrix D in Equation 1 by multiplication with an appropriate rotation matrix operating on I and equivalently by correcting the element-wise square of D in Equation 4. The reason for this small offset is very probably the mechanics of the polarizer wheel construction. The stepper motor positions rotating the polarizer are discretized in steps of 2.5 • with a repeatability of each polarizer position better than 30 (Howard et al., 2008). It seems that the −1 • offset we found yields the absolute angular reference of the 2.5 • position grid of the polarizer.
We consider the agreement between the expected PDFs and the observed distributions from different test regions to be quite satisfactory in general. However, our comparison has a weakness. The intensity-variance relation we use has some uncertainties as it has been determined from the observed and routinely processed data rather than from a laboratory calibration which yields a genuine PTC. This effective intensity-variance relation is the basis for estimating the variance parameter σ P by which the observed polarized brightness S P has to be normalized for the comparison with the PDFs. Ideally, the variance would have been measured from a series of homogeneously illuminated images. Many space-bourne coronagraphs allow one to take these images through frost glass filters. But space missions seldom afford to spend dozens of such images for this purpose. Our approach to measuring the variance from partially processed data is not ideal, but our results show that it may be adequate.

Appendix A: Derivation of the PDF of Equation 14 for the Polarized Brightness
With the isotropic approximation = σ −2 P 1 the integrand in Equations 7 and 13 simplifies considerably to P 2D S P | I , σ −2 P 1 = 1 2πσ 2 P exp − 1 2σ 2 P (S P − S P ) 2 .
Using the Stokes space angle φ = 2α instead of polarization angle and similarly φ = 2α, P iso (S P ) = S P 2πσ 2 P φ+π φ−π exp − S 2 P + S where we used the integral representation of the modified Bessel function I 0 of the first kind (Abramowitz and Stegun, 1964, 9.6.16). For S P → 0 we can set I 0 → 1 so that the χ 2 2 distribution for two degrees of freedom results For large S P , we find that P iso is negligible except for S P ≈ S P . For large arguments of I 0 we can then use the asymptotic approximation I 0 S P S P σ 2 P → σ 2 P 2π S P S P exp S P S P σ 2 P , which yields P iso (S P ) → 1 √ 2πσ P S P S P exp − (S P − S P ) 2 2σ 2 P ≈ 1 √ 2πσ P exp − (S P − S P max ) 2 2σ 2 P .
Here S P max has to be chosen so that the maximum of P (S P ) is at the correct position. This gives S P max = S P 2 1 + 1 + 2σ 2 P S 2 P ≈ S P 1 + σ 2 P 2S 2 P .

Appendix B: Derivation of the PDF of Equation 19 for the Polarization Angle
After splitting S P and S P into magnitude and unit direction as defined in the main text, we write the negative exponent of P 2D in Equation 7  Note that e e > 0 because is positive definite but e e in S 0 may become negative when α − α is large.

Appendix D: Details of the χ 2 Tests
For the calculation of the χ 2 variable we essentially follow Fisher (1934) and Barlow (1989). Given N pixels in the test region 1 to 3 with measured normalized polarized brightness values x i , i = 1, N, we bin the data in classes j = 1, M according to their xvalue. Alternatively, x i is the polarization angle error of the N pixels in region 4 and 5 with approximately the same S P /σ P . The random number of n i entries in each class with x-values ∈ [x j − dx/2, x j + dx/2], should then be comparable to the expected number μ j = N P (x j ) dx where dx is the bin width of the class and P (x) the expected PDF. If all x i are samples from the same distribution then the numbers n j have a multinomial distribution with mean and approximate variance of μ j and for large N χ 2 = j =1,M (n j − μ j ) 2 μ j approaches a χ 2 k -distribution for k degrees of freedom. Here, k = M − c − r where c are the number of parameters in P (x) adapted to fit the measured distribution. For the expected PDF we chose Equation 14 which has essentially one parameter, S P /σ P . It was not directly fitted to minimize the discrepancy between measurement and PDF but determined from the original pixel intensities I and J from which in the end also the x i were calculated. We therefore still set c = 1. Moreover, we ignored in the summation for χ 2 all classes j for which n j ≤ 5. They also do not count for the degrees of freedom and their number r must therefore also be subtracted to obtain k.
A key test index is the p-value (called χ 2 -probability in Barlow, 1989, p. 150) defined as the probability that the value of χ 2 is larger than calculated under the hypothesis that the x i are a random realization of PDF.
where χ 2 is the value summed up above. Small values for this probability indicate a bad fit. A dividing line is commonly drawn at a p-value of 0.05.
The p-value has been criticized because it seems to vary strongly with different samples from the same distribution (e.g., Halsey et al., 2015). As alternative, the effect size has been proposed as test index but here the definitions vary considerably. We here adopt the definition of (Cohen, 1988, p. 216 to 223). The square of this effect size is already mentioned in (Fisher, 1934, p. 82) as Lexis' measure of dispersion Q.
It is largely independent of the sample size N and also of M because dx ∝ 1/M. Values about 0.1 to 0.3 and larger are supposed to suggest a difference between the PDF and the measured distribution, i.e., some manipulation on the measured samples show an "effect" as its distribution deviates sufficiently from the standard. Reversely, small effect sizes are taken as evidence that the observed sample was drawn from a distribution described by PDF.