# Hybrid Evolutionary Optimization Approach for Oilfield Well Control
Optimization

Ajitabh Kumar
Visage Technology

###### Abstract

Oilfield production optimization is challenging due to subsurface model
complexity and associated non-linearity, large number of control parameters,
large number of production scenarios, and subsurface uncertainties.
Optimization involves time-consuming reservoir simulation studies to compare
different production scenarios and settings. This paper presents efficacy of
two hybrid evolutionary optimization approaches for well control optimization
of a waterflooding operation, and demonstrates their application using Olympus
benchmark. A simpler, weighted sum of cumulative fluid (WCF) is used as
objective function first, which is then replaced by net present value (NPV) of
discounted cash-flow for comparison. Two popular evolutionary optimization
algorithms, genetic algorithm (GA) and particle swarm optimization (PSO), are
first used in standalone mode to solve well control optimization problem.
Next, both GA and PSO methods are used with another popular optimization
algorithm, covariance matrix adaptation – evolution strategy (CMA-ES), in
hybrid mode. Hybrid optimization run is made by transferring the resulting
population from one algorithm to the next as its starting population for
further improvement. Approximately four thousand simulation runs are needed
for standalone GA and PSO methods to converge, while six thousand runs are
needed in case of two hybrid optimization modes (GA-CMA-ES and PSO-CMA-ES). To
reduce turn-around time, commercial cloud computing is used and simulation
workload is distributed using parallel programming. GA and PSO algorithms have
a good balance between exploratory and exploitative properties, thus are able
identify regions of interest. CMA-ES algorithm is able to further refine the
solution using its excellent exploitative properties. Thus, GA or PSO with
CMA-ES in hybrid mode yields better optimization result as compared to
standalone GA or PSO algorithms.

##  1 Introduction

Oilfield planning and operations related decisions such as well placement and
production settings have direct impact on cost and revenue for any given
project. Hence, field development planning is done to maximize fluid recovery
or net present value of cash flow, and is often accomplished using reservoir
simulation studies. Simulation results are used to compare the performance of
various well configuration and controls. Optimization algorithms are used in a
workflow to search for well settings which maximize return on investment [1].
Search involves running thousands of reservoir simulation run for a range of
possible solutions and identifying a few solutions of interest for further
deliberations.

Optimization algorithms used for such automated workflow are mainly divided in
two categories: point based and population based. In this work, population
based evolutionary optimization approach is used for oilfield well control
optimization problem as these algorithms are relatively simple to program and
more effective for a large scale, highly non-linear, and multi-modal problems.
Moreover, these algorithms do not require gradient calculations and can also
be used for discrete problems. Two popular evolutionary optimization
algorithms, genetic algorithm (GA) and particle swarm algorithm (PSO), are
used in standalone mode for solving oil well control optimization problem.
Next two optimization methods, GA and PSO, are used with another popular
optimization method, covariance matrix adaptation – evolution strategy (CMA-
ES), in hybrid mode. CMA-ES is also a gradient-free evolutionary optimization
algorithm and has been successfully used in different problem settings [6].

Computational requirement of such oilfield well design and control study is
quite high, and hence commercial cloud computing is used to handle workload.
Parallel programming is used to distribute computational workload among nodes
of a large cluster. Two different computational setup are used, one in which a
single cluster with a large number of cores is used, and another in which
multiple mid-sized clusters are used with a common file server.

The main contributions of this work are:

  * •

Hybrid evolutionary algorithms, GA-CMA-ES and PSO-CMA-ES, are proposed for
oilfield well control optimization problem.

  * •

Hybrid algorithms are tested using well control optimization study of large
Olympus benchmark case.

  * •

Robust optimization is done using an ensemble of 10 simulation model
realizations, and results are compared with that from a single realization
case.

  * •

Parallel distribution of computational workload on cloud is done in two ways,
one between nodes of large cluster and another between nodes of multiple mid-
sized clusters with a common file server.

##  2 Related Works

Udy et al. summarized various approaches for oilfield production optimization
and compared algorithms including adjoint methods, particle swarm
optimization, simulated annealing and genetic algorithms [1]. Tanaka et al.
compared three different optimization approaches GA, PSO and ensemble-
optimization (EnOpt) for Olympus benchmark case [3]. Harb et al. used hybrid
evolutionary optimization method, black hole particle swarm optimization
(BHPSO), for simultaneously optimizing well count, location, type, and
trajectory [4]. Fonesca et al. used a hybrid of ensemble optimization EnOpt
and gradient free covariance matrix adaptation (CMA) method for well control
problem [5], while Bouzarkouna et al. used CMA-ES method for optimal well
placement under uncertainty [2].

##  3 Approach

Olympus benchmark case is a synthetic simulation model based on a typical
North Sea offshore oilfield, and is used in this work to compare the
effectiveness of different optimization algorithms [11]. Field is a 9km by 3km
by 50m in size with minor faults and limited aquifer support. Subsurface
uncertainty is captured by 50 realizations with different permeability,
porosity, net-to-gross and initial water saturation, of which only 10 random
realizations are used in this work to reduce computational requirement. GA and
PSO optimization algorithms are first tested for ensemble of 10 realizations,
and then going forward a representative realization is used for the remaining
studies. A realization in such studies stands for a random simulation model,
while ensemble of realizations is often used to capture subsurface
uncertainty.

Simulation model has 11 production wells and 7 injection wells which are all
on bottom-hole pressure (BHP) control (Figure 1). For the optimization study,
each well’s BHP control setting is changed independently at four different
pre-defined times during the 20 year simulation period. Thus, there are total
72 different control variables defined as continuous real parameter which can
take any value within a given range.

Two different objective functions are used, one as the weighted sum of
cumulative fluid production/injection while second as the net present value of
discounted cash-flow [1, 14]. This is done to have a rigorous analysis by
encoding varying degree of non-linearity among the two objective functions.
Weighted sum of cumulative fluid (WCF) is defined as:

$WCF=Q_{op}-0.1*(Q_{wp}+Q_{wi})$ (1)

where $WCF$ is weighted sum, $Q_{op}$ is total oil produced, $Q_{wp}$ is total
water produced, and $Q_{wi}$ is total water injected. Net present value (NPV)
of discounted cash-flow incorporates additional degree of non-linearity and is
defined as:

$NPV=\sum_{f=1}^{N_{f}}\sum_{\tau=1}^{N_{\tau}}\frac{Q_{f,\tau}C_{f}}{(1+d_{f})^{\tau-1}}$
(2)

where $N_{f}$ is the number of fluid component, $N_{\tau}$ is number of time
intervals, $Q_{f,\tau}$ is the volume component produced/injected in time
interval $\tau$, $C_{f}$ is the price/cost while $d_{f}$ is discount factor
for component $f$.

Figure 1: Olympus simulation model overview: initial oil saturation and well
locations

Parallelization is easy for evolutionary optimization algorithms as objective
function calculation for population samples can be done in parallel [3]. This
calculation includes making the reservoir simulation run for the population
sample (or simulation case), and then determining the required value. Message-
Passing Interface (MPI) is used in this study to distribute workload among
nodes and thus speed-up calculations. For a large cluster setting, master node
assigns simulation cases to different worker nodes of cluster and then
accumulates simulation results. A number of small clusters can also be used to
parallalize workload if a large cluster is not available or is not efficient.
For such a setting, master node on main cluster uses a common file server
accessible by smaller clusters which in turn make simulation runs.

##  4 Algorithm

Three different gradient-free, evolutionary optimization algorithms used in
this study are genetic algorithm (GA), particle swarm optimization (PSO), and
covariance matrix adaptation-evolution strategy (CMA-ES).

Genetic algorithm (GA) has been a popular optimization algorithm and uses
three basic steps for population refinement, namely selection, crossover and
mutation [9]. For a real-coded genetic algorithm (RCGA) used in this study,
real variables themselves are used for population refinement and objective
function calculation, as compared to binary encoding done for binary-coded
genetic algorithm (BCGA) [10]. In this work, RCGA is used to obviate the need
of encoding and decoding control variables. Crossover and mutation steps of GA
help the algorithm in exploring complex and multi-modal solution space
effectively.

Particle swarm optimization (PSO) uses collective intelligence based on
particle position and velocity. It is the exchange of information among
particles that leads to the emergence of global patterns [7, 8, 13]. Sample
update is done using inertia, personal-best, and global-best values along with
random components. It is essential to properly tune hyperparameters for both
GA and PSO algorithms to find good solutions and help optimization run
converge. Based on domain expertise, the algorithms can be made more
exploratory or more exploitative by tuning hyperparameters. The
hyperparameters can also be varied during an optimization run to have more
exploration in the beginning and more exploitation towards the end of the
optimization run.

Covariance matrix adaptation - evolution strategy (CMA-ES), on the other hand,
is a stochastic population-based optimization approach where population is
updated using weighted mean of samples from the previous generation and a
random draw from a multivariate normal distribution [5, 6]. Next, parameters
of multivariate normal distribution are updated using objective function
evaluations. CMA-ES attempts to learn a second-order model of the underlying
objective function and follows a natural gradient approximation of the
expected fitness, similar to the approximation of the inverse Hessian matrix
in Quasi-Newton methods [6, 18]. For multi-modal functions, where local optima
cannot be interpreted as perturbations to an underlying convex (unimodal)
topology, its performance can strongly decrease.

Table 1: Solving Rastrigin test function using GA and PSO Dimensions | Population size | Iterations | Best cost (GA) | Best cost (PSO)
---|---|---|---|---
2 | 40 | 100 | 0.99 | 7.7e-10
50 | 40 | 100 | 172.1 | 322.0
50 | 100 | 200 | 97.8 | 158.7
Figure 2: Cost history for solving Rastrigin test function using PSO algirithm

In continuous domains, the volume of search space grows exponentially with
dimensions, and is known as the curse of dimensionality. This phenomenon
causes a population size that is capable of adequate search space coverage in
low dimensional search spaces to degrade into sparse coverage in high
dimensional search spaces. GA and PSO algorithms work well for a problem with
smaller number of control variables (or dimensions), but have problem
converging when the number of control variables increase [12, 17]. To
demonstrate this, both GA and PSO are used for solving Rastrigin function with
varied number of dimensions. Rastrigin function has global minimum at zeros
for any given number of dimensions and the associated cost is also zero. Cost
is calculated as:

$f=10d+\sum_{i=1}^{d}[x_{i}^{2}-10\cos(2\pi x_{i})]$ (3)

