1. Introduction

The 1940s mark to some degree a turning point not only in electrodynamics, but also in particle physics, with the creation of two distinct matrix calculi for the modeling of polarization phenomena [1]. The “Jones calculus” and “Stokes–Mueller calculus”—named after some of their inventors—were met with an enthusiastic response, as they overcome the “cumbersome” character of the more traditional algebraic and trigonometric methods used at the time. The analysis of both calculi rapidly became the subject of various publications, with the contemporaries of Clark Jones and Hans Mueller meticulously comparing the phenomenological Stokes–Mueller approach with the more theoretical Jones one [2–7].

Arguably, half a century later, these matrix representations have gained a renewed momentum with the development of “polarized Monte Carlo (MC)” programs for the modeling of polarized light propagation. Such programs, which, unlike scalar MC programs, acknowledge the vectorial nature of light, gained popularity in the 1990s as part of the ongoing research in astrophysics and atmospheric physics [8–10], before being reappropriated by various fields, including tissue optics [11–15].

These programs stochastically generate a large number of individual photon paths, where each path is constructed by means of successively sampled scattering events. In between two scattering events, the light propagation is represented by a photon path segment that is defined by a certain direction and a certain polarization state. The choice of that state’s representation is certainly a matter of some importance, and it is interesting to notice that most groups opt for the Stokes vector—hence dealing with intensities—whereas only a few authors prefer to use the Jones one—involving amplitudes [16,17]. The popularity of the Stokes formalism could mainly be explained by the recurrent claim which assumes that the Jones formalism does not allow for the treatment of depolarizing effects [12,18]. While the latter claim holds true for analytical treatments of the light propagation, it becomes irrelevant in the framework of stochastic and numerical simulations of the MC type. Thus, it is apparent, within the modern context of MC simulations, that the Jones and the Stokes–Mueller formalisms both deserve an appropriate reinterpretation. In this contribution, we lay the groundwork of our comparative study by discussing from a theoretical standpoint the significance and characteristics of both matrix methods when modeling light propagation in dielectric materials. In a second paper, we present the results yielded by two distinct MC programs (one of Jones-type, the other of Stokes–Mueller-type), where polarimetric imaging experiments on various suspensions of polystyrene spheres were simulated [19].

This paper is organized as follows: in Section 2, we begin by briefly revisiting the model of a polarized light beam, as it occurs in a MC simulation, and recall the Jones vector, for it being a representation of the polarization state that is close to the underlying physics. In Section 3, we focus on the statistical aspects of polarization measurement. (For a comprehensive treatment of the topic, see [20].) We introduce the Stokes vector from a “measurable quantity” viewpoint, recalling how the Stokes vector of the light backscattered or transmitted from a sample can be constructed from a series of measurements. The subjects of Section 4 are the Jones and Mueller matrices operating on the Jones and Stokes vectors, respectively. We interpret the Mueller matrix in terms of observables and discuss the role of the fluctuations in the observed system. In Section 5, we detail the modeling of the polarization analysis involved in polarimetric imaging experiments, since this key issue is, unfortunately, scarcely reported in the literature. In Section 6, we discuss the manifestation of fluctuations/“depolarizing effects” within both types of MC simulations.

2. Polarized Light

The term “polarization” merely expresses the fact that the electromagnetic field is a vector field. Here, we are mainly interested in dielectrics, and therefore in the electric field vector $\mathbf{E}$. Moreover, we are concerned with light fields, i.e., with radiative transfer of energy (and momentum) in transversal electromagnetic waves. A simple yet realistic example of a light field is a quasi-monochromatic quasi-plane beam wave is given by

Beam waves are useful for wave-optical modeling of photon path segments sampled in a MC simulation [21]: each segment corresponds to a beam wave with a propagation vector $\mathbf{k}$ and transversal polarization vector $\widehat{\mathbf{e}}$. The total light field can be thought of as a superposition of such beam waves. However, the beam profile $X(\mathbf{r})$ turns out to be of little importance in a MC simulation (Section 7.5, again in Ref. [21]). Therefore, one can conceptually reduce the beam to a “polarized geometrical ray,” providing that one chooses to neglect interference effects.

Polarization of a transversal beam wave with a known $\mathbf{k}$ is completely specified by two complex numbers, which can be arranged into the so-called Jones vector

