Polarimeters characterize the polarization properties of materials. They find application in, for instance, optical samples [2] , cancer non-invasive screening tools [3] in clinics, hyper-numerical-aperture lithography [4] [5] [6] [7] [8] where controlled polarization enhances the contrast and thus enabling smaller structures to be written on the wafer.
Inherited from standard interferometry [9] , the doublepass configuration detecting the phase shift between its two arms has been developed for sensing applications such as dilatometric measurement [10] and pH monitoring [11] . In polarimetry, a double-pass layout enables angle-resolved measurements, whereby the polarization response of a sample for a range of polar angles and all azimuths can be measured in a synchronous approach. This simplifies the measurement setup and saves time compared to, otherwise, an apparatus with a function of rotating a solid angle over a certain range. Since the light is transmitted through the sample being tested twice, each ray nominally picks up the same polarization properties in both the outgoing and return paths. Given the same apparatus errors outside the sample being tested, the double-pass configuration offers double the sensitivity of the polarization properties. The interferometric merit of the double-pass, on the other hand, is utilized in aligning the optical components in the angle-resolved polarimetry. While experimental demonstrations have validated the concept of double-pass polarimetry in angle-resolved polarization measurements [2] , repeatability analysis to tolerances of the double-pass polarimeter has not been * yuym3133@smee.com.cn; all authors contributed equally. studied systematically. The present work attempts to fill this gap by providing a detailed sensitivity analysis of the polarimeter repeatability.
An example of the operation of a double-pass polarimeter includes the calibration and measurement procedures. In the calibration apparatus as illustrated in Fig. 1a , a coherent laser illuminates a polarizer (P1) and a quarterwave plate (Q1) successively before being reflected by a non-polarizing beam splitter (BS). The coherent laser, the polarization components P1 and Q1, together with the reflective path of the BS, form the polarization state generator (PSG). The laser beam then passes through the calibration sample in the forward and reverse directions with the help of a mirror. The change of the polarization state of the beam caused by the calibration sample and the mirror is analyzed by the polarization state analyzer (PSA) and readout from the CCD. The PSA consists of the transmission path of the BS, the quarter-wave plate Q2 and the polarizer P2. The goal of the calibration setup is to characterize the polarization properties of the PSG and PSA accurately using calibration samples and the eigenvalue calibration method (ECM) [1] . The polarization properties of the calibration samples can be extracted using the same setup. In the measurement procedure the calibration samples and the mirror are subsequently replaced with an objective lens, the sample under test (SUT) and a hemispherical mirror as shown in Fig. 1b . The focus of the laser beam from the objective is aligned to coincide with the center of the curvature of the hemispherical mirror, to ensure that the beam is reflected back along the incoming optical path. The SUT is placed away from the focus for the laser beam to cover its pupil.
In this work, we break down the angle-resolved measurement of a SUT into 3 steps. In step 1, the trans- Step 2
Step 1
Step mittance amplitude for the two orthogonal polarization eigenstates and the retardance of the calibration samples are extracted from the calibration apparatus by comparing the intensities with and without the calibration samples.
Step 2 is an algorithmic procedure used to obtain the polarization properties of the PSG and PSA. This algorithm depends not only on the measured intensities with and without the calibration samples, but also on the polarization properties extracted in step 1. In step 3, the measurement setup employs the calibrated PSG and PSA to measure the polarization of the objective, SUT and hemispherical mirror together. Because the PSG and PSA are maintained unchanged during calibration and measurement, they cause no systematic change in the errors in measuring the SUT. Tolerance analysis of the components affecting the polarization from steps 1 to 3 results in the polarization measurement sensitivities.
We characterize the polarization of a non-depolarizing sample in terms of its diattenuation and retardance, which quantify the transmission amplitude difference between the two orthogonal brightest and darkest axes and the phase difference between the two orthogonal fastest and slowest axes, respectively. For non-depolarizing samples, the Jones matrix representation of the polarization is all that is required and is simpler than the Mueller matrix representation, in that the Jones matrix uses fewer parameters, only 4 complex elements compared to 16 real elements for the Mueller matrix. The diattenuation and retardance across the pupil can be expanded in terms of the orientation Zernike polynomials (OZP) based on the Jones matrix formalism [12] [13] [14] , and the RMS of the coefficients quantifies the diattenuation or retardance across the entire pupil by analogy with Zernike polynomials for wavefront expansion. By inputting tolerances of available commercial products into the numerical model, we predict a sensitivity of 0.5% RMS OZP (a unit stands for the RMS of the corresponding OZP coefficients) for a diattenuation pupil, equivalent to a pupil with a mean diattenuation of 1%. Likewise, the prediction of the sensitivity for a retardance pupil is 0.3 • RMS OZP corresponding to a pupil with a mean retardance of 0.6 • .
This numerical tool takes the systematic and stochastic errors of each component in the system for both the calibration and measurement as inputs, and derives the sensitivities of diattenuation and retardance to errors in the measured values in a bottom-up approach. Whereas double-pass polarimeters can find application for characterizing incident-angle dependent variable attenuators [15] , wide-view-angle polarizers and retarders [6, 8, 16] in lithographic equipment, this numerical tolerance analysis paves the way for predicting the sensitivity of the polarization properties for those optical components. Furthermore, this numerical tool can help to improve sensitivity via error-budgeting [17] . Depending on the relative contribution of each tolerance error, targeted hardware could be replaced to improve the sensitivity.
Step 1: Determining the properties of calibration samples
Classical calibration procedures usually rely on standard samples with well-known properties [18] or similar devices with higher accuracy. The former approach requires strict sample fabrication, while the later one limits the accuracy of the polarimeters to be calibrated to, roughly, that of the calibrating polarimeter. The ECM developed by Compain et al. [1] largely relaxes the requirement for special calibration samples, and is able to extract the polarization properties of the calibration sample from the polarimeter itself, hence nominally guaranteeing measurement accuracy. The ECM uses linear dichroic polarizers and retarders with retardation far from 180 • [1, 19] . These polarization elements need to be homogeneous [20] . That is their eigen polarization states of polarizing elements are orthogonal. Here we extend the ECM to double-pass polarimetry. Due to the flat mirror in the double-pass layout sketched in Fig. 1a , wave plates with retardance of 90 • are excluded from use as calibration samples. A dichroic polarizer and a 1/6-wave plate are selected as calibration samples in this work.
Intensities modulated by the PSG and PSA are recorded. The calibration sample is first retracted from the optical path in the setup in Fig. 1a , leaving only the mirror. This results in the intensity matrix i 0 
Here matrix a is the calculated PSA matrix from the intensity measurement. It is constructed from the 1 st to the u th configuration of the PSA, using the first row of the Mueller matrices of the PSA. The calculated PSG matrix w is formed by v different configurations of Stokes vectors. The middle term m mirror on the right hand side (RHS) of Eq. (1) is the measured Mueller matrix for the mirror. We then insert the dichroic polarizer and the 1/6-wave plate separately to obtain the intensity matrices
in which the subscript stands for the i th calibration sample. Matrices m f i and m b i can be further decomposed to
is the rotation matrix corresponding to azimuthal rotation angle θ of the calibration samples. The superscript f and b denote that the light passes through the calibration sample in a forward path and a backward path after reflection from the mirror, respectively. The measured Mueller matrix of the dichroic polarization elements m i , with zero azimuthal angle, can be expressed as [1] 
in which t X and t Y are the measured transmittance amplitudes of the sample along the two orthogonal directions, X and Y. We define the Z direction of the coordinate to be aligned with the ray propagation direction, the X direction to be pointing inside, and the Y direction to be pointing upwards at the start of the beam near the laser as demonstrated in Fig. 1b . The measured retardance difference between the X and Y directions is φ.
The quotient matrix c i is defined as the product of the inverse of the intensity matrix i 0 and the matrix i i , which gives
The Mueller matrix of the mirror m mirror is in the form of Eq. (3), where the non-identity of the reflectance is expressed by the transmittance amplitudes t X and t Y , and the retardance φ of the mirror in Eq. (3) is taken to be the sum of 180 • and noises. The last relation (≈) becomes an equality when no noise is present in the measured intensity matrices i 0 , i 1 or in the control of the azimuthal angle θ of the calibration samples. To ensure the uniqueness of the solutions, the full rank of the PSG's w matrix is required for the inversion in Eq. (4) and the PSA's matrix a has the same requirement. The true combinations of the azimuthal angles of the polarizing elements in the PSG are chosen to maximize the absolute value of the determinant of the true PSG matrix W in order to minimize the inversion error of W in calculating the Mueller matrix for the calibration sample. Here, we use the convention that matrices with uncapitalized and capitalized letters symbolize the measured (or calculated) values and the actual (or true) values, respectively. The true PSA matrix A is optimized in the same way. The coherent laser source beam in the PSG is modeled as a linearly polarized electrical field of E in = [1; 1]/ √ 2. The PSG uses 4 configurations in our simulations for convenience in performing the inversions, i.e., set v = 4 in Eq. (1). Each configuration is obtained by varying the azimuthal angles of the polarizer P1 and the quarter-wave plate Q1. The true values of the polarization properties of the PSG and PSA used in the simulation are summarized in Tab. I.
The maximum absolute value of the determinant of the PSG is optimized to |detW |=0.58 and that of the PSA is |detA|=0.06. The reflectance and transmittance amplitudes of the BS are idealized to be √ 0.5 in this modeling. Note that the Mueller matrix of the BS only affects the optimization of the azimuthal angle configurations for P1, Q1, Q2 and P2. It has no influence on 
the calibration error for the PSG ∆W = w − W or that for the PSA ∆A = a − A in step 2. In the experiments, the Mueller matrix for the transmission and reflection paths through the BS could be measured in advance using a single pass polarimeter in transmission [19, 21] and reflection [1] to ensure the calibration accuracy. The quotient matrix c i in Eq. (4) is similar, in the linear algebra sense, to the square of the measured Mueller matrix of the calibration sample [m i ] 2 given that matrix [R(θ)w] is invertible. Therefore, the quotient matrix c i and [m i ] 2 share the same eigenvalues. While the transmittance amplitudes can be calculated from the two real eigenvalues λ 1 and λ 2 , as t X = 4 √ λ 1 and t Y = 4 √ λ 2 , the retardance of the calibration sample φ is a function of the two complex eigenvalues λ 3 and λ 4 , as
Error sources depending on the measurement time scale are categorized into stochastic noise and systematic errors. Characteristic time scales are the total measurement time for intensities without the calibration samples i 0 , those with the calibration samples i i and the sampleswitch time in between. Stochastic noise with a time scale shorter than the total measurement time, comes from the laser source, the CCD, vibration of the rotatory positioners and the mechanical mounts of the optical elements. Each pixel of the CCD has a fluctuation of ±0.3% in the measured intensity which is modeled as statistically the same, and comes primarily from the repeatability of the laser source [22] and the random spatial non-uniformity in the CCD [23] . Both cross-talk between neighboring pixels and electrical shot noise contribute to the spatial non-uniformity. Cross-talk is simulated via the correlation length of these noise sources across the CCD. The correlation length is taken to be 1 pixel for simplicity, i.e. no cross-talk is assumed. For longer correlation lengths, filtering algorithms may be applied to reduce the noise influence. The impact of electrical shot noise on the signal to noise ratio decreases as the number of photons increases (assuming the photon-to-electron conversion rate of 1). By carefully selecting the measurement conditions so that the CCD is near saturation, electrical shot noise buried in the signal controlled by the power of the laser and the integration time of the CCD can have less than 1/10 of the influence on diattenuation and retardance caused by the quantization noise due to the analogueto-digital conversion (ADC) of the CCD. Electrical shot noise can therefore be safely neglected under the assumption of near CCD saturation, 10 14 photons per pixel in the model. The stochastic vibration of the rotatory positioners attached to polarization elements P1, Q1, Q2, and P2 in axial direction is taken to be 0.01 • . This follows from the Thorlabs' motorized rotator K10CR1 [24] specifications. Tilted variation of the PSG and PSA on the other hand is allocated to the polarization properties of the mirror and the BS in addition to the stochastic noise of the retardance and reflectance across the mirror. The pre-measurement of the mirror can be performed by a single pass polarimeter in reflection mode using analysis [1] similar to this step to obtain those stochastic noise. The difference lies in that for the double-pass layout the light probes a SUT in both forward and backward directions, while in the single pass polarimeter the light incidents on a SUT (mirror here) only once. To calibrate the mirror under normal incidence, an additional BS is required to deflect the beam from reflection in the single pass polarimeter and should be calibrated in advance. The tolerance types for stochastic noise and their values are summarized in Tab. II. Elements of the measured intensity matrix in Eq. (1) equal to the true values plus errors, i x,y 0 = I x,y 0 + ∆I x,y 0 (where x = 1, 2, ..., u; y = 1, 2, ..., v). The measured azimuthal angle θ = Θ + ∆Θ in Eq. (4) is the sum of true value Θ and precision ∆Θ of rotatory positioners. The reflectance error of the mirror ∆R X/Y = r X/Y − R X/Y , the retardance error ∆Φ = φ − Φ and the transmittance amplitude error ∆T X/Y = t X/Y − T X/Y all follow the same convention.
Although the PSG and PSA nominally have systematic errors, with the settings of the PSG and PSA being the same between the calibration (Fig. 1a) and measurement (Fig. 1b) setups there is no systematic change in the errors for the PSG and PSA matrices w and a. Therefore, the systematic error introduced in the modeling comes from the mirror in steps 1-2, the objective and the hemispherical mirror in step 3. The CCD is a common element in both the calibration and measurement layout, nevertheless, information loss in the process of ADC cannot be calibrated out. Hence the systematic error from the CCD must be included in all 3 steps. Systematic errors are listed in Tab. III.
With the stochastic and systematic errors of each component of the polarimeter listed above, we simulate both the stochastic noise and systematic errors in the properties of the calibration samples using a bottom-up approach. To reduce the rotational asymmetrical noise such as the tilt angle of the calibration samples, we rotate the calibration sample azimuthally and take the average over the -90 • to 90 • range. Figure 2a shows the calibrated transmittance amplitudes t X , t Y and retardance φ as a function of the azimuthal angle of the 1/6-wave plate sample. The average measurement value is displayed as a red line with the true values of the transmittance amplitudes being T X = 0.98, T Y = 0.97 and that of the true retardance being Φ = 60 • . The stochastic noise is defined as the difference between the average measurement over all azimuths and the true value. A 1000-trial simulation in Fig. 2b indicates that the calibration stochastic noise for transmittance is ∆T X < ±0.0002, ∆T Y < ±0.0002 and that for retardance is ∆Φ < ±0.012 • as shown in Tab. II. For the calibration sample polarizer, this step is sufficient to determine the transmittance amplitude of the bright transmission axis, but not the dark axis. Using two 40 dB polarizers in series could ensure a stochastic error ∆T Y < ±0.0002. The calibrated polarization properties of the 1/6-wave plate and polarizer, together with their stochastic noise determines the calibration accuracy for the PSG w and PSA a matrices in step 2.
In this subsection, we calculate the PSG w and PSA a matrices as well as expand the working range of the azimuthal angles of the calibration samples from those used in the past [25] . To compute the PSG matrix w, w in Eq. (4) is first replaced with an unknown matrix x. This results in m DP i x − xc i = 0, where the sample matrix for the double-pass polarimeter is defined as
The linear operator h i is a 16×16 matrix, the elements of which are detailed in Eqs. (A3) and (A4) in Appendix A. The calibrated PSG matrix w is then the non-zero solution to Eq. 
In this way, all 16 calculated eigenvalues λ(1) < λ(2) < . . . < λ(16) of k must be positive and real. The eigenvector with eigenvalue closest to 0 is the calculated PSG w matrix, after the 16×1 eigenvector being reshaped into [19] , and the first sample is a polarizer with azimuthal angle Θ = 0 • . To find suitable combinations of azimuthal angles that guarantee calibration accuracy, we plot the error estimator log[λ(2)/λ(1)] as a function of the azimuthal angles for calibration sample 2 (a polarizer with different azimuthal angle to sample 1) and sample 3 (a 1/6-wave plate) in Fig. 3a . Both of their azimuthal angles are varied from -90 • to 90 • .
Stochastic noise contributions to the error estimator include the rotational repeatability of the calibration samples as limited by mechanical positioners, intensity fluctuations, stochastic polarization noise of the calibration samples in step 1 and that of the mirror. Systematic errors come from the quantization error of the ADC, and the polarization properties of the mirror. Values of these errors are given in Tabs. II and III. Each combination of azimuthal angles in Fig. 3a is averaged over 100 trails to reduce the influence from stochastic noise. The larger the value of the error estimator, the closer the smallest eigenvalue of k in Eq. (6) is to 0, and consequently the more accurate the calculated PSG w matrix will be. We observe that the error estimator is relatively small, log[λ(2)/λ(1)] < 9, when the 1/6-wave plate has the same azimuthal angle (Θ ≈ 0 • , the middle horizontal reddish line in Fig. 3a) as that of the first polarizer. It is likely that the lack of calibration accuracy is due to the azimuthal angle overlap of the two orthogonal eigenstates of the first polarizer and the 1/6-wave plate, blurring the precision of the eigenvalue-based calibration method.
We further calculated the calibration error between the calibrated and the true PSG matrices ∆W and that of the PSA error matrix ∆A. The calculated PSG w matrix is normalized by its transmission before the comparison with the true PSG matrix W , because the eigenvector of Eq. (6) can be scaled with any real number. As the calculated PSA matrix a is derived from the measured intensity using Eq. (1), it will give an inverse scaling factor to the calculated PSG matrix w if the normalization is not done. Consequently, normalization of transmission only serves for obtaining the error for the PSG ∆W and the PSA ∆A. The PSG w and PSA a matrices without normalization will not affect the measurement accuracy of a SUT in step 3. The logarithm of the error of the PSG ∆W and the PSA ∆A as a function of the azimuthal angles of the 1/6-wave plate and the polarizer are plotted in Fig. 3b and Fig. 3c , respectively. The first element (1,1) of the 4×4 error matrices can be chosen without loss of generality. The other 15 elements of the error matrices ∆W and ∆A share roughly the same calibration error. The cross areas in the middle of the error matrices for the PSG and PSA display a relatively worse accuracy, and are aligned with the error estimator map, log[λ(2)/λ(1)] in Fig. 3a . As a result, the requirement for alignment of the calibration samples can be relaxed to all the yellowish areas in Fig. 3a , corresponding to the error estimator log[λ(2)/λ(1)] > 10. Former experimental observations reveal the calibration accuracy of the PSG and PSA matrices, where an average of a standard deviation over all 16 Mueller matrix elements is employed for quantification [19] . In those experiments, the averaged standard deviation is 6.7×10 −4 for the PSG matrix and 6.0×10 −4 for the PSA matrix over 38 calibrations. We simulate the pixel-based PSG and PSA matrices for 10000 trails, and obtain the averaged standard deviation of 5.9×10 −4 for the calibrated PSG and that of 3.6×10 −4 for the calibrated PSA, which is in line with the experiments, verifying our tolerance analysis for the calibration.
The alignment of the objective to the center of the hemispherical mirror can be monitored by adding an interferometer arm to form an interference pattern on the CCD. This added arm would extend horizontally from the laser and the PSG, and have a mirror at the end.
The simulation flow leading to the prediction of the sensitivities for the angle-resolved measurements is depicted in Fig. 4a . The simulation uses a generated Jones matrix covering the whole pupil (in short Jones pupil matrix) as the true Jones pupil matrix of a SUT J true . It is synthesized by the RMS of the coefficients of up to order 72 in an expansion using the OZP for diattenuation and retardance [12] [13] [14] .
The Jones pupil matrix is converted to a Mueller pupil matrix to be compatible with the Mueller matrix description M true of the PSG and PSA in step 1-2. The true PSG matrix W , PSA matrix A, stochastic noise from the hemispherical mirror (whose values are listed in Tab. II), combined with the systematic errors from the objective and hemispherical mirror (whose values are listed in Tab. III), result in the true intensity I true . Objectives usually contain multiple lenses to ensure a specific image quality over the field of view. As polarization relies on the order of the components the light passes through, the polarization of the light traveling through the objective in the forward direction, from collimated space to the focus as sketched in Fig. 4b tom, differs from the light transmitted by the objective in the backward or return direction. As an example, we choose a Japanese patent 61 2925 860129 in the CODE V database [26] with a half incident angle of 25.4 • to investigate its polarization properties. Without applying optical coatings to the objective, we trace the polarization of the objective in both the forward and backward directions. The backward beams exiting the object have maximum deviation angles of 0.29 • along the periphery due to the imperfect wavefront of the objective. Systematic errors of diattenuation and retardance of the objective considering the retrace error in the backward direction are listed in Tab. III. Random intensity noise at each pixel and information loss from the ADC are added to the intensity as error sources to form the measurement intensity I measure . The calibrated PSG matrix w, PSA matrix a, and idealized Mueller matrix [1 0 0 0; 0 1 0 0; 0 0 -1 0;0 0 0 -1] are employed to calculate the Mueller matrix of the SUT in the forward path M measure . The Mueller pupil matrix is converted to the Jones pupil matrix afterwards. This procedure removes the information about depolarization contained in the Mueller matrix to obtain the Jones matrix. Depolarization in the measured Mueller matrix M measure comes from overlap of incoherent electromagnetic fields [3] . To convert the Mueller matrix with limited depolarization to the Jones matrix J measure , the non-depolarization condition for the conversion trace(M T M ) = 4m 2 11 [27] is approximated as |trace(M T M ) − 4m 2 11 | < 0.01. For measurements that meet this condition, the Jones matrix can be derived from the Mueller matrix via expressions given in Ref. [28] .
The RMS of the OZP coefficients for either diattenuation or retardance is a single number used to quantify the goodness of a Jones pupil via the relative transmittance amplitude difference between the brightest and darkest axes or retardance delay between the fastest and slowest axes across the pupil of a SUT, respectively. Mathematical details of the OZP can be found in Appendix B. For 72 terms the highest power in the radial direction of the OZP is 10, corresponding to the highest radial power of the 36 th term of the fringe Zernike polynomials [26] .
Though the true reflectance amplitude of the hemispherical mirror is not unity, only the difference between the reflectance in the X and Y directions affects the diattenuation and retardance of the pupil. This is because the measured Jones pupil matrix J measure is further decomposed into a product of apodization, a partial polarizer, a retarder and two other physically meaningful matrices [5] , and only the diattenuation pupil in the partial polarizer and the retardance pupil in the retarder will be further expanded by the OZP. Writing the reflectance in the X and Y directions of the hemispherical mirror as r X = r Y + ∆R XY , the average of the reflectance in the X and Y direction contributes only to the apodization of the SUT. The difference of the reflectance amplitudes ∆R XY will be counted in the first term of the OZP expansion (see Eqs. (B6) and (B7) for the mathematics). Since rotating the hemispherical mirror azimuthally for 90 • swaps the reflectance values r X and r Y , taking the average of the fitting coefficients to the OZP expansion, measured with 0 • and 90 • hemispherical mirror rotation, improves the accuracy of the OZP coefficients for the diattenuation and retardance pupils.
We decompose both J true and J measure into an OZP description of retardance and diattenuation, using the first 72 terms. The RMS of the coefficients are calculated as RMS = 72 j=1 coe j /(j + 1), with coe j denoting the j th OZP coefficient. Comparison of the true Jones pupil matrix for the SUT and the measured value is made by running the simulation through the flow in Fig. 4a for 100 trials. The repeatability in terms of the RMS of the OZP coefficients replaces the mean value in the standard variance [29] with the true value, defined as
Before predicting the sensitivity presented by RMS OZP in the Jones matrix description, we apply our tolerance analysis to the SUT in terms of the Mueller matrix in the measurement procedure similar to that reported experimental observations in Ref. [1] . We compare the true Mueller matrix M true with the measured Mueller matrix M measure in the simulation flow as sketched in Fig. 4a . Both of the two matrices are normalized to their (1, 1) elements, so that the relative error of the (2, 2), (3, 3) and (4, 4) elements of the matrices can be calculated under the condition of a non-identity Mueller matrix of the mirror. Off-diagonal elements of the Mueller matrices M true and M measure are small due to the generated weak polarization properties of the SUT, leading to unphysically large relative errors, and thus they are safely disregarded in the comparison. We obtain a maximum 0.4% over all three Mueller matrix pupils, in good agreement with the 0.5% in the reported experiment.
Sensitivity is defined in terms of a boundary. In Fig. 5 , the boundary where repeatability equals the true value is the line with a slope of 1 through the origin (0,0). Away from the gray shadow areas, the repeatability (i.e. the measurement uncertainties) are smaller than the true values. The sensitivity of the diattenuation pupil depends on the corresponding retardance. Larger retardance leads to better sensitivity for diattenuation in general. The same phenomenon applies to the sensitivity of retardance as well. It is likely that the measurement is more sensitive when the SUT exhibits strong polarization properties, and the retardance and diattenuation are not decoupled in calculating the repeatability of either of them. To reduce the sensitivity from a set of values to a single value, we quantify the sensitivity of diattenuation with an additional requirement: the corresponding retardance of the pupil should be of the same order of magnitude as the diattenuation. This results a 0.5% RMS OZP sensitivity for diattenuation. With the same requirement, the predicted sensitivity for retardance is 0.3 • RMS OZP.
Visualization of the pupils for the true diattenuation and retardance when their sensitivities are reached (labeled with black arrows in Fig. 5) , i.e. 0.5% RMS OZP for diattenuation and 0.3 • RMS OZP for retardance, is shown in Fig. 6a and Fig. 6d , respectively. The mean of the measured pupils of diattenuation and retardance comes from the measured Jones pupil J measure , and a decomposition of the measured Jones pupil in terms of diattenuation and retardance thereafter. Reconstruction of the pupils for diattenuation and retardance is based on the first 72 terms of the OZP expansion, where each pixel of the pupils for diattenuation and retardance is averaged over 100 trials. The sensitivity of diattenuation shows an average of 1% over all pixels of a pupil displayed in Fig. 6b with repeatability around 1/3 of that displayed in Fig. 6c . For the sensitivity of retardance, the average Fig. 4a . Black arrows point to the pupils that meet our definition of sensitivity, and polarization properties of these pupils are visualized in Fig. 6. of the pupil is 0.6 • as shown in Fig. 6d with a repeatability around 1/3 of that as well, as shown in Fig. 6f . Directional lines on the diattenuation pupils denote azimuthal angles for the partial polarizer, while they denote those for the retarder on the pupils of retardance. The azimuthal angle pupil reconstructed from the OZP coefficients may have a 90 • shift, due to the limitation of the inverse trigonometric functions described in Eqs. (B14) and (B15) in Appendix B. Horizontal lines represent an azimuthal angle of 0 • , while vertical lines represent 90 • . White lines with directions other than vertical or horizontal represent error, the larger the error of the direction, the farther away the direction of the white line is from either vertical or horizontal. When treating the 90 • shift to be error-free, the repeatability of the azimuthal angle is 5 • for diattenuation and 3 • for retardance averaged across the pupil, lower than 1/10 of the mean values of that for diattenuation and retardance. Concentric circles with different radii on the pupil correspond to different incident angles of the laser beam away from the focal plane of the objective. In Fig. 6 , from inside to outside the concentric circles correspond to angles at the objective of 5 • , 15 • and 25 • . Pupils of diattenuation and retardance provide a visualization of the azimuthal anisotropy and polarization response of a refractive sample under non-normal incidence.
In conclusion, we have performed a detailed tolerance analysis of the calibration and measurement procedures for a double-pass polarimeter, and have predicted the sensitivity of the polarimeter to systematic errors and stochastic noise. The eigenvalue calibration method ECM [1] is used in the polarimeter calibration, resulting in the Mueller and Stokes description of the PSG and PSA characteristic matrices. The Mueller pupil matrix of an arbitrary non-depolarizing SUT is predicted before it is converted to a Jones pupil matrix. Our tolerance model for the calibration of the PSG and PSA, as well as the measurement of the Mueller matrix pupil are consistent with previous experimental observations [1, 19] . Thanks to the Jones pupil decomposition and the OZP expansions of diattenuation and retardance, the whole pupil of the SUT can be described by two values, diattenuation and retardance in terms of the RMS of the OZP coefficients. The sensitivity prediction for diattenuation is 0.5% and that for retardance is 0.3 • . The double-pass polarimeter offers a platform to measure angle-resolved SUTs, revealing the azimuthal inhomogeneity of retardance and diattenuation. The ECM, tolerance analysis and the subsequent conversion of the measured Mueller pupil matrix of the SUT to a Jones pupil matrix in terms of the OZP expansions to predict sensitivities and visualize retardance and diattenuation pupils can also be applied to a single pass polarimeter. Though the incident angle would not be resolved in the single pass polarimeter, without a BS and a mirror fewer noise sources are included. The singe pass polarimeter can achieve better sensitivity of diattenuation and retardance as well as resolve the small inhomogeneity of the pupil under normal incidence.
The authors are thankful to Zejiang Meng for introducing the ECM algorithm, Vladimir Nikishkin for coding assistance, and Wei Wang for the initial contact with vendors for the specifications of the components used in the modeling as well as his exuberant personality. Last but not the least, we acknowledge an extended vacation due to COVID-19 outbreak.
We have derived a simplified theory of error propagation for double-pass polarimetry to cross-check our numerical simulations with tolerances. By employing perturbation theory to the first order, we theoretically calculate the error of the PSG matrix ∆W given the measurables without any calibration samples for i 0 and with the calibration sample for i i .
The noise propagation of the quotient matrix c i combining Eqs. (A2) and (4) to the first order results in the expression
Letting the unknown x = W + ∆W , the linear operator h i (x) in Eq. (5) is expanded to the first order with the quotient matrix from Eq. (A2), as
Factoring the first term on the RHS in Eq. (A3), we operate on the elements of the matrices. It follows that by applying the relation for the least square fit ∆W p,q = δ p,F δ q,G ∆W F,G , where δ is the Kronecker delta, p, q, F and G are summed from 1 to 4, we have
The single indices µ and ν label all possible combinations of F, G and p, q. The last two terms on the RHS in Eq. (A3) are influenced by the intensity with the calibration sample ∆(AM b i M mirror M f i W ) and without it ∆(AM mirror W ). Hence, the intensity error is defined
. Assuming G i,µ,ν is invertible, Eq. (A3) can be simplified to
This expresses the linear relationship between one element of the PSG error matrix ∆W ν and the sum of the stochastic noise of the calculated calibration sample ∆(M DP i ) times the true PSG matrix W and the noise of the measured intensity ∆I i,µ .
To verify the validity of our numerical tool for the tolerance analysis, we simulate the PSG error ∆W matrix as a variation of the intensity error. The stochastic noise of the calibration sample ∆(M DP i ) in Eq. (A6) is idealized to be 0. Modeling results show that the error across the pupil of one element of the 4×4 PSG matrix ∆W 3,1(β=9) increases linearly with the intensity noise as expected as shown in Fig. 7 . The intensity noise normalized by the intensity of the PSA and PSG varies from ≈ I(r, ω) + j coe j OZ j (r, ω) (B6) J ret (φ, β, r, ω) = cos φ(r, w) 2 I(r, ω)
cos 2β(r, ω) sin 2β(r, ω) sin 2β(r, ω) − cos 2β(r, ω)
≈ cos φ(r, ω) 2 I(r, ω) + i j coe j OZ j (r, ω),
where the approximation sin φ(r,ω) 2 ≈ φ(r,w) 2 is used. The term OZ j (r, ω) is further decoupled into a position (r) dependent term and an orientor matrix depending on the azimuths ω, as 
where n indexes the highest power in radial direction and =0, 1 for the 2 orientor matrix in Eq. (B5). The order label j represents combinations of the OZP indices m, n, with the relation n−m = 2l, l = 0, 1, ...n, n ∈ Z + . Corresponding relation between j and m, n, up to the first 16 terms of the OZP is displayed in Tab. IV. The OZP expansion of diattenuation and retardance are approximations. To test the accuracy of these approximations, we use the first 72 orders of the OZP. A to-be OZP expanded and reconstructed diattenuation pupil consists of 4 pupils of elements, among which two pupils are independent. We label the upper-left element in the matrix as Jinput dia xx = d 2 cos 2γ and that in the upper-right as Jinput dia xy = d 2 sin 2γ. Similarly, two independent matrix elements for retardance are Jinput ret xx = sin φ 2 cos 2β and Jinput ret xy = sin φ 2 sin 2β. As shown in Fig. 8 , the pupils of the two independent elements in diattenuation or retardance matrices are compared between the reconstruction from the coefficients to the OZP and the inputs. Two independent elements of the matrix j coe j OZ j (r, ω) are reconstructed from the OZP coefficients, as
Jreconst xx (r, ω) = Fig. 6a . b The input pupil for retardance is from Fig. 6d . The difference between the reconstructed and the input original pupil based on the first 72 orders of the OZP is around 1/10 of either of them.
Differences from the input Jinput and reconstructed Jreconst Jones matrix pupils are an order of magnitude less than either the input Jinput or the reconstructed Jreconst. Diattenuation and retardance pupils used for comparison comes from those in Fig. 6a and Fig. 6d , respectively. This difference is around 1/3 of that between the mean of the measured pupil and the repeatability pupil for either diattenuation (in Figs. 6b-6c) or retardance (in Figs. 6e-6f) , demonstrating that the OZP expansion well represents diattenuation and retardance pupils. Therefore, errors contributed from the reconstruction with the first 72 orders of the OZP is negligible in calculating the sensitivities for diattenuation and retardance. The reconstruction of the diattenuation and retardance pupils as well as the direction of the partial polarizer and the retarder from the two independent elements Jreconst xx and Jreconst xy is given by 