where $d$ is the number of dimensions, $x_{i}$ is the value along dimension
$i$, and $f$ is the cost. As the dimensionality increases, both GA and PSO
require larger population size and iterations, but still do not manage to
converge to the global minimum (Table 1, Figure 2). While CMA-ES is not very
effective for multi-modal problems, it has been found to have excellent
convergence properties for unimodal problems. Thus once GA and PSO algorithms
have established the region of interest, CMA-ES can be used to achieve
convergence or further improvement. In general, effectiveness of any
optimization algorithm depends on the problem class under consideration, and
is commonly known as the no free lunch theorem for optimization [19].

In this work, GA and PSO are each first used with ensemble of 10 realizations
as well as a single representative realization to establish their efficacy
with such optimization problems. Going forward, optimization studies are done
only for single representative model in order to cut down on computational
cost. For using two optimization methods in a hybrid mode, one optimization
algorithm is first allowed to run for a number of generations, and then second
algorithm takes over. Two hybrid approaches, GA-CMA-ES and PSO-CMA-ES, are
tested and compared with standalone GA and PSO algorithms [14, 15, 16, 18].

##  5 Results

Olympus model has 11 producer wells and 7 injector wells, all on bottom-hole
pressure (BHP) constraint (Figure 1). Constraints can be changed independently
for each well at four predefined times giving a total of 72 control variables.
Well control settings for an ensemble of 10 realizations as well as a single
representative realization are first optimized using GA, and then steps are
repeated using PSO. Table 2 shows related parameter values for the
optimization runs. A discount factor of 0.08 was used for all fluid components
in NPV calculation. For ensemble case, all 10 realizations (or simulation
models) are simulated for any given population sample and results are
averaged. As can be seen in Figures 4 and 4, both GA and PSO algorithms are
able to optimize well control settings and yield higher weighted cumulative
fluid as compared to the default case. As sampling is randomized in both the
algorithms, every optimization run will proceed in unique fashion. Thus
multiple optimization runs are required for comparing different algorithms.
Going forward, only single representative reservoir model is optimized to cut
down on computational requirement.

Table 2: Optimization study parameters Population size | 40 | Oil price | $ 40 per barrel
---|---|---|---
Iterations (GA or PSO) | 100 | Produced water cost | $ 4 per barrel
Iterations (Hybrid) | 150 | Injected water cost | $ 2 per barrel

Figure 3: GA and PSO based ensemble optimization run with WCF objective

Figure 4: GA and PSO based single realization optimization run with WCF
objective

To further validate algorithm efficacy, objective function is next changed to
net present value (NPV) of discounted cash-flow, and the results for GA and
PSO based optimization runs are shown in Figure 5. Again algorithms seem to be
effective in handling such simulation based oil production optimization
problem. As compared to weighted sum of cumulative fluid (WCF) based
optimization, NPV based optimization yields muted increase over default/base
case because of the discounting.

Hybrid optimization run is made by transferring the resulting population from
one algorithm to the next as its starting population for further improvement.
Total 6000 simulation runs are made for such studies, out of which first 4000
runs (or 100 generations) are done using primary algorithm followed by 2000
runs (or 50 generations) using secondary algorithm. Two different hybrid
optimization settings are: GA followed by CMA-ES (GA-CMA-ES) and PSO followed
by CMA-ES (PSO-CMA-ES). Once region of interest is found using GA or PSO, CMA-
ES algorithm further refines the solution as it is a better algorithm for
local exploitation. Figures 7 and 7 show that CMA-ES is able to further
improve upon the solution, and thus GA-CMA-ES and PSO-CMA-ES approach of
hybrid optimization is better as compared to standalone GA or PSO. There is a
marked jump in objective function value after 4000 runs when CMA-ES takes over
from GA or PSO. It needs to be mentioned that CMA-ES does not have good
exploratory mechanism for multi-modal solution space, and hence it does not
perform as good as GA or PSO in standalone mode. Domain expertise is crucial
for success of such optimization studies as their design would vary depending
on underlying solution space and objective function.

Figure 5: GA and PSO based single realization optimization run with NPV
objective

Figure 6: Hybrid GA-CMA-ES and PSO-CMA-ES based single realization
optimization run with WCF objective

Figure 7: Hybrid GA-CMA-ES and PSO-CMA-ES based single realization
optimization run with NPV objective

Population based algorithms can be easily programmed in parallel mode with
objective function for each sample being calculated using a separate node. In
the calculation of objective function for a given population sample, compute
intensive step is making the reservoir simulation run. Either a large cluster
can be used to divide workload among cores, or multiple mid-size clusters can
be used with a common file server. By reading from and writing to a common
file server (Figure 8), multiple clusters can be used to carry-out compute
intensive simulation runs, thus reducing turn-around time.

##  6 Conclusion

Evolutionary algorithms are efficient at solving oilfield well control
optimization problem. Two such algorithms are used in hybrid mode to explore
multi-modal solution space and attain convergence. GA-CMA-ES and PSO-CMA-ES
hybrid optimization approach perform better as compared to using GA or PSO for
such optimization problem. For the hybrid optimization, GA or PSO algorithm is
first used to explore solution space and identify region of interest, and then
CMA-ES is used to further refine the solution. Domain expertise is crucial for
success of such optimization studies as their design would vary depending on
solution space and objective function. Cloud computing and parallel
programming are used to handle the computational workload. Two different
cluster settings, one with single large cluster and other with multiple mid-
sized clusters, have been used to handle high computational workload.

## Acknowledgments and Disclosure of Funding

Research was made possible by generous cloud credits provided by Amazon Web
Services, DigitalOcean, Google Cloud and IBM Cloud. Olympus benchmark
reservoir model was made available by the Netherlands Organization for Applied
Scientific Research (TNO).

Figure 8: Compute resource setup with multiple clusters and common file-server

## References

[1] Udy, J., Hansen, B., Maddux, S., Petersen, D., Heilner, S., Stevens, K.,
Lignell, D. & Hedengren, J.D. (2017) Review of field development optimization
of waterflooding, EOR, and well placement focusing on history matching and
optimization algorithms. Processes 5(34).

[2] Bouzarkouna, Z., Ding, D.Y. & Auger, A. (2012) Well placement optimization
under uncertainty with CMA-ES using the neighborhood. $13^{th}$ European
Conference on the Mathematics of Oil Recovery.

[3] Tanaka, S., Wang, Z., Dehghani, K., He, J., Velusamy, B. & Wen, X. (2018)
Large scale field development optimization using high performance parallel
simulation and cloud computing technology. SPE Annual Technical Conference and
Exhibition.

[4] Harb, A., Kassem, H. & Ghorayeb, K. (2019) Black hole particle swarm
optimization for well placement optimization. Computational Geosciences, 24
(6):1979-2000.

[5] Fonseca, R.M., Leeuwenburgh, O., Van den Hof, P.M.J. & Jansen, J.D. (2013)
Improving the ensemble optimization method through covariance matrix
adaptation (CMA-EnOpt). SPE Reservoir Simulation Symposium.

[6] Hansen, N. (2016) The CMA evolution strategy: a tutorial. arXiv preprint
arXiv:1604.00772.

[7] Wang, D., Tan, D. & Liu, L. (2018) Particle swarm optimization algorithm:
an overview. Soft Computing, 22:387-408.

[8] Sengupta, S., Basak, S. & Peters, R.A. (2019) Particle swarm optimization:
a survey of historical and recent developments with hybridization
perspectives. Machine Learning and Knowledge Extraction, 1:157-191.

[9] Katoch, S., Chauhan, S.S. & Kumar, V. (2021) A review on genetic
algorithm: past, present, and future. Multimedia Tools and Applications,
80:8091-8126.

[10] Chuang, Y.C., Chen, C.T. & Hwang, C. (2015) A real-coded genetic
algorithm with a direction-based crossover operator. Information Sciences,
305:320-348.

[11] Fonseca, R.M., Rossa, E.D., Emerick, A.A., Hanea, R.G. & Jansen, J.D.
(2020) Introduction to the special issue: overview of OLYMPUS optimization
benchmark challenge. Computational Geosciences, 24 (6):1933-1941.

[12] Chen, S., Montgomery, J. & Bolufé-Röhler, A. (2015) Measuring the curse
of dimensionality and its effects on particle swarm optimization and
differential evolution. Appl Intell, 42:514-526.

[13] Zhang, Y., Wang, S. & Ji, G. (2015) A comprehensive survey on particle
swarm optimization algorithm and its applications. Mathematical Problems in
Engineering, 2015.

[14] Baumann, E.J.M, Dale, S.I. & Bellout, M.C. (2020) FieldOpt: A powerful
and effective programming framework tailored for field development
optimization. Computers and Geosciences, 135.

[15] Xu, P., Luo, W., Lin, X., Qiao, Y. & Zhu, T. (2019) Hybrid of PSO and
CMA-ES for global optimization. IEEE Congress on Evolutionary Computation.

[16] Santos, J., Ferreira, A. & Flintsch, G. (2019) An adaptive hybrid genetic
algorithm for pavement management. International Journal of Pavement
Engineering, 20 (3):266-286.

[17] Clerc, M. & Kennedy, J. (2002) The particle swarm - explosion, stability,
and convergence in a multidimensional complex space. IEEE Transactions on
Evolutionary Computation, 6 (1):58-73.

[18] Muller, C.L., Baumgartner, B. & Sbalzarini, I.F. (2009) Particle swarm
CMA evolution strategy for the optimization of multi-funnel landscapes. IEEE
Congress on Evolutionary Computation, pp. 2685-2692.

[19] Wolpert, D.H., & Macready, W.G. (1997) No free lunch theorems for
optimization. IEEE Transactions on Evolutionary Computation, 1 (1):67-82.


# Three-Polarizer Treatment of Linear Polarization in Coronagraphs and
Heliospheric Imagers

Craig E. DeForest Southwest Research Institute
1050 Walnut Street, Suite 300
Boulder, CO 80302, USA  Daniel B. Seaton Southwest Research Institute
1050 Walnut Street, Suite 300
Boulder, CO 80302, USA  Matthew J. West Southwest Research Institute
1050 Walnut Street, Suite 300
Boulder, CO 80302, USA Craig E. DeForest<EMAIL_ADDRESS>
(Received October 4, 2021)

###### Abstract

Linearly polarized light has been used to view the solar corona for over 150
years. While the familiar Stokes representation for polarimetry is complete,
it is best matched to a laboratory setting and therefore is not the most
convenient representation either for coronal instrument design or for coronal
data analysis. Over the last 100 years of development of coronagraphs and
heliospheric imagers, various representations have been used both for direct
measurement and analysis. These systems include famous representations such as
the ($B$, $pB$) system that is analogous to the Stokes system in solar
observing coordinates, and also internal representations such as in-instrument
Stokes parameters with fixed or variable “vertical” direction, and brightness
values through a particular polarizing optic or set thereof. Many polarimetric
instruments currently use a symmetric three-polarizer measurement and
representation system, which we refer to as “($M,~{}Z,~{}P$)”, to derive the
($B$, $pB$) or Stokes parameters. We present a symmetric derivation of
($B,~{}pB$) and Stokes parameters from ($M,~{}Z,~{}P$), analyze the noise
properties of ($M,~{}Z,~{}P$) in the context of instrument design, develop
($M,~{}Z,~{}P$) as a useful intermediate system for data analysis including
background subtraction, and draw a helpful analogy between linear polarimetric
systems and the large existing body of work on photometric colorimetry.