Dirac’s bra-ket notation $|E{〉}_{\mathbf{k}}$ emphasizes the relationship of the Jones vector with state vectors of quantum theory [22]. The quantities $E_x=|E_x|\exp i{\varphi }_x$ and $E_y=|E_y|\exp i{\varphi }_y$ do not necessarily mean local field components at some position in a light wave. In general, they should be understood as state amplitudes, which fully characterize the mode of the light field (e.g., a beam wave) that is selected for detection by an optical receiver. (A special case would be a point-like receiver, e.g., a single molecule, which would indeed probe the local quantities.) The norm $〈E|E〉={|E_x|}^2+{|E_y|}^2=P$ represents the received power, rather than the “intensity,” as often denoted with $I$.

Any polarization state can be expressed as a linear combination $|E〉=E_a|e_a〉+E_b|e_b〉$ of two orthonormal states $|e_a〉$ and $|e_b〉$ such that $〈e_a|e_a〉=〈e_b|e_b〉=1$ and $〈e_b|e_a〉=0$. Three particularly important pairs are indicated in Table 1.

Table 1. Fundamental Pairs of Polarization States

View Table

3. Measuring Polarization: Stokes Vector

In order to understand the differences between the “Jones” and the “Stokes” approaches to polarization treatment, it is important to bear in mind that neither the amplitudes $E_x(t)$, $E_y(t)$, nor the power $P(t)$ can be measured directly. In optics, the only measurable quantity is the energy, i.e., the number of photons ${\int }_{T}P(t)\mathrm{d}t$ accumulated in a detector of finite area $A$ during a finite integration time $T$. Subsequently, one can divide by $T$ to obtain the mean power $\overline{P}=〈P{〉}_T={〈〈E|E〉〉}_T$, and divide by the detector area $A$ to get the mean irradiance $\overline{I}=\overline{P}/A={〈〈E|E〉〉}_{TA}$. Thus, measuring polarization consists in measuring the mean power $\overline{P}$ or the mean irradiance $\overline{I}$ transmitted through polarization filters, “smearing” thereby any kind of fluctuations. The fluctuations can occur in time, or in space, across the area of a detector exposed to many overlapping partial beam waves whose interference gives rise to the so-called speckles. Now, we mainly focus on fluctuations in time, assuming a “single-mode” detector that selects a single partial beam, or a detector whose size is much smaller than the mean size of a speckle, i.e., smaller than the area of coherence. The generalization including spatial fluctuations should be obvious.

Two types of filters are needed to determine the polarization of a beam wave: a linear polarizer (often a Polaroid foil) with variable orientation and a pair of circular polarization filters realized by combining a linear polarizer with a wave plate [26–28]. Transmission of the state amplitudes through a filter is expressed as $|E_{\mathrm{tr}}〉=\mathbf{F}|E_{\mathrm{in}}〉$, where $\mathbf{F}=|f〉〈f|$ is the projection operator and $|f〉$ is the Jones vector characterizing the filter. This $|f〉$ is normalized such that $〈f|f〉=1$. Power transmitted through the filter is $P_{\mathrm{tr}}={|〈f|E_{\mathrm{in}}〉|}^2$. A long time ago, Sir George Gabriel Stokes realized that all the information on polarization that is accessible from such measurements of $P_{\mathrm{tr}}$ can be summarized in four parameters [30]:

The physical significance of this formal definition is depicted in Fig. 1, where, by the same occasion, the conventions used in the present work are specified. The Stokes parameters can be determined by measuring the mean power $\overline{P}$, or the mean intensity $\overline{I}$, transmitted through three fundamental pairs of orthogonal filters: $I$ stands for the beam’s observed “total intensity,” $Q$ gives the prevalence of $|LX〉$ polarization over $|LY〉$ polarization, $U$ gives the prevalence of $|L+〉$ polarization over $|L-〉$ polarization, and $V$ gives the prevalence of $|C+〉$ polarization over $|C-〉$ polarization. (Note that the six measurements suggested in Fig. 1 are not linearly independent: four independent measurements would be sufficient to determine the Stokes parameters.)

 

Fig. 1. Graphical representation of the four Stokes parameters as determined from measurements with the three fundamental filter pairs defined in Table 1. The arrow $\mathbf{e}(t)$ drawn in the $x–y$ plane of the transversal coordinate system represents the real part of the complex polarization vector of the light transmitted through a filter: $\mathbf{e}(t,z)=\mathfrak{R}[{\widehat{\mathbf{e}}}_f\exp (-i\omega t+ikz)]$ at $z=0$. The pictograms defining $I$, $Q$, $U$, and $V$ represent the motion of the tip of this arrow in the $x–y$ plane, as seen when looking from the point of view of the receiver, against the beam propagation. To avoid confusion that arises from the “point of view” definition, we define right-circular polarization $|C-〉$ through the chirality: the loci of the tips of $\mathbf{e}(t,z)$, at any instant $t$, form a right-handed spiral that is pushed (not rotated) through the $x–y$ plane.

Download Full Size | PDF

Chandrasekhar [29] assembled those four parameters into a vector $\mathrm{S⃗}=(I,Q,U,V)$, which is nowadays known as the Stokes vector, and which is the central quantity in a “Stokes–Mueller” type of MC simulation. However, working with Stokes vectors requires caution: a Stokes vector is a rather abstract quantity. Like a Jones vector, it is defined only in the chosen transverse coordinate system, but unlike a Jones vector, it cannot be directly associated with the components of an ordinary vector in 3D space. Moreover, there is an unfortunate ambiguity of the definition, which causes confusion [34] affecting the sign of the elements $Q$, $U$, and mainly $V$. All this turns the comparing and debugging of Stokes–Mueller-type MC programs into a tedious and headache-causing exercise.

We are now in the position to appreciate the differences between “Stokes” and “Jones.” While a Jones vector is composed of two complex numbers, a Stokes vector comprises four real numbers. Thus, the amount of information contained in the two representations of the polarization state is the same. However, the type of information is different: a Jones vector provides the full information about an instantaneous polarization state, including the common phase $\varphi$, in which one can also include the optical oscillations $\omega t$. Therefore, Jones vectors can be used to describe interference effects through the “coherent” addition of amplitudes. On the other hand, an instantaneous Jones vector does not contain any information on the fluctuations of the amplitude and therefore, it “cannot describe unpolarized light.” Moreover, a Jones vector, representing a state amplitude, is not directly measurable.

The measurable quantity is the Stokes vector. In the measurement, the information on the phase $\varphi$ is lost, but what is gained is a part of the information on fluctuations that is contained in the quadratic averages $〈{|E_x|}^2{〉}_T$, $〈{|E_y|}^2{〉}_T$ and $〈E_x{E}_y^{*}{〉}_T=〈|E_x||E_y|\exp i({\varphi }_x-{\varphi }_y){〉}_T$. The components $Q$, $U$, and $V$ measure the degree of correlation, also called “degree of coherence,” of the fluctuating components of the polarization vector [35]. A light field is said to be “completely unpolarized” if the amplitude components are completely uncorrelated so that $〈{|E_x|}^2{〉}_T=〈{|E_y|}^2{〉}_T$ and $〈E_x{E}_y^{*}〉=〈E_x〉〈E_y^{*}〉=0$. In this case, $\mathrm{S⃗}=(I,0,0,0)$. On the other hand, in the absence of fluctuations, i.e., when the amplitude components remain fully correlated during the experimental integration time $T$, the averaging brackets become superfluous:

This set of equations can be solved for $|E_x|$, $|E_y|$ and ${\varphi }_x-{\varphi }_y$. Thus, the corresponding Stokes and Jones vectors are equivalent, except that the common phase $\varphi$ is lost in $\mathrm{S⃗}$: only the phase difference ${\varphi }_x-{\varphi }_y$ can be recovered from a Stokes vector. One easily verifies that in the absence of fluctuations $Q^2+U^2+V^2=I^2$, whereas in general $Q^2+U^2+V^2\le I^2$. Thus, the quantity $\mathrm{\Pi }=\sqrt{Q^2+U^2+V^2}/I$ can be used to quantify the degree of polarization, so that $\mathrm{\Pi }=0$ means completely unpolarized light and $\mathrm{\Pi }=1$ fully polarized light, a so-called “pure polarization state.”