Solar K Corona (2042) — Heliosphere (711) — Polarimetric Instrumentation
(1277) — Polarimetry (1278)

††journal: Astrophys. J.

##

> Light is a thing that cannot be reproduced but must be represented by
> something else: by color. (Paul Cézanne)

##  1 Introduction

The solar corona is linearly polarized (Arago, 1843). That property has been
exploited over nearly a century of coronal observations: both for background
removal in coronagraphs (Lyot, 1930) and for 3D analysis (Poland & Munro,
1976; de Koning & Pizzo, 2011; DeForest et al., 2017).

The degree of coronal polarization is conventionally reported via two
parameters: a $B$ (“unpolarized brightness”) parameter and its counterpart
$pB$ (“polarized brightness”); and both parameters may be mapped over an image
plane to create separate “$B$ images” and “$pB$ images”. Because the corona is
a distributed object, it is best tracked via its radiance (delivered optical
power per unit area, per unit solid angle), and therefore “brightness” and
“radiance” are synonyms. The K-corona is visible primarily via Thomson
scattering and therefore it is polarized perpendicular to the plane containing
the observer, the Sun, and the scattering point. At each point in a two
dimensional image plane, the K corona is thus polarized perpendicular to a
line extending through the imaged center of the Sun and the given point
(Figure 1a); this direction is not only predicted by theory (e.g., Billings,
1966; Howard & Tappin, 2009), but also readily verified through direct
measurement (e.g., Filippov et al., 1994). $pB$ may thus be calculated as

$pB=B_{T}-B_{R},$ (1)

where $B_{T}$ is the radiance observed through a linear polarizer oriented
tangentially to a solar-concentric circle passing through the image point of
interest, and $B_{R}$ is the radiance observed through a linear polarizer
oriented radially to the Sun through the same point (e.g., Minnaert, 1930;
Altschuler & Perry, 1972). Whether Equation 1 is accidental or fundamental to
the definition of $pB$ is a matter of ongoing historical confusion and even
dispute. This is explored in an Appendix to this article, which introduces a
${}^{\circ}pB$ and ${}^{\perp}pB$ to distinguish historical usage. In this
work we treat Equation 1 as fundamental; this treatment implies that $pB$ is
similar to the Stokes $Q$ parameter (e.g., Hecht & Zajac, 1974), but in a
coordinate system that is rotated relative to the instrument.

The ($B,~{}pB$) representation of coronal polarization is convenient because
it matches the observing geometry, and in fact many coronagraphs have been
constructed to measure directly that element of linear polarization (e.g.,
Altschuler & Perry, 1972). However, because the direction of “radial” varies
in the image plane, this adds complexity to the instrument itself. Recent
spaceborne coronagraphs have instead measured the full linear polarization
state of incident light in the instrument frame (e.g., Brueckner et al., 1995;
Howard et al., 2008). This requires capturing the first three Stokes
parameters ($I,~{}Q,~{}U$) at each location of the image plane, either
directly or indirectly via a representative measurement.

Figure 1: Four panels describe linear polarization analysis in a coronal
context. (a) An object near the Sun has position angle $\alpha$ and
polarization angle $\beta$; $B_{R}$ and $B_{T}$ describe brightnesses through
radially and tangentially aligned polarizers. (b) Stokes $Q$ and $U$ describe
polarization in the “+” and “$\times$” directions relative to the instrument
(or solar North). (c) $pB$ and $pB^{\prime}$ are Stokes parameter analogs in
the solar observing reference frame. (d) Observing polarization through three
polarizers mutually separated by $\pi/3$ radians is sufficient to capture the
polarization state ($I,Q,U$) or, equivalently, ($B,pB,pB^{\prime}$). Although
Thomson-scattered light is polarized in the $B_{T}$ direction, we show the
polarization vector slightly misaligned, to emphasize the general case.

The Stokes ($I,~{}Q,~{}U$) parameters are defined in terms of exposures
through individual crossed polarizers in a laboratory, and the textbook
approach to measurement involves four exposures through crossed polarizers in
the fixed instrument reference frame: horizontal/vertical for Q and diagonal
for U (Figure 1b; Hecht & Zajac, 1974). However, only three independent
measurements are necessary to determine the linear polarization states of a
beam of light. This was the basis of the triplet polarizing camera used by
Öhman (1947) to observe the 1945 eclipse, and Öhman attributed the technique
to an earlier analysis by Fesenkov (1935). The technique is described briefly
in Chapter 4 of Billings (1966), and by Newkirk et al. (1970); these authors
both used the Stokes parameters as an intermediary system, rather than
treating the triplet analysis as primary. Polarization triplet analysis has
been and is still used routinely with coronal images from Skylab (Poland &
Munro, 1976), SOHO/LASCO (Brueckner et al., 1995) and STEREO/SECCHI (Howard et
al., 2008), and is planned for other in-development instruments.

Despite its common use, full descriptions of the polarizer triplet technique,
including instrumental effects, are elusive in the current literature; and
those that exist are somewhat unsatisfying because, in translating through the
Stokes system, the analytic work loses the symmetry of the polarizer triplet
system. In this article, we derive formulae to find directly the
($B,~{}pB,~{}pB^{\prime}$) or Stokes ($I,~{}Q,~{}U$) representations of linear
polarization from polarizer triplet data; present an analytic noise analysis
for the most likely sources of noise or systematic error in a real instrument
using a polarizer triplet; introduce the use of an ($M,~{}Z,~{}P$) system of
“virtual polarizer triplets” (Figure 1d) to represent the linear polarization
state of light in frames other than the observing frame; and draw a surprising
but helpful analogy between the major systems of linear polarimetry and
corresponding systems for representing color.

##  2 Definitions

Table 1 defines several quantities, including radiances, polarization
parameters, and angles, that are relevant to a polarizing coronagraph or
heliospheric imager. These quantities are used throughout Section 3, and are
collected here for reference. Figure 1 illustrates the angles and associated
angles of the principal quantities in Table 1.

Quantity | Expression | Definition
---|---|---
$B_{\\{T,R\\}}$ | – | Radiance through an ideal polarizer oriented
tangentially or radially relative to the Sun
$B$ | $B_{T}+B_{R}$ | “Unpolarized brightness” (radiance)
$pB$ | $B_{T}-B_{R}$ | Coronal “Polarized brightness”
$\alpha$ | – | Solar position angle of an image point
$\beta$ | – | Direction of polarization
$\theta$ | – | Polarizer angle (also subscripted as $\theta_{i}$)
$\theta_{i}$ | $\theta+\\{-1,0,1\\}\pi/3$ | One of three angles in an ($M,~{}Z,~{}P$) triplet
$\phi$ | – | Second polarizer angle
$B_{\theta}$ | – | Radiance through a polarizer at angle $\theta$
$B_{i}$ | – | Radiance through a polarizer at angle $\theta_{i}$
$B_{\\{|,-,/,\backslash\\}}$ | – | Radiance through an ideal polarizer
oriented at $n\pi/4$ ($n\in\mathbb{Z})$
$I$ | $B_{|}+B_{-}$ | Stokes I – synonym for $B$ (sum of
$B_{\backslash}+B_{/}$ | radiance through any two perp. polarizers)
$Q$ | $B_{|}-B_{-}$ | Stokes Q (in instrument frame)
$U$ | $B_{\backslash}-B_{/}$ | Stokes U (in instrument frame)
$p$ | $\sqrt{Q^{2}+U^{2}}/I$ | Polarization fraction (note: for
Thomson scattering, $pB=(p)(B)$)
$S_{i}$ | $\sin\left[2\left(\theta_{i}-\alpha\right)\right]$ | Convenient abbreviation for sine expression
$C_{i}$ | $\cos\left[2\left(\theta_{i}-\alpha\right)\right]$ | Convenient abbreviation for cosine expression
Table 1: Useful quantities in a linear polarimetric coronagraph. In this
context, Stokes $I$ is a radiance, not an intensity. Lists in curly braces
enumerate options for the relevant expression.

##  3 Linear Polarization from Three Polarizers

We start by deriving a formula for $pB$ and $B$ at position angle $\alpha$,
with coronal light polarized along the direction $\alpha+\pi/2$. This further
leads to expressions for extracting the direction of polarization, and the
full linear portion of the Stokes representation. Since $\alpha$ is the angle
of a radial line, the polarization direction $\beta$ is just $\alpha+\pi/2$
(tangential) for Thomson scattered light in the corona.

Solving the definitional expressions for $B$ and $pB$ in Table 1, for $B_{T}$
and $B_{R}$:

$B_{T}=\frac{B+pB}{2}\hskip 21.68121ptand\hskip
21.68121ptB_{R}=\frac{B-pB}{2}$ (2)

Admitting the light through a polarizer at angle $\theta$ projects the
electric field to the new direction; intensity is the square of the electric
field amplitude, so:

$B_{\theta}=B_{T}\left[\sin\left(\theta-\alpha\right)\right]^{2}+B_{R}\left[\cos\left(\theta-\alpha\right)\right]^{2}.$
(3)

Substituting and applying the double-angle cosine formula:

$B_{\theta}=\frac{1}{2}\left\\{B-pB\cos\left[2\left(\theta-\alpha\right)\right]\right\\},$
(4)

which is the projection equation for partially polarized light. Solving for
$pB$ in terms of $B$ and $B_{\theta}$ yields:

$pB=\frac{B-2B_{\theta}}{\cos\left[2\left(\theta-\alpha\right)\right]}~{},$
(5)

which resolves $pB$ in terms of the radiance through a single arbitrary
polarizer and also the total radiance with no polarizers in the beam. Equation
5 is problematic, because the denominator is small when $(\theta-\alpha)$ is
near $\pm\pi/4$. If multiple polarizers are available at different angles, one
can regularize Equation 5 by taking an average, and weighting each term by the
square of the offending cosine. Adopting the $C_{i}$ abbreviation (Table 1):

$pB=\frac{\sum_{i}\left\\{\left(B-2B_{i}\right)C_{i}\right\\}}{\sum_{i}C_{i}^{2}}$
(6)