It is important to realize that the measured degree of polarization depends on the length of the characteristic fluctuation time scale, i.e., on the coherence time ${\tau }_c$, with respect to the integration time $T_i$ of an individual measurement. If $T_i\ll {\tau }_c$, then one always measures completely polarized light, $\mathrm{\Pi }=1$. Suppose that we perform a large number of measurements, so that each is done with a short integration time $T_i\ll {\tau }_c$, but the whole ensemble covers a period $T\gg {\tau }_c$. Each of the measurements samples an instantaneous realization of the fluctuating polarization state, i.e., the individual members of the ensemble of Stokes vectors ${\mathrm{S⃗}}_i$ are randomly, yet fully polarized; $\mathrm{\Pi }=1$ for each of them. However, the sum $\mathrm{S⃗}=\sum _i{\mathrm{S⃗}}_i$ will yield partially polarized light with $\mathrm{\Pi }<1$, or in the extreme case, a completely unpolarized Stokes vector $\mathrm{S⃗}=(I,0,0,0)$, i.e., $\mathrm{\Pi }=0$.

We conclude by generalizing the findings: the so-called unpolarized or partially polarized light should be understood as an ensemble of “independent light streams” [30] with random polarization states. A Stokes vector captures the essential second-order averages, i.e., the correlations of the field amplitudes. Therefore, a Stokes vector can be regarded as being equivalent to such an ensemble. It should be noted that this picture of independent light streams suits particularly well the ensembles of photon paths that are generated within polarized MC simulations.

4. Jones Calculus and Mueller Calculus

In the 1940s, several researchers, inspired perhaps by the success of the matrix formulation of quantum mechanics, were compelled to adapt the matrix calculus concepts to classical optics. Clark Jones, researcher with Polaroid Corporation and Harvard scholar, presumably weary of the vector calculus in electrodynamics, was prompted to design an efficient formalism for the analysis of light propagation through optical elements such as, for example, sheets of Polaroid foil [37]. Jones decided to work on the amplitude level. In his formalism, an optical element is represented by a $2\times 2$ complex matrix operator $\mathbf{O}$, which we now call the Jones matrix. The Jones matrix acts on an input polarization state $|E^{\mathrm{in}}〉$ represented by a 2D complex column vector, namely, the Jones vector from Eq. (2), and produces an output state $|E^{\mathrm{out}}〉$ such that

The time argument in the elements of the Jones matrix in Eq. (11) is a modern addition: we anticipate fluctuations in the optical system under consideration. For example, the axis of a badly fastened analyzer may jitter because of a heavy truck passing below the lab window. Jones was primarily concerned with the propagation of a light beam through optical elements aligned in a straight line on an optical bench, where such disturbances are rare. Later, however, he expanded his calculus to include scattering in heterogeneous materials. Thereby, he realized an apparently serious deficiency of his calculus. In the fifth of a series of eight papers on the topic, Jones writes that his calculus “is completely unable to treat systems which serve to depolarize the light which passes through them.” Which is of course correct: a single instantaneous Jones matrix acting on a pure polarization state cannot produce anything else than a pure polarization state with $\mathrm{\Pi }=1$. In other words, a single Jones matrix does not depolarize. Jones then pursues by stating that the latter restriction “can be avoided altogether by the use of a more powerful calculus which Professor Hans Mueller has developed” [37].

Mueller (and others [38]) realized that in frequently encountered situations, one does not care for interference and therefore, propagation of polarized light can be treated in terms of the observables, i.e., Stokes vectors $\mathrm{S⃗}$. The transformation of the polarization state upon transmission through the optical system is effected by a $4\times 4$ real matrix operator, which we shall call Perrin–Mueller matrix (PM), denoted as $M$: ${\mathrm{S⃗}}_{\text{out}}=\mathbf{M}{\mathrm{S⃗}}_{\mathrm{in}}$. To highlight the role of internal fluctuations in the optical system under consideration, we introduce the PM matrix in two steps. First we note that any Jones matrix can be transformed into an instantaneous PM matrix. This can be done formally by expanding the Jones matrix in terms of Pauli matrices (see, e.g., Section 4.6 in [21]); final formulas can be found in any book on polarimetry/ellipsometry, e.g., [42]. Here, we instead define the PM matrix in terms of the observables, in the same spirit as the elements of the Stokes vector in Fig. 1:

Here, ${\mathrm{M⃗}}_1(t)$, ${\mathrm{M⃗}}_2(t)$, ${\mathrm{M⃗}}_3(t)$, and ${\mathrm{M⃗}}_4(t)$ denote the four columns of the PM matrix. The pairs ${\mathrm{S⃗}}_{LX}(t)$, ${\mathrm{S⃗}}_{LY}(t)$ and ${\mathrm{S⃗}}_{L+}(t)$, ${\mathrm{S⃗}}_{L-}(t)$ and ${\mathrm{S⃗}}_{C+}(t)$, ${\mathrm{S⃗}}_{C-}(t)$ are the instantaneous output Stokes vectors as would be measured when transmitting the $3\times 2$ fundamental pure polarization states defined in Table 1 through the optical system that is frozen in a certain instant of time $t$. (In the first column of the PM matrix, ${\mathrm{S⃗}}_{+}(t)$, ${\mathrm{S⃗}}_{-}(t)$ means any of the three pairs.) Notice that the PM matrix is actually a shorthand notation for the three pairs: all six Stokes vectors can be obtained from the four columns of the PM matrix.

In the second step, we explicitly average over the internal fluctuations during the measurement period $T$:

Having understood the structure of the PM matrix in terms of the observables, one can attempt the interpretation of its elements. The first column ${\mathrm{M⃗}}_1$ is easy: it represents the response of the optical system to unpolarized input, i.e., its capability to repolarize the impinging light stream (recall polarization of blue sky). An unpolarized light stream can be modeled as an incoherent superposition of two light streams with mutually orthogonal polarizations. Obviously, the first column is the Stokes vector ${\mathrm{M⃗}}_1=\mathbf{M}{\mathrm{s⃗}}_{\mathrm{in}}$ that results from the unpolarized input Stokes vector ${\mathrm{s⃗}}_{\mathrm{in}}=(1,0,0,0)$.

The remaining three columns reflect the symmetries of the system with respect to the rotation of impinging polarization. Their interpretation is more difficult, since symmetries are affected not only by birefringence and chirality, but also by the geometry of the light path through the sample. We recommend to investigate instead of the PM matrix the three pairs of fundamental Stokes vectors, which can be easily derived from it. We shall return to this point in the second part, in the context of polarization imaging.

5. Measuring Polarization: the Role of the Optical Receiver

Before proceeding with the comparison between MC simulations of Jones type and Stokes–Mueller type, we must clarify an important issue concerning the optical setup needed for polarization analysis. Polarimetry requires a careful optical design. One cannot just send the light stream “somehow” through a filter, because real polarizing elements, such as Polaroid foils and/or wave plates, are sensitive to the incidence angle of the incoming light stream. Therefore, collimation and nearly normal incidence are required. On the other hand, a high angular aperture may be required for high-resolution imaging polarimetry. Both requirements are satisfied by the optics of a telecentric system, as shown in the left part of Fig. 2: the objective collects the backscattered light through a high angular aperture and the tube lens refocuses the collimated beams on the imager. Polarization analysis is done in the collimated space.

 

Fig. 2. Left: optical setup for polarization imaging. Polarization analysis is done in the collimated space of the telecentric lens pair. Right: geometrical representation of polarization analysis of the light stream received by one pixel of the imager. The collimating objective lens is replaced by the refractive plane ${\sigma }_D$, which deflects a ray (dotted line) exiting from the observed pixel area on the sample surface in the direction ${\widehat{\mathbf{n}}}_D$. The polarization analyzer is represented by its filter state ${\widehat{\mathbf{e}}}_f$, here, for example, ${\widehat{\mathbf{e}}}_f=|L+〉$. Further explanations in text.

Download Full Size | PDF

Most papers on MC simulations of Stokes–Mueller type treat the detection issue rather generously, giving the impression that the Stokes vector, which is somehow attached to the detected photon, only needs to be “transformed” into the detector’s coordinate system. Unfortunately, there is no general consensus on the nature of this “transformation” [12,14,15,44], but this causes little concern, since all programs yield qualitatively similar results. However, for quantitative comparisons between simulation and experiments, or between two simulations, an accurate modeling of the detection procedure is paramount. Therefore, we present here our detection model in some detail. Thereby, we take the opportunity to point out the difference between an electric field vector $\mathbf{e}$ and a Jones vector $|e〉$ and identify a potential source of confusion in published MC simulations.

Our model of the detection process is sketched in the right-hand part of Fig. 2. Here, we outline only the ray optical approximation, that is, we regard the partial beams which represent segments of a photon path, as polarized rays. (The wave optical background is given in Section 7.5.6 of Ref. [21].) The ray representing a photon’s path emerges from the sample in the direction ${\widehat{\mathbf{k}}}_{\text{out}}$ with the corresponding electric field vector ${\widehat{\mathbf{e}}}_{\text{out}}$. The collimating lens (objective) is modeled as a refractive plane ${\sigma }_D$ whose normal ${\widehat{\mathbf{n}}}_D$ is aligned the detector’s optical axis. The refraction of the ray is a “deterministic” scattering event which turns the ray from the propagation direction ${\widehat{\mathbf{k}}}_{\text{out}}$ into the direction ${\widehat{\mathbf{k}}}_s={\widehat{\mathbf{n}}}_D$. The vectors ${\widehat{\mathbf{k}}}_{\text{out}}$ and ${\widehat{\mathbf{n}}}_D$ span the scattering plane, which coincides with a meridian plane, if ${\widehat{\mathbf{n}}}_D$ is aligned with the $z$ axis of the lab system. The normal to the scattering plane defines the $x$ axis of the meridian system: ${\widehat{\mathbf{x}}}_m={\widehat{\mathbf{k}}}_{\text{out}}\times {\widehat{\mathbf{n}}}_D/|{\widehat{\mathbf{k}}}_{\text{out}}\times {\widehat{\mathbf{n}}}_D|$. The polarization vector ${\widehat{\mathbf{e}}}_{\text{out}}$ can be decomposed into its perpendicular and parallel components: ${\widehat{\mathbf{e}}}_{\text{out}}={\mathbf{e}}_{\perp }+{\mathbf{e}}_{\parallel }$, whereby ${\mathbf{e}}_{\perp }\parallel {\widehat{\mathbf{e}}}_{xm}$ and ${\mathbf{e}}_{\parallel }$ is in the scattering plane. Refraction causes the parallel component ${\mathbf{e}}_{\parallel }$ to be rotated by an angle $\theta$ around ${\widehat{\mathbf{x}}}_m$ into the plane ${\sigma }_D$. The perpendicular component ${\mathbf{e}}_{\perp }$ remains unchanged. The polarization vector of the refracted ray is

Admittedly, the geometrical analysis above, including the projection concept, is somewhat unusual and therefore we translate the whole procedure into Jones formalism. We start with a photon whose Jones vector $〈e{|}_{\text{out}}=(e_x^{*},e_y^{*}{)}_{\text{out}}$ is given in a local transversal coordinate system ${\widehat{\mathbf{x}}}_{\text{out}}$, ${\widehat{\mathbf{y}}}_{\text{out}}$, ${\widehat{\mathbf{z}}}_{\text{out}}={\widehat{\mathbf{k}}}_{\text{out}}$. The state is normalized, i.e., $〈e|e〉=1$. The first step is to transform the Jones vector into the meridian coordinate system $m$ specified by ${\widehat{\mathbf{x}}}_m$, ${\widehat{\mathbf{y}}}_{mk}$, ${\widehat{\mathbf{z}}}_m={\widehat{\mathbf{k}}}_{\text{out}}$. The $m$ system is rotated with respect to the out system by an angle ${\varphi }_m$ around the common $z$ axis ${\widehat{\mathbf{k}}}_{\text{out}}$. Thus, $|e{〉}_m=\mathbf{R}({\varphi }_m)|e{〉}_{\text{out}}$, where $\mathbf{R}({\varphi }_m)$ is the corresponding 2D rotational matrix (in Fig. 2, $\varphi$ rotations are shown only schematically). Notice that we only transformed the coordinates into a new coordinate system. The polarization vector ${\widehat{\mathbf{e}}}_{\text{out}}$ remains unchanged; only the components of the Jones vector are transformed. Next we perform the refraction, regarding it as a deterministic scattering event represented by a scattering matrix $\mathbf{S}$: $|e{〉}_{mD}=\mathbf{S}|e{〉}_m$. Refraction rotates the $m$ system around ${\widehat{\mathbf{x}}}_m$ together with the ray’s field vector embedded therein. This is a true rotation, not only a coordinate transformation. Thus, ${\widehat{\mathbf{e}}}_{\text{out}}$ is rotated into ${\widehat{\mathbf{e}}}_D$ that lies in the refractive plane ${\sigma }_D$, but the Jones vector remains unchanged. Thus, the $\mathbf{S}$-matrix is trivial, namely, $\mathbf{S}=\left(\begin{array}{cc}1& 0\\ 0& 1\end{array}\right)$.