which is numerically stable if the $\theta_{i}$ are not all separated by
intervals of $n\pi/2$. Choosing the $(M,Z,P)$111(for “Minus, Zero, Plus”)
basis of polarizers at $\theta-\pi/3$, $\theta$, and $\theta+\pi/3$, the
denominator sum evaluates to 3/2 and the $B$ terms sum to zero, so:

$pB=-\frac{4}{3}\sum_{i\in\\{M,Z,P\\}}\left(B_{i}C_{i}\right)$ (7)

Meanwhile, solving Equation 4 for $B$ yields:

$B=2B_{\theta}+C\,pB,$ (8)

where the $C$ is not subscripted because here it appears outside of a sum.
Equation 8 is not particularly useful by itself – but averaging over the three
polarizer positions of the $(M,Z,P)$ system eliminates the cosine through a
trigonometric identity, yielding:

$B=\frac{2}{3}\sum_{i\in\\{M,Z,P\\}}B_{i}$ (9)

Equations 7 and 9 give $B$ and $pB$ in closed and symmetrical form, given
radiance data through three polarizers with relative angles of $(-\pi/3)$,
$(0)$, and $(+\pi/3)$ radians (60$\degree$ separation) relative to a baseline
angle $\theta$, with the assumption that the direction of linear polarization
is perpendicular to $\alpha$.

Note that, although $\alpha$ is defined (in Table 1) as position angle around
the Sun, nothing in Equations 7 and 9 requires any particular value of
$\alpha$ – or, for that matter, any particular orientation of the main
polarization direction $\theta$. In fact, in the special case where
$\alpha=\pi/2$, $B$ and $pB$ are just the Stokes $I$ and $Q$ parameters in the
solar reference frame. Working by analogy to the Stokes parameters, we can
substitute $\alpha\rightarrow\alpha+\pi/4$ into Equation 7; this yields a
similar quantity that reduces to Stokes U in the same circumstance:

$pB^{\prime}=-\frac{4}{3}\sum_{i\in\\{M,Z,P\\}}\left(B_{i}S_{i}\right)$ (10)

where the new quantity $pB^{\prime}$ bears the same relationship to $pB$ as
Stokes $U$ does to Stokes $Q$. In systems where $pB^{\prime}$ is important, we
must generalize Equation 4 to include $pB^{\prime}$:

$B_{\theta}=\frac{1}{2}\left\\{B-C\,pB-S\,pB^{\prime}\right\\}$ (11)

where the difference between Equation 4 and Equation 11 is that the former
assumes tangential polarization while the latter treats arbitrary linear
polarization. Repeating the derivation of Equation 7 using Equation 11 instead
of 4 yields the same result, because the $C_{i}$ weighted summation in
Equation 7 eliminates the $S_{i}\,pB^{\prime}$ terms from Equation 11.

The $B$, $pB$, $pB^{\prime}$ system is an analog of the Stokes $I$, $Q$, $U$
system, but rotating with $\alpha$ around the Sun rather than being fixed in
the instrument frame of reference; the familiar Stokes parameters may be
recovered anywhere by substituting $\alpha=\pi/2$ into Equations 7, 9, and 10.
In the special case where the polarizer triplet is aligned with the
instrument, one may also set $\theta=0$ and arrive, after many cancellations,
at the simplified expressions:

$Q=\frac{2}{3}\left(2B_{Z}-B_{M}-B_{P}\right)$ (12)

and

$U=\frac{2}{\sqrt{3}}\left(B_{P}-B_{M}\right),$ (13)

which (together with the identity $B=I$) define the Stokes $I,Q,U$ triplet
relative to the reference Z polarizer in the $M,Z,P$ system, with a minimum of
calculation. Equations 12 and 13 reproduce the derivations presented by Öhman
(1947) and Billings (1966).

The remaining important system to represent linear polarization is
$(B,\theta,p)$, which uses the overall brightness $B$, direction of
polarization $\theta$, and normalized degree of polarization $p$. The
direction of polarization is available from $pB$ and $pB^{\prime}$ (or,
equivalently, from Stokes $Q$ and $U$). Differentiating Equation 11 to find
the maximum value of $B_{\theta}$:

$\left.\frac{\partial
B_{\theta}}{\partial\theta}\right|_{\theta=\theta_{max}}=\left\\{S\,pB-C\,pB^{\prime}\right\\}=0~{},$
(14)

so:

$\theta_{max}=\frac{1}{2}\arctan\left(\frac{pB^{\prime}}{pB}\right)+\frac{\pi}{2}+\alpha~{},$
(15)

where the four-quadrant arctan is implied and selects maxima rather than
minima, and the $\pi/2$ arises from canceling negative signs in the numerator
and denominator of the arctan (adding $\pi$ to the result of the four-quadrant
arctan). Making use of Equation 15, substituting the arctan into Equation 11,
and recognizing two trigonometric cancellations yields the expected formula
for $p$:

$p\equiv\frac{2B_{\theta,max}-B}{B}=\frac{\sqrt{pB^{2}+pB^{\prime 2}}}{B}.$
(16)

Equations 15 and 16 recover the textbook formulas for $\theta_{max}$ and $p$
from the Stokes formalism, in the context of the $(B,pB,pB^{\prime})$ system
instead of $(I,Q,U)$.

###  3.1 Error Sources in Three-Polarizer Polarimetry

Polarimetry is affected by multiple error and noise sources that are specific
to the measurement. Here we derive expressions for calculating the error or
noise associated with three-polarizer measurements of $B$ and $pB$, given
important photometric or mechanical tolerances of the instrument. The three
major potential sources of error are photometric noise, polarizer misalignment
errors, and finite polarizer effectiveness. We consider noise as an
approximately Gaussian-distributed random variable, with mean value of zero,
added to the true value of each quantity.

To keep track of which quantities represent physical truth (with no noise nor
measurement error) and which represent inference from measurement (with noise
and other error sources included), we introduce an overbar to indicate
directly or indirectly measured parameters:

$\bar{X}\equiv X+\Delta X~{},$ (17)

where $X$ is an observable parameter such as $B$ or $pB$, $\bar{X}$ is the
observed value, and $\Delta X$ is a noise term.

Furthermore, we take the indices $i$ and $j$ to run over the ($M,~{}Z,~{}P$)
polarizer angles when mentioned in a summation, in keeping with Equations 7,
9, and 10.

####  3.1.1 Photometric errors

Propagating photometric error and noise through multi-polarizer measurements
is straightforward. Photometry is most frequently limited by Poisson-
distributed photon shot noise or similar uncorrelated noise sources.
Systematic errors can also contribute to total error; we neglect these terms,
which is equivalent to treating systematics in each channel as randomly
distributed across channels. To each value $B_{i}$ we add a noise term $\Delta
B_{i}$ that we treat as a sample of a random variable. To propagate noise from
the three samples, we take the partial derivative of Equations 7, 9, and 10
with respect to an arbitrary polarizer brightness, then sum the $\Delta B_{i}$
terms in quadrature. This approach works because Equations 7, 9, and 10 are
linear in $B_{i}$; and the $B_{i}$ noise samples are also taken to be
uncorrelated. Equation 9 is straightforward to differentiate:

$\frac{\partial\,\bar{B}}{\partial\bar{B}_{i}}=\frac{2}{3}.$ (18)

In the most common case for image detectors, there is both a signal-
independent noise component and a Poisson-statistics photon counting noise
component, adding in quadrature. In that case, we have

$\Delta B_{i}=\sqrt{(\Delta_{c}B)^{2}+B_{0}B_{i}},$ (19)

where $\Delta_{c}B$ is the signal-independent noise term (such as detector
dark noise or read noise), $B_{0}$ is an instrument-specific constant of
proportionality relating radiance units to detection quanta at an individual
pixel or detector, and $\sqrt{B_{0}B_{i}}$ is the Poisson noise term that
arises from counting statistics. Merging the noise terms explicitly using
Equation 9 and quadrature summmation yields

$\Delta
B=\frac{2}{3}\sqrt{\sum_{i}\left[\left(\Delta_{c}B\right)^{2}+B_{0}B_{i}\right]}~{},$
(20)

or, substituting from Equation 11,

$\Delta
B=\frac{2}{3}\sqrt{3(\Delta_{c}B)^{2}+\frac{B_{0}}{2}\sum_{i}\left(B-C_{i}pB-
S_{i}pB^{\prime}\right)}~{}.$ (21)

The $C_{i}$ terms and the $S_{i}$ terms sum to zero, leaving

$\Delta B=\frac{2}{\sqrt{3}}\sqrt{(\Delta_{c}B)^{2}+B_{0}B/2}~{},$ (22)

which, in the special case of unpolarized light (so that $B_{i}=B/2$ for each
$i$) would reduce, for each $i$, to

$\Delta B=\frac{2}{\sqrt{3}}\Delta B_{i}~{},$ (23)

which arises from Equation 22 by substituting $B_{i}$ for $B/2$ and noticing
that the sum matches the right-hand side of Equation 19.

Because Equation 22 is independent of $pB$ and $pB^{\prime}$ for polarizer
triplets, the photometric $\Delta B$ may be calculated for the non-polarized
case and remains valid in the polarized case, if the primary noise sources
scale independently of the signal and/or with Poisson counting statistics. The
form of Equation 23 highlights that the photometric noise scales as expected
from individual samples being averaged together, shrinking by a factor of
$3^{1/2}$ when the polarizer triplet exposures are merged.

Turning to $pB$, the derivative of Equation 7 is just:

$\frac{\partial\,p\bar{B}}{\partial\bar{B}_{i}}=-\frac{4}{3}C_{i}.$ (24)

The quadrature sum is:

$\Delta pB=\frac{4}{3}\sqrt{\sum_{i}C_{i}^{2}\Delta
B_{i}^{2}}=\frac{4}{3}\sqrt{\frac{3}{2}\left(\Delta_{c}B\right)^{2}+B_{0}\sum_{i}\left(C_{i}^{2}B_{i}\right)}$
(25)

Substituting with Equation 11 to eliminate the $B_{i}$ term yields

$\Delta
pB=\frac{4}{3}\sqrt{\frac{3}{2}\left(\Delta_{c}B\right)^{2}+\frac{3}{2}\frac{B_{0}B}{2}-\frac{B_{0}}{2}\left(pB\sum_{i}C_{i}^{3}+pB^{\prime}\sum_{i}C_{i}^{2}S_{i}\right)}$
(26)

where the summation terms oscillate in $\theta-\alpha$, with individual
amplitude $3/4$ and frequency 6 per full circle in $\theta$. It is useful to
bound the right-hand-side of Equation 26. The maximum value of the square root
occurs when the oscillating term has its minimum (most negative) possible
value. That occurs at the extremal value of $pB$ and/or $pB^{\prime}$, such
that $pB^{2}+pB^{\prime 2}=B^{2}$. Setting $pB=B$ and $pB^{\prime}=0$, and
selecting the minimum value of the $C_{i}^{3}$ sum, yields