For polarization analysis we must express $|e{〉}_{mD}$ in the detector’s coordinate system. If the detector system is aligned with the lab system so that ${\widehat{\mathbf{n}}}_D={\widehat{\mathbf{z}}}_{\text{lab}}$, then the transformation is a rotation of $|e{〉}_{mD}$ by ${\varphi }_D$ around ${\widehat{\mathbf{n}}}_D$: $|e{〉}_D=\mathbf{R}({\varphi }_D)|e{〉}_{md}=\mathbf{R}({\varphi }_D)\mathbf{R}({\varphi }_m)|e{〉}_{\text{out}}$. (The polarization vector ${\widehat{\mathbf{e}}}_D$ itself remains unchanged.) Finally, polarization analysis is achieved by applying a filter operator $\mathbf{F}=|f〉〈f|$, where $|f〉$ is any one of the six polarization states listed in Table 1. Putting all the steps together, the state transmitted through the filter is

Here, ${\mathbf{R}}_M$ are the PM matrix versions of the 2D transformation matrices $\mathbf{R}$. The filter operation does not need to be carried out explicitly, since it is implicit in the definition of the Stokes vector. The detection scheme from Eq. (20) is implemented in the program by Ramella-Roman and co-workers [14,15]. On the other hand, from the paper of Jaillon and Saint-Jalmes [44], one deduces a slightly different version: ${\mathrm{S⃗}}_D={\mathbf{R}}_M({\varphi }_D){\mathbf{R}}_M(\theta ){\mathbf{R}}_M({\varphi }_m){\mathrm{S⃗}}_{\text{out}}$. There is a $\theta$-transformation matrix in place of the trivial unit scattering matrix. It can be shown that, in this way, one projects ${\widehat{\mathbf{e}}}_{\text{out}}$ on the detector plane ${\sigma }_D$, instead of rotating into it by refraction. This scheme is equivalent to neglecting the second term in Eq. (18), which gives the transmission probability $P_f={|{\widehat{\mathbf{e}}}_{\text{out}}·{\widehat{\mathbf{e}}}_f^{*}|}^2$. Clearly, the neglected term becomes small when the aperture of observation is restricted to below 0.1. The scheme of Jaillon et al. is actually a low aperture approximation [45].

6. Jones and Stokes–Mueller in MC Simulations

In a multiply scattering sample, the ensemble of PM matrices is “created” in two ways, based on two kinds of randomness. The first kind is the randomness of the quantum world. Loosely speaking, each photon is sent to the detector through a randomly selected path, which converts a pure input state into an output state with a full yet random polarization. After the incoherent addition of the individual random outcomes, one ends up with a partially polarized Stokes vector $\mathrm{\Pi }<1$. Precisely this kind of ensemble averaging is implicitly done in a MC simulation, be it of “Jones” or “Stokes–Mueller” type.

The second kind of randomness results from thermal fluctuations. Thermal fluctuations are certainly present in any practically relevant sample (tissues or tissue phantoms), but rarely affect the polarization state (Section 7.3 in [21]). A cause for depolarization would be scattering from dielectrically anisotropic objects, e.g., collagen needles, whose orientation exhibits thermal fluctuations. An explicit expression for the averaged PM matrix of anisotropic dipole scatterers can be found in [48]. Here, one may recognize a potential advantage of the Stokes–Mueller approach in MC simulations: computation time may be saved by using a preaveraged PM matrix. On the other hand, in most cases the averaging would anyway have to be done numerically and such numerical averaging can directly be implemented in the MC algorithm. To the best of our knowledge, no MC simulation with depolarizing scatterers has been done so far, and therefore, we focus on the first kind of randomness, namely, the sampling of random paths.

Consider a chosen photon path that encompasses $n$ scattering events. Mathematically, this path is given by a chain of Jones matrices:

 

Fig. 3. Geometry of a scattering event. Notice here our convention: the $x$ axis of the $i$th coordinate system is perpendicular to the scattering plane ${\sigma }_i$.

Download Full Size | PDF

We assume that the input state is normalized, i.e., $〈e_{\mathrm{in}}|e_{\mathrm{in}}〉=1$. The probability to sample this particular chosen path is $P_{\mathrm{path}}=〈E_{\mathrm{out}}|E_{\mathrm{out}}〉=|E_x^{\text{out}}{|}^2+|E_y^{\text{out}}{|}^2$. In other words, in a MC simulation, this particular path will be sampled with an occurrence frequency that is proportional to $P_{\mathrm{path}}$.

All this can be translated one-to-one into the Stokes–Mueller language, i.e., each of the matrices can be transformed into the corresponding PM matrices that operate on Stokes vectors:

To simulate unpolarized light [25,49] we apply the so-called “mixed sampling”: each photon path’s input state is sampled over a random ensemble covering all possible states. An alternative to this would be the “orthogonal sampling” based on Stokes’ equivalence theorem: we propagate two orthogonal polarization states, for example, $〈e_X|=(1,0)$ and $〈e_Y|=(0,1)$, which we chose with equal probability 1/2. The frequency of occurrence of our particular path in the unpolarized simulation is $P_O=(P_X+P_Y)/2$. The same of course applies with the Stokes–Mueller formalism, when sampling the paths with the two orthogonal inputs ${\mathrm{s⃗}}_X=(1,1,0,0)$ and ${\mathrm{s⃗}}_Y=(1,-1,0,0)$, namely, $P_O=(I_{\text{out}}^X+I_{\text{out}}^Y)/2$.

Since the Stokes–Mueller formalism is incoherent, i.e., linear in intensity, we can add any two measurements along the chosen path, for example, the two measurements with the orthogonal polarized inputs ${\mathrm{s⃗}}_X$ and ${\mathrm{s⃗}}_Y$:

Polarization sampled Jones MC gives the same mean as a Stokes–Mueller MC with unpolarized input, but one can expect subtle statistical differences that concern second-order intensity statistics, i.e., “fluctuations” of the intensity and related quantities that occur in averages of squared Stokes vectors (e.g., $〈E_x{E}_x^{*}〉$ or $〈E_x{E}_y^{*}〉$). It is possible that Stokes–Mueller MC carried out with unpolarized input ${\mathrm{s⃗}}_O$ improves the convergence of a MC simulation in a similar fashion as the “photon packet” MC greatly improves the statistics of absorption [50]. On the other hand, natural unpolarized light is certainly best modeled with the mixed sampling. Certainly, one does not get the same unpolarized stream when alternating two orthogonal states and when varying the polarization in a similar way as we do it with mixed sampling. But the differences do not occur in the means, only in the fluctuations.

7. Conclusion

In this paper, we have presented a concise description of the physical concepts that govern the propagation of polarized light in MC simulations. Our findings suggest the following: (i) The modeling of polarized light detection plays a decisive role in the simulation of polarimetric experiments, and should, as such, not be overlooked. Here, we have chosen to base our model on the geometrical approximation, for it being closer to the various approaches reported in the literature. (ii) Both “Jones” and “Stokes–Mueller” are equally suitable to model the propagation of not only pure states but also partially polarized/totally unpolarized states. This is because MC simulations involve ensemble sampling. As long as interference effects are not of interest, both formalisms are expected to be perfectly equivalent in the first-order statistics. Second-order statistical differences may be expected in the case of unpolarized light.

Admittedly, each formalism has its advantages and limitations. The Jones vector provides a more intuitive description of the light polarization, and further allows for the simulation of interferences, thanks to the phase information it carries. On the other hand, the preaveraged Stokes–Mueller notation can be beneficial when there are fluctuations in a sample, as faster convergence can be achieved. However, as already mentioned, the Stokes vector constitutes a rather abstract quantity: consequently, the implementation of a MC program of Stokes–Mueller type might be more prone to errors.

We do not argue that a polarized MC program should be restricted to either Jones or Stokes–Mueller. Both methods can be combined, or novel models can be developed to satisfy individual purposes [51]. Yet, it is doubtful that Jones and Stokes have been exploited to the maximum of their capacities. The simulation programs we have at our disposal supply a fertile terrain of investigation, provided that they are reliable: this constitutes another topic of concern, which we treat in the second paper dedicated to our comparative study [19], where we quantitatively compare two independent polarized MC programs. The latter comparison also serves as an empirical support to the theoretical considerations presented in this first paper.