$\Delta
pB\leq\frac{2\sqrt{2}}{\sqrt{3}}\sqrt{\left(\Delta_{c}B\right)^{2}+\frac{3}{4}B_{0}B}~{},$
(27)

which can be simplified further. Adding $(\Delta_{c}B)^{2}/2$ to the interior
is allowed in the inequality, and allows a much simpler expression at the cost
of slightly expanding the bound:

$\Delta pB\leq 2\sqrt{\left(\Delta_{c}B\right)^{2}+B_{0}B/2}=\sqrt{3}\,\Delta
B~{}.$ (28)

Systems that have nonzero values of $pB^{\prime}$ have the same upper bound
for $\Delta pB$, because when the total relative polarization is maximized,
mixing $pB$ and $pB^{\prime}$ is equivalent to rotating the polarizers
(changing $\theta$), which does not affect the minimization of the summation
terms across $\theta-\alpha$.

The upper bound for $pB^{\prime}$ clearly must be the same as for $pB$,
because the two quantities are related by a rotation – and neither $\theta$
nor $\alpha$ appear in Equation 28. We therefore immediately write:

$\Delta pB^{\prime}\leq\sqrt{3}\,\Delta B~{}.$ (29)

In sum: although the photometric noise level in the derived value of $pB$
(and, by extension, other related quantities $pB^{\prime}$, $Q$, and $U$) from
Equation 7 varies with overall polarizer angle $\theta$, a simple upper bound
exists that is independent of $\theta$ and (equivalently) position angle
$\alpha$. Photometry in polarized brightness parameters is worse than the
unpolarized overall photometry, by a factor of up to $\sqrt{3}$.

####  3.1.2 Polarizer Misalignments

Instrument polarizers are mechanical devices and subject to alignment
tolerance. The triplet $(M,Z,P)$ system (introduced between Equation 6 and
Equations 7 and 9) relies on a trigonometric identity to simplify the weighted
average in Equation 6. Errors in polarizer angle for a particular exposure
translate directly to errors in $B$ and $pB$ independent of photometric noise.
Even assuming either perfect alignment or proper angular calibration (using,
e.g, Equation 6 or perturbation analysis on Equation 7), individual polarized
intensity images will be polarized at an angle that differs slightly from
nominal, due to noise in the alignment process. Following the same methodology
as in Section 3.1.1, we characterize polarimetric response to this type of
noise by treating noise-associated error in polarizer angle as a random
variable. We then propagate the noise from the foundational equations to
Equations 7, 9, and 10, by partially differentiating.

We start by characterizing $B_{\theta}$’s $\theta$-dependence. From Equation
11,

$\frac{\partial B_{i}}{\partial\theta_{i}}=S_{i}\,pB-C_{i}\,pB^{\prime}$ (30)

and, from Equation 9 (and making use of the same identity as in Equation 7),

$\frac{\partial\bar{B}}{\partial\theta_{i}}=\frac{2}{3}\left\\{S_{i}\,pB-
C_{i}\,pB^{\prime}\right\\}~{}.$ (31)

Treating three $\Delta\theta_{i}$’s in quadrature (and making use of the same
identity as in Equation 26):

$\Delta B=\frac{\sqrt{2}}{\sqrt{3}}\sqrt{pB^{2}+pB^{\prime
2}}\Delta\theta~{}.$ (32)

Equation 32 describes a noise term in $B$ whose magnitude is dependent on the
total degree of linear polarization regardless of direction.

Applying the same process to $pB$, we have from Equation 6 and Equation 30:

$\begin{split}\frac{\partial
p\bar{B}}{\partial\theta_{i}}=&\frac{1}{\sum_{k}C_{k}^{2}}\left\\{\sum_{j}(B-2B_{j})(-2S_{i}\delta_{ij})\right\\}\\\
&+\frac{1}{\sum_{k}C_{k}^{2}}\left\\{\sum_{j}(2/3)(S_{j}\,pB-
C_{j}\,pB^{\prime})C_{j}\right\\}\\\
&+\frac{1}{\sum_{k}C_{k}^{2}}\left\\{\sum_{j}\left(-2(S_{j}\,pB-
C_{j}\,pB^{\prime})\delta_{ij}C_{j}\right)\right\\}\\\
&-\frac{2(C_{i})(-2S_{i})}{\left[\sum_{k}C_{k}^{2}\right]^{2}}\sum_{j}\left\\{(B-2B_{j})C_{j}\right\\}~{},\end{split}$
(33)

where $\delta_{ij}$ is the Kronecker $\delta$. Substituting and, where
possible, eliminating those terms that sum to zero:

$\frac{\partial
p\bar{B}}{\partial\theta_{i}}=\frac{2}{3}\left\\{-2(B-2B_{i})S_{i}-2(S_{i}\,pB-
C_{i}\,pB^{\prime})C_{i}+4C_{i}S_{i}\,pB\right\\}.$ (34)

Solving Equation 11 for the quantity $B-2B_{\theta}$, and substituting in to
simplify Equation 34 further:

$\frac{\partial
p\bar{B}}{\partial_{\theta_{i}}}=\frac{2}{3}\left\\{-2\left(C_{i}\,pB+S_{i}\,pB^{\prime}\right)S_{i}-2\left(S_{i}\,pB-
C_{i}\,pB^{\prime}\right)C_{i}+4C_{i}S_{i}\,pB\right\\}~{},$ (35)

or, gathering terms and applying the double-angle identity:

$\frac{\partial
p\bar{B}}{\partial{\theta_{i}}}=\frac{4}{3}\cos\left[4\left(\theta_{i}-\alpha\right)\right]pB^{\prime}~{}.$
(36)

Equation 36 describes a noise coupling term in $p\bar{B}$ that is proportional
to $pB^{\prime}$. Therefore, because $pB^{\prime}$ is negligible in the bright
inner corona, small errors in polarizer angle have negligible effect on the
$pB$ calculations in that field. This surprising result helps explain why
triplet polarization has performed well in historical coronagraph data: in the
inner few apparent solar radii of the corona, such systems are insensitive to
small polarizer misalignments. The corresponding error term is as curiously
silent as a dog that does nothing in the night (Doyle, 1892).

Applying three $\Delta\theta_{i}$’s in quadrature yields:

$\Delta pB=\frac{2\sqrt{2}}{\sqrt{3}}\,pB^{\prime}\,\Delta\theta~{}.$ (37)

Arguing from symmetry, we can also immediately write:

$\Delta pB^{\prime}=\frac{2\sqrt{2}}{\sqrt{3}}\,pB\,\Delta\theta~{}.$ (38)

In sum: polarizer misalignment in principle produces photometric error in
derived polarization parameters. The resulting error in unpolarized brightness
$B$ depends on total polarization and not, to first approximation order, on
its direction. The corresponding error in $pB$ calculation depends only on the
cross mode polarization $pB^{\prime}$, and vice versa.

####  3.1.3 Polarizer Effectiveness

The analysis leading to Equations 7, 9, and 10 assumed perfect polarizers. In
practice, real polarizers are not perfectly effective and instead have a
finite extinction coefficient. A real instrument with physical polarizers will
produce observed polarizer brightnesses $\bar{B}_{\theta}$ with:

$\bar{B}_{\theta}=B_{\theta}+\epsilon B~{},$ (39)

where $\epsilon$ is a small leakage coefficient (and of course an overall
efficiency factor is removed via photometric calibration). If $\epsilon$ is
known, it can be calibrated out. Equation 39 propagates into Equation 7 as:

$pB=-\frac{4}{3}\sum_{i}\left[\left(\bar{B}_{i}-\epsilon_{i}B\right)C_{i}\right]~{}.$
(40)

In the case where all the $\epsilon_{i}$ terms are equal, Equation 40 reduces
to Equation 7; i.e. finite polarizer extinction effects do not affect $pB$
measurement nor, by symmetry, $pB^{\prime}$ measurement, provided that they
are constant across polarizers.

Likewise, substituting $B=\bar{B}_{i}-\epsilon_{i}B$ into Equation 9, we
obtain:

$B=\frac{2}{3}\frac{\sum_{i}\bar{B}_{i}}{\left(1+\frac{2}{3}\sum_{i}\epsilon_{i}\right)}~{},$
(41)

so that, if all $\epsilon_{i}$ terms are equal, then the calculated $B_{i}$
value must be scaled by a factor of $(1-2\epsilon)$, to first approximation
order, relative to the naïve perfect-polarizer formula in Equation 9.

In sum: the polarizer extinction coefficient $\epsilon^{-1}$ does not affect
inferred values of $pB$ nor of $pB^{\prime}$, provided that it is constant
across polarizer positions. The extinction coefficient does affect the
inferred value of $B$, which can be corrected provided that $\epsilon$ is
known. As an example, polarizers with $\epsilon^{-1}=200$ would incur a $1\%$
error in the ratio $pB/B$ from this source alone, if not corrected using
Equation 41.

##  4 Virtual Polarizer Triplets

So far we have established the theory of measuring the ($B$, $pB$,
$pB^{\prime}$) coronal polarization parameters or, equivalently, the ($I$,
$Q$, $U$) Stokes parameters, using polarizer triplets. But the triplet
formulation is useful not only as a means of measurement, but also as a means
to represent the linear polarization state of light. A “virtual polarizer
triplet” formulation represents light as three brightness parameters through
three ideal polarizers separated by $\pi/3$ radians ($60^{\circ}$) each: the
($M$, $Z$, $P$) representation.

In the context of coronal imaging, virtual polarizer triplets have two
distinct advantages over the Stokes representation. First, the virtual triplet
is a symmetric representation: the three brightness channels all have similar
properties and none is treated preferentially, unlike the Stokes and $pB$
systems in which Stokes $I$ (or $B$) is treated differently from the other
parameters. Second, each parameter in the virtual triplet representation is
positive-definite; this is important in the context of coronal measurements,
because it improves background subtraction when the background itself is
polarized.

Background subtraction for both coronagraphs and heliospheric imagers requires
accumulating and modeling a “minimum background” with varying degrees of
sophistication depending on field of view (e.g., Brueckner et al., 1995;
Howard et al., 2008). For measurements close to the Sun, the $pB^{\prime}$
component is negligible and the background is nearly unpolarized (e.g. Wlérick
& Axtell, 1957). At solar elongations of a few degrees, the polarization of
the F corona becomes significant (Weinberg & Hahn, 1980) and its polarization
must be treated differently. At solar elongations of 10∘ or more, the
starfield itself is an important source of background light and must be
eliminated independently (DeForest et al., 2011). Adding insult to injury, the
starfield is itself slightly linearly polarized at up to the few-percent level
by the interstellar medium (e.g., Mignani et al., 2019). The positive definite
nature of virtual polarizer images allows the use of existing background
estimation methods, which apply to unpolarized radiance, to be used on each of
the three polarizer channels independently. This allows the use of existing
techniques for removal of these background sources and instrument stray light
independently in each polarizer channel, to preserve the polarization as well
as the overall brightness of the background signal (DeForest et al., 2017).
Each virtual polarizer channel is itself “just” a radiance channel containing
a fixed subset of the overall brightness observed by the instrument; and
existing techniques for finding a running-minimum brightness apply.

Because the starfield and the F corona are fixed in different coordinate
systems, removing them requires operating not only in different pixel
coordinates (DeForest et al., 2011), but also in different polarization
coordinates. In particular, a linear polarization signal from a particular
star will maintain a fixed direction in celestial coordinates but have
variable direction of polarization in solar observing coordinates. Therefore,
coronal polarization values accumulated in an instrument coordinate system
must be transformed into solar coordinates to identify and remove the F
corona; then into celestial coordinates to identify and remove the starfield;
then back into solar coordinates for further analysis. This echoes the pixel
resampling that is also necessary to eliminate those noise sources.

In Section 3 we derived how to convert brightness values measured through a
triplet of polarizers, into either Stokes parameters (Equations 9, 12 and 13)
or ($B$, $pB$, $pB^{\prime}$) (Equations 9, 7, and 10, respectively). The
inverse transform from either the Stokes or ($B,~{}pB,~{}pB^{\prime}$) system
to a virtual polarizer triplet at angle $\theta$ is given in Equation 11.
Here, we derive how to convert directly between $M$, $Z$, $P$ virtual triplet
representations with different values of $\theta$.

Considering a polarizer at angle $\phi$, and substituting Equations 7, 9, and
10 into Equation 11, gives:

$B_{\phi}=\frac{1}{2}\left\\{\frac{2}{3}\left[\sum_{i}B_{i}\right]+\cos\left[2\left(\phi-\alpha\right)\right]\left[\frac{4}{3}\sum_{i}\left(B_{i}C_{i}\right)\right]+\sin\left[2\left(\phi-\alpha\right)\right]\left[\frac{4}{3}\sum_{i}\left(B_{i}S_{i}\right)\right]\right\\}~{}.$
(42)

Gathering terms across the finite sums gives

$B_{\phi}=\frac{1}{3}\sum_{i}\left\\{B_{i}+2B_{i}~{}\left(\cos\left[2\left(\phi-\alpha\right)\right]C_{i}+\sin\left[2\left(\phi-\alpha\right)\right]S_{i}\right)\right\\}~{},$
(43)

which simplifies to

$B_{\phi}=\frac{1}{3}\sum_{i}\left\\{B_{i}\left(1+2\cos\left[2\left(\phi-\theta_{i}\right)\right]\right)\right\\}~{},$
(44)

and thus to

$B_{\phi}=\frac{1}{3}\sum_{i\in\\{M,Z,P\\}}B_{i}\left[4\left(\cos\left[\phi-\theta_{i}\right]\right)^{2}-1\right]~{}.$
(45)

Equation 45 yields the predicted brightness through a polarizer at angle
$\phi$, given the brightness through an ($M,~{}Z,~{}P$) triplet of three
linear polarizers at angle $\theta$, i.e. a set of three linear polarizers at
$(\theta-\pi/3)$, $(\theta)$, and $(\theta+\pi/3)$.

Important properties of Equation 45 include: (a) $B_{\phi}$ is by construction
nonnegative for physical values of the $B_{i}$’s; (b) setting
$\phi=\theta_{i}$ for any of the $i$’s recovers the relation $B_{\phi}=B_{i}$;
and (c) the equation is numerically stable, with no potential poles or
singularities.

Applying Equation 45 to create predicted brightnesses for an ($M,~{}Z,~{}P$)
triplet of polarizers at $(\phi-\pi/3$), $(\phi)$, and $(\phi+\pi/3)$
constitutes a change of basis, fully describing the polarization state not via
direct brightnesses transmitted through a triplet of physical polarizers, but
via predicted brightnesses through a triplet of ideal “virtual polarizers” at
an arbitrary angle $\phi$. Provided only that both $\theta$ and $\phi$ are
known, the transformation via Equation 45 preserves the full polarization
state recorded by the original triplet.

##  5 Representations of Polarization and of Color

There is a strong analogy between representations of color and representations
of linear polarization.

The normal human visual system has three separate color channels at long,
medium, and short wavelengths, conventionally labeled “red”, “green”, and
“blue”, comprising a 3-D space (e.g., Young, 1802; Maxwell, 1857; Feynman et
al., 1963). The components are conventionally labeled ($r,~{}g,~{}b$), with
each representing radiance of light with a particular spectral characteristic
(a “primary color” of the system, which forms a basis vector for the 3-space).
In practice, a large number of slightly different ($r,~{}g,~{}b$) systems
exist, using slightly different primary colors. All such spaces share the
property of being positive-definite, like the radiance values discussed in
Section 3; this forms a positive-definite “gamut” that encompasses a subset of
all possible colors (Figure 2a). Representing colors outside the gamut would
require negative brightnesses of the primary colors, which is not physically
possible. Because of peculiarities of the human eye, no three primary color
vectors can both have positive-definite spectra (which may be reproduced in
the laboratory) and also encompass all humanly perceptible colors. The CIE
1931 ($X,~{}Y,~{}Z$) color space (Smith & Guild, 1931) is a standard
($r,~{}g,~{}b$) system that uses non-physical primary colors with negative
spectral intensities at some wavelengths, to provide a positive-definite gamut
that encompasses all of human vision, at the cost of requiring conversion to a
physical ($r,~{}g,~{}b$) system before the color can be rendered for viewing.
(Note that the $Z$ in ($X,~{}Y,~{}Z$) is distinct from the $Z$ in
($M,~{}Z,~{}P$): the former represents a blue-like fictional primary color,
while the latter represents zero polarizer offset from a reference angle
$\theta$).

Figure 2: Systems for representing color parameterize a 3-space that is
analogous to the 3-space needed to represent linear polarization. (a)
($r,~{}g,~{}b$) uses conventional Cartesian coordinates to represent three
separate radiances in different primary colors. (b) ($H,~{}S,~{}L$) or
($H,~{}C,~{}L$) use conical or cylindrical coordinates, respectively, around a
“white line” balancing $r$, $g$, and $b$, with $L$ and $C$ having units of
radiance and $S\equiv C/L$. (c) ($Y,~{}Cb,~{}Cr$) uses Cartesian coordinates
rotated to align with the white line.

The ($r,~{}g,~{}b$) systems are far from the only representations of color
space. Two other large classes are equally important.

Hue systems are based on the Munsell (1912) color wheel, which parameterizes
the space as hue, “chroma”, and “luminance” ($H,~{}C,~{}L$). These systems are
based around a cylindrical projection of ($r,~{}g,~{}b$) space. Luminance is a
radiance of white light (distance from the origin along a defined “white line”
in colorimetric space, so that $r$, $g$, and $b$ have a given fixed ratio);
chroma is a radial distance in the plane perpendicular to the white line at
particular radiance (the “chroma plane”; and hue gives the angle at which the
radius is to be drawn in that plane (Figure 2b). A more common variant of the
($H,~{}C,~{}L$) system is ($H,~{}S,~{}L$); this system replaces $C$ (chroma)
with a relative value $S$ (“saturation”) defined via $S\equiv C/L$, and L
(luminance) with a closely related quantity “lightness”.

“Opponent” color systems work similarly to ($H,~{}C,~{}L$) systems: they
represent color with luminance (usually $Y$ from the ($X,Y,Z$) system) and 2-D
signed chroma values in the chroma plane. The ($Y,~{}Cb,~{}Cr$) system (Union,
2017) has an unsigned luminance signal indicating distance from the origin
along the white line, and two signed “chrominance” signals (with units of
radiance) that describe location in the perpendicular plane (Figure 2c). Here,
$Y$ represents luminance and is normally a linear combination of the
($r,~{}g,~{}b$) signals: $Y\equiv K_{r}\,r+K_{g}\,g+K_{b}\,b$, with
$\sum_{i}K_{i}^{2}=1$. For convenience in this cursory treatment, we ignore
this scaling, effectively setting $Y=L$ and $K_{i}=1/\sqrt{3}$ for all the
$K$s. The ($Y,Cb,Cr$) and related color-opponent systems are important both
for video coding and because they describe well the perceptual aspects of the
human visual system, as can be verified by direct visual experiment
(Churchland, 2005).

Converting from $(r,~{}g,~{}b)$ to $(Y,~{}Cb,~{}Cr)$ is a simple linear
transformation. With
$\hat{Y}\cdot\hat{r}=\hat{Y}\cdot\hat{g}=\hat{Y}\cdot\hat{b}$ by construction,
the realization that the chroma-plane-projected
($\hat{r}^{\prime},~{}\hat{g}^{\prime},~{}\hat{b}^{\prime}$) vectors are
separated by $\pi/3$, and direct evaluation of $\cos(\pi/3)$ and
$\sin(\pi/3)$, the coefficients may be written down by simple inspection:

$Y=\frac{2}{3}\left[r+g+b\right]~{},$ (46)
$Cb=\frac{1}{\sqrt{3}}\left[2b-r-g\right]~{},\textrm{and}$ (47)
$Cr=\frac{2}{\sqrt{3}}\left(r-g\right).$ (48)

The $Cb$ and $Cr$ formulae describe perpendicular vectors in the chroma plane.
The $Cb$ formula arises because of the symmetry of the triangle formed by the
projected ($\hat{r}^{\prime},~{}\hat{g}^{\prime},~{}\hat{b}^{\prime}$)
directions in the chroma plane, and the fact that $\hat{Cb}$ is parallel to
$\hat{b}^{\prime}$. The $Cr$ formula is simpler because $\hat{Cr}$ is parallel
to the $\hat{g}^{\prime}-\hat{r}^{\prime}$ line and perpendicular to
$\hat{b}^{\prime}$.

Equations 46, 47, and 48 echo Equations 9, 12, and 13 respectively,
establishing the analogy between ($r,~{}g,~{}b$) colors and polarizer triplet
brightnesses — and, similarly, between ($Y,Cb,Cr$) colors and the linear
Stokes parameters.

Moving to the $(H,~{}S,~{}L)$ system, measuring the hue angle relative to
$\hat{b}^{\prime}$ immediately gives:

$H=\arctan\left[Cr/Cb\right]~{},$ (49)

where the four-quadrant arctan is implied, and

$S=\frac{\sqrt{Cr^{2}+Cb^{2}}}{Y}~{}.$ (50)

Again, this echoes directly the structure of Equations 15 and 16, converting
from the Stokes formulation to the ($B,\theta,p$) system, with $Y$, $H$, and
$S$ corresponding to $B$, $2\theta^{\prime}$, and $p$ respectively (where
$\theta^{\prime}$ has the appropriate reference angle subtracted for either
the $(I,~{}Q,~{}U)$ or $(B,~{}pB,~{}pB^{\prime})$ system).

Two useful insights come immediately from the connection between polarization
and color spaces. The first is an understanding of why polarizer triplets are
convenient and numerically stable representations of polarization space: they
are simply 3-vectors in the space, forming an orthonormal basis. Polarizer
triplet (or virtual polarizer triplet) brightness parameters have a higher
degree of mutual symmetry than do Stokes parameters (which are themselves a
different orthonormal basis of the same space), but as orthornormal bases they
are easily analyzed, manipulated, and interconverted.

The second useful insight is the notion of a “physical gamut” – the range of
values that are both representable, and also not physically impossible, within
a given representation space. Because radiances of light cannot be negative,
and the spectra of the human eye’s photosensitive pigments (treated as very-
high-dimensional vectors) are not perpendicular, no physically realizable
primary colors can represent all of human color in a tricolor system. The
specific ($X,~{}Y,~{}Z$) system of ($r,~{}g,~{}b$) primary colors (Smith &
Guild, 1931) uses impossible (non-physical) spectra to achieve a gamut that
can represent any possible color of light, at a cost of requiring conversion
to a specific physical colorimetric system before color can be rendered in the
physical world.

A peculiarity of ($X,~{}Y,~{}Z$) is that it can also represent physically
impossible colors that cannot exist at all, even in principle, because they
cannot be rendered with any positive-definite spectrum of light. This is also
the case for ($M,~{}Z,~{}P$) systems for representing polarization: some of
the possible vectors yield values of $pB$ and $pB^{\prime}$ (or, equivalently,
Stokes $Q$ and $U$) that are outside the Stokes inequality $B^{2}\geq
pB^{2}+pB^{\prime 2}$ which determines physically possible polarization
states. In the ($X,~{}Y,~{}Z$) system, the physically realizable colors occupy
an irregularly shaped locus in the chroma plane, driven by the specific
absorption spectra in the human eye. In an ($M,~{}Z,~{}P$) system, the
physically realizable polarizations (that satisfy the Stokes inequality) form
a circle in the “polarization plane” described by ($pB,~{}pB^{\prime}$) or
($Q,~{}U$); the circle is a cross-section of the Poincaré sphere (Born & Wolf,
1999), neglecting the third Stokes perturbation parameter $V$.

Happily, the circle is inscribed in the triangle formed by the ($M,~{}Z,~{}P$)
gamut on the polarization plane in polarimetric space, as can be verified by
direct substitution into Equation 11. This triangle is exactly analogous to
the triangle in the chroma plane formed by the ($r,~{}g,~{}b$) gamut, which
may be seen in Figure 2, panels (b) and (c).

The presence of the “Poincaré circle” physical gamut inside the triangular
($M,Z,P$) gamut ensures that all valid polarization states may be represented
by any suitable ($M,Z,P$) triplet of polarizers (as in Equation 45), because
the circle maps to itself under rotations about the center point.

Many impossible polarizations (which do not satisfy the Stokes inequality and
are therefore outside the Poincaré circle) may be represented with a
particular polarizer triplet. Substituting such triplet values into Equation
45 may yield negative polarizer radiances in the new system. As a trivial
example, a point with $B_{\theta,M}=B_{\theta,P}=0$ and $B_{\theta,Z}=1$
yields negative values of $B_{\phi}$ for nearly all values of $\phi-\theta$,
because no physically realizable polarization state can produce that
combination of radiances through an ideal linear polarizer triplet.
Contrariwise, all physically realizable states exist within the Poincaré
circle, which is inscribed in all possible polarizer-triplet gamuts on the
polarization plane, and therefore those states all yield positive values from
Equation 45 regardless of the value of $\phi-\theta$.

In sum: the principal systems of representing linear polarization are
mathematically analogous to the principal systems of representing color. This
yields both intuition, via the connection to ($r,~{}g,~{}b$) color systems,
and insights into why the polarizer triplet system works well for representing
linear polarization of light.

##  6 Conclusions

Linear polarization is an important measurement of the solar corona and
heliosphere. In the solar case, the Stokes parameters (I,Q,U), which are
defined in the instrument frame of reference, are not the most convenient
measure of polarization; instead, because of the peculiarities of polarization
in the visible solar corona, a large body of past and current literature
treats the polarized brightness as a sort of “virtual Stokes parameter”
oriented relative to the image-plane solar radial direction rather than in a
particular instrument direction, with its Stokes complement $pB^{\prime}$
generally ignored. Building on that practical and implicit treatment, we have
presented the derivation of the Stokes-equivalent $B$, $pB$, and $pB^{\prime}$
parameters directly from a triplet of brightness values in three polarizers.
While the three-polarizer approach to measuring coronal polarization was
described (in Russian) by Fesenkov (1935) and more briefly (in English) by
Öhman (1947) and Billings (1966), and is now commonly used in many
instruments, the derivation has generally used an asymmetric transition
through the Stokes parameters and/or been given cursory treatment.

Instruments making use of a polarizer triplet are subject to particular
uncertainties in measurement that are peculiar to the polarizer triplet
system. A literature search found no complete analytic treatment of these in-
principle knowable noise sources. Accordingly, we have presented the effects
of three major error or noise sources in three-polarizer instruments, on
polarimetry derived using the formulae presented here. These error sources
are: per-channel photometric error; polarizer orientation error; and finite
polarizer extinction coefficient. The most surprising result is that polarizer
orientation error does not affect the $pB$ measurement in first approximation
order, in the inner few solar radii of the solar corona where total
polarization is in the tangential ($pB$) direction. Our noise results supply
indicative limits relating these major noise sources in polarization
measurements, to noise in the calculated quantities $B$, $pB$, and $pB’$; but
we did not consider non-Gaussian probability distributions of the input data,
nor noise sources with significant covariance across measurements. Further,
specific definitions of derived parameters can affect noise terms in subtle
ways for derived quantities that are not linearly related to the source
quantities. For example, even “well-behaved” zero-mean, Gaussian-distributed
photometric noise (such as we considered here) produces a net offset (non-zero
mean perturbation) in determination of the total polarized brightness,
${}^{\circ}pB\equiv(Q^{2}+U^{2})^{0.5}=(p)(B)$, due to geometrical effects in
parameter space (Inhester et al., 2021)

The three-polarizer representation of polarization is symmetric and readily
transformed to different orientations. Using “virtual polarizer triplets”, to
represent polarization data symmetrically in multiple frames, conveniently
maintains the positive-definite and symmetric properties of polarizer triplet
measurements without corresponding physical polarizers. The mathematical
transformation between different three-polarizer systems is readily evaluated
and numerically stable, making virtual polarizer triplet analysis a useful
approach for converting polarization data between different coordinate
systems.

Possible applications for a virtual polarizer triplet representation include
regularizing polarimetry from a coronagraph or heliospheric imager that orbits
the Sun or Earth, and subtracting multiple types of data background that are
polarized in different coordinate systems relative to the instrument. Three
such overlapping backgrounds are potentially-polarized instrument stray light
(fixed in the instrument frame), wide field F corona (fixed in the solar
frame), and the starfield (fixed in the celestial frame).

There is a direct mathematical analogy between representations of linear
polarization, and representations of visual color; this analogy may help guide
intuition for the analytic derivations presented here and potentially other
applications of polarization space. The ($M,~{}Z,~{}P$),
($B,~{}pB,~{}pB^{\prime}$) and Stokes, and ($B,~{}\theta,~{}p$) systems of
representing linear polarization are seen to be direct analogs of the ($r$,
$g$, $b$), ($Y$, $Cb$, $Cr$), and ($H$, $S$, $L$) systems for representing
color. With this understanding, it becomes clear why the triplet polarizer
system works well for representing and manipulating polarization values. The
triplet brightness values are immediately seen to form an orthonormal basis of
the polarization space; and physically allowable polarization states are seen
to be representable independent of reference angle, because of the circular
geometry of the Poincaré sphere.

The authors thank A. Caspi, M. Beasley, and C. Lowder for useful discussion
and review of the article. The analysis and discussion were improved by
insights from the anonymous referee. This work was funded through PUNCH, a
NASA Small Explorer mission, via NASA Contract No. 80GSFC18C0014.

## References

  * Altschuler & Perry (1972) Altschuler, M. D., & Perry, R. M. 1972, Sol. Phys., 23, 410, doi: 10.1007/BF00148104
  * Arago (1843) Arago, F. 1843, Ann Philos Discov Month Report Progr Sci Art, 1, 209
  * Billings (1966) Billings, D. E. 1966, A guide to the solar corona (Academic Press)
  * Born & Wolf (1999) Born, M., & Wolf, E. 1999, Principles of Optics, 7th Ed. (Cambridge University Press), 31–33
  * Brueckner et al. (1995) Brueckner, G. E., Howard, R. A., Koomen, M. J., et al. 1995, Sol. Phys., 162, 357, doi: 10.1007/BF00733434
  * Churchland (2005) Churchland, P. 2005, Philosophical Psychology, 18, 527
  * Crifo et al. (1983) Crifo, F., Picat, J. P., & Cailloux, M. 1983, Sol. Phys., 83, 143, doi: 10.1007/BF00148249
  * Dai et al. (2014) Dai, X., Wang, H., Huang, X., Du, Z., & He, H. 2014, ApJ, 780, 141, doi: 10.1088/0004-637X/780/2/141
  * de Koning & Pizzo (2011) de Koning, C. A., & Pizzo, V. J. 2011, Space Weather, 9, 03001, doi: 10.1029/2010SW000595
  * DeForest et al. (2017) DeForest, C. E., de Koning, C. A., & Elliott, H. A. 2017, ApJ, 850, 130, doi: 10.3847/1538-4357/aa94ca
  * DeForest et al. (2011) DeForest, C. E., Howard, T. A., & Tappin, S. J. 2011, ApJ, 738, 103, doi: 10.1088/0004-637X/738/1/103
  * DeForest et al. (2013) —. 2013, ApJ, 765, 44, doi: 10.1088/0004-637X/765/1/44
  * Dere et al. (2005) Dere, K. P., Wang, D., & Howard, R. 2005, ApJ, 620, L119, doi: 10.1086/428834
  * Doyle (1892) Doyle, A. C. 1892, The Strand Mag., 4, 645
  * Fesenkov (1935) Fesenkov, V. G. 1935, Russ. Astron. J., XII, 310
  * Feynman et al. (1963) Feynman, R. P., Leighton, R. B., & Sands, M. 1963, The Feynman Lectures on Physics, Vol. I (California Institute of Technology), 35–1—35–10
  * Filippov et al. (1994) Filippov, B. P., Molodensky, M. M., & Koutchmy, S. 1994, in IAU Colloq. 144: Solar Coronal Structures, ed. V. Rusin, P. Heinzel, & J.-C. Vial, 601–604
  * Frazin et al. (2012) Frazin, R. A., Vásquez, A. M., Thompson, W. T., et al. 2012, Sol. Phys., 280, 273, doi: 10.1007/s11207-012-0028-3
  * Guhathakurta et al. (1999) Guhathakurta, M., Fludra, A., Gibson, S. E., Biesecker, D., & Fisher, R. 1999, J. Geophys. Res., 104, 9801, doi: 10.1029/1998JA900082
  * Hayes et al. (2001) Hayes, A. P., Vourlidas, A., & Howard, R. A. 2001, ApJ, 548, 1081, doi: 10.1086/319029
  * Hecht & Zajac (1974) Hecht, E., & Zajac, A. 1974, Optics (Addison-Wesley)
  * Heiles (2000) Heiles, C. 2000, 119, 923, doi: 10.1086/301236
  * Howard et al. (2008) Howard, R. A., Moses, J. D., Vourlidas, A., et al. 2008, Space Sci. Rev., 136, 67, doi: 10.1007/s11214-008-9341-4
  * Howard & Tappin (2009) Howard, T. A., & Tappin, S. J. 2009, Space Sci. Rev., 147, 31, doi: 10.1007/s11214-009-9542-5
  * Inhester et al. (2021) Inhester, B., Mierla, M., Shestov, S., & Zhukov, A. N. 2021, Sol. Phys., 296, 72, doi: 10.1007/s11207-021-01815-3
  * Kohl et al. (1995) Kohl, J. L., Esser, R., Gardner, L. D., et al. 1995, Sol. Phys., 162, 313, doi: 10.1007/BF00733433
  * Koomen et al. (1975) Koomen, M. J., Detwiler, C. R., Brueckner, G. E., Cooper, H. W., & Tousey, R. 1975, Applied Optics, 14, 743, doi: 10.1364/AO.14.000743
  * Leinert et al. (1981) Leinert, C., Richter, I., Pitz, E., & Planck, B. 1981, A&A, 103, 177
  * Lyot (1930) Lyot, B. 1930, Bulletin Astronomique, 6, 305
  * Maxwell (1857) Maxwell, J. C. 1857, Trans. R. Soc. Edinburgh, 21, 275
  * Mignani et al. (2019) Mignani, R., Shearer, A., Słowikowska, A., & Zane, S. 2019, Astronomical Polarisation from the Infrared to Gamma Rays (Springer)
  * Minnaert (1930) Minnaert, M. 1930, ZAp, 1, 209
  * Moran & Davila (2004) Moran, T. G., & Davila, J. M. 2004, Science, 305, 66, doi: 10.1126/science.1098937
  * Munro & Jackson (1977) Munro, R. H., & Jackson, B. V. 1977, ApJ, 213, 874, doi: 10.1086/155220
  * Munsell (1912) Munsell, A. H. 1912, Am. J. Psych., 23, 236
  * Newkirk et al. (1970) Newkirk, Gordon, J., Dupree, R. G., & Schmahl, E. J. 1970, Sol. Phys., 15, 15, doi: 10.1007/BF00149469
  * Ofman et al. (2000) Ofman, L., Romoli, M., Poletto, G., Noci, G., & Kohl, J. L. 2000, ApJ, 529, 592, doi: 10.1086/308252
  * Öhman (1947) Öhman, Y. 1947, Stockholms Observatoriums Annaler, 15, 2.1
  * Poland & Munro (1976) Poland, A. I., & Munro, R. H. 1976, ApJ, 209, 927, doi: 10.1086/154791
  * Reginald et al. (2017) Reginald, N. L., Gopalswamy, N., Yashiro, S., Gong, Q., & Guhathakurta, M. 2017, Journal of Astronomical Telescopes, Instruments, and Systems, 3, 014001, doi: 10.1117/1.JATIS.3.1.014001
  * Saito (1965) Saito, K. 1965, PASJ, 17, 1
  * Smith & Guild (1931) Smith, T., & Guild, J. 1931, Transactions of the Optical Society, 33, 73, doi: 10.1088/1475-4878/33/3/301
  * Union (2017) Union, I. T. 2017, Recommendation ITU-T H.273 “Coding-independent code points for video signal type identification” (International Telecommunications Union)
  * Vorobiev & Ninkov (2014) Vorobiev, D., & Ninkov, Z. 2014, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 9143, Space Telescopes and Instrumentation 2014: Optical, Infrared, and Millimeter Wave, ed. J. Oschmann, Jacobus M., M. Clampin, G. G. Fazio, & H. A. MacEwen, 91435F, doi: 10.1117/12.2059201
  * Weinberg & Hahn (1980) Weinberg, J. L., & Hahn, R. C. 1980, in Solid Particles in the Solar System, Vol. 90, 19–22
  * Wlérick & Axtell (1957) Wlérick, G., & Axtell, J. 1957, ApJ, 126, 253, doi: 10.1086/146396
  * Young (1911) Young, R. K. 1911, Lick Observatory Bulletin, 205, 166, doi: 10.5479/ADS/bib/1911LicOB.6.166Y
  * Young (1802) Young, T. 1802, Phil. Trans. R. Soc. London, 92, 12

##  7 Appendix: definitions of pB

Historically, the fact that the corona is tangentially polarized has led to
ambiguity in the term $pB$, which depending on context is used to refer either
to total polarization or to only a certain component of polarization. In the
former case, authors define and use a quantity like
${}^{\circ}pB\equiv(Q^{2}+U^{2})^{1/2}$ (where $Q$ and $U$ are the relevant
Stokes parameters and Stokes $V$ is neglected). ${}^{\circ}pB$ is thus the
total polarized brightness regardless of direction, and is analogous to $C$ in
the ($H,C,L$) color system in Section 5. In the latter case, authors define
and use a quantity like ${}^{\perp}pB\equiv B_{T}-B_{R}$, i.e. the difference
between the radiance observed through a tangentially-aligned and a radially-
aligned linear polarizer. ${}^{\perp}pB$ is effectively a variant on the
Stokes parameter $Q$, oriented radially in a solar image plane (Figures 1a and
1c).

The distinction between ${}^{\circ}pB$ and ${}^{\perp}pB$ is nearly moot in
lower coronal studies, where observed light happens to be polarized
perpendicular to the focal-plane solar radial direction, though ${}^{\circ}pB$
responds counter-intuitively to photometric noise (Inhester et al., 2021); but
the two quantities generalize differently in cases where the overall
polarization doesn’t happen to be tangential – for example, in wide field
coronagraphs or polarizing heliospheric imagers, where non-coronal polarized
light sources are significant sources of background radiance. In addition to
instrumental stray light, two such sources are the zodiacal light (Leinert et
al., 1981), and the starfield itself (Heiles, 2000).

Use of the ${}^{\circ}pB$ formalism is widespread in the coronal literature.
Examples include Öhman (1947) and Billings (1966), who define an $I_{p}$
similarly to our ${}^{\circ}pB$. Munro & Jackson (1977) used ${}^{\circ}pB$ in
their analysis of the corona, apparently to avoid having to deal with angle
explicitly (based on their Equation 3, which treats Thomson scattering). More
recently, the SOHO/UVCS (Kohl et al., 1995), SOHO/LASCO (Brueckner et al.,
1995), and STEREO/SECCHI (Howard et al., 2008) analysis pipelines calculated
${}^{\circ}pB$ explicitly (e.g., Ofman et al., 2000; Dere et al., 2005); and
sources as recent as Vorobiev & Ninkov (2014) and Reginald et al. (2017) use
the formalism in describing new polarimetric techniques using multiplexed
single-frame polarimetric detectors.

The ${}^{\perp}pB$ formalism is equally widespread through the coronal
literature, and seems to extend even earlier than ${}^{\circ}pB$. Minnaert
(1930), in his seminal work on coronal polarimetry, follows Young (1911) in
defining the quantity $p\equiv(J_{t}-J_{r})/(J_{t}+J_{r}+2A)$, in which $p$ is
relative polarization, $J_{t}$ and $J_{r}$ are radially and tangentially
polarized elements of the visual corona, and $A$ is an unpolarized background
brightness. We can immediately recognize that the quantity $(J_{t}+J_{r}+2A)$
is more commonly abbreviated $B$, and that Minnaert’s $p(J_{t}+J_{r}+2A)$ is
an expression for ${}^{\perp}pB$. The HAO K coronameter (Altschuler & Perry,
1972) was specifically built to exploit the tangential direction of coronal
polarization, and ${}^{\perp}pB$ appears explicitly (as $pB$) in their
analysis. Additional important references using ${}^{\perp}pB$ span six
decades and include Saito (1965), Koomen et al. (1975), Crifo et al. (1983),
Hayes et al. (2001), Moran & Davila (2004), Howard & Tappin (2009), de Koning
& Pizzo (2011), DeForest et al. (2013), Dai et al. (2014), and DeForest et al.
(2017). All of these sources define a ${}^{\perp}pB$ and refer to it as
“polarized brightness”; most refer to the quantity as $pB$ although some use
$I_{p}$ instead (note that $I_{p}$ is used in these sources as a synonym for
${}^{\perp}pB$ and in other sources, mentioned above, as a synonym for
${}^{\circ}pB$).

The distinction between ${}^{\circ}pB$ and ${}^{\perp}pB$ has been so blurred,
that some authors might even use both in the same paper, without
acknowledgement of the distinction. For example, both Guhathakurta et al.
(1999) and Frazin et al. (2012) compare $pB$ values across different
instruments, some of which generate ${}^{\circ}pB$ data products and some of
which generate ${}^{\perp}pB$ data products.

The polarized brightness parameter $pB$ is well entrenched in the literature,
and the two meanings ${}^{\circ}pB$ and ${}^{\perp}pB$ are close enough to not
differ greatly in the inner corona. However, in wide-field coronagraphs and
polarizing heliospheric imagers, other polarized light sources come into play.
In these contexts the two parameters differ greatly in meaning; past muddiness
requires more clarity from present and future authors, who ought to define
$pB$ explicitly (as many do already) wherever it is used.