## Polarimetric data reduction: a Bayesian approach

Optics Express, Vol. 15, Issue 1, pp. 83-96 (2007)

http://dx.doi.org/10.1364/OE.15.000083

Acrobat PDF (664 KB)

### Abstract

In this paper, we introduce a general Bayesian approach to estimate polarization parameters in the Stokes imaging framework. We demonstrate that this new approach yields a neat solution to the polarimetric data reduction problem that preserves the physical admissibility constraints and provides a robust clustering of Stokes images in regard to image noises. The proposed approach is extensively evaluated by using synthetic simulated data and applied to cluster and retrieves the Stokes image issuing from a set of real measurements.

© 2007 Optical Society of America

## 1. Introduction

1. D. Miyazaki, M. Saito, Y. Sato, and K. Ikeuchi, “Determining surface orientations of transparent objects based on polarization degrees in visible and infrared wavelengths,” J. Opt. Soc. Am. A **19**,687–694 (2002). [CrossRef]

2. D. Miyazaki, M. Kagesawa, and K. Ikeuchi, “Transparent surface modeling from a pair of polarization images,” IEEE Trans. PAMI **26**,920–932 (2004). [CrossRef]

3. J. M. Bueno and P. Artal, “Double-pass imaging polarimetry in the human eye,” Opt. Letters. **24**64–66 (1999) [CrossRef]

4. S. D. Giattina, et al., “Assessment of coronary plaque collagen with polarization sensitive optical coherence tomography (PS-OCT),” Int. J. Cardiol. **107**,400–409 (2006). [CrossRef] [PubMed]

*N*(

*N*≥ 4 ) independent basis states. The complete set of the

*N*measurements yields a matrix equation, which relates the out-coming Stokes vector

**S**from the sample to the measured raw intensity data vector

**I**for each pixel. This matrix relation has to be inverted properly. This implies that an efficient calibration procedure of polarimetric imaging systems must be employed in order to extract the desired polarization-encoded images effectively. Many interesting studies in regard to this problem have been published in the recent literature, see for example [6–9

9. S. N. Savenkov, “Optimization and structuring of the instrument matrix for polarimetric measurements,” Opt. Eng. **41**,965–972 (2002). [CrossRef]

## 2. Problem statement

### 2.1 Observation model in the Stokes imaging polarimetry framework

*j*,

_{y}*j*). This can be summarized by the following equation:

_{x}**I**

*is the acquired image,*

_{m}**P**is the “polarization measurement matrix” (PMM) that depends on a parameter vector (η), while Sin is the Stokes vector to be estimated.

**I**

*(*

_{m}*j*,

_{y}*j*) will denote the intensity vector measured at the (

_{x}*j*,

_{y}*j*) pixel location while

_{x}**S**

*(*

_{in}*j*,

_{y}*j*) represents the corresponding Stokes vector.

_{x}**I**

*, and the Stokes image will be arranged in data cube structures such that:*

_{m}**P**

^{#}of the PMM to obtain a minimum-norm, least-square estimate

**S**̂

*for the Stokes vector at each pixel location:*

_{in}**P**

^{#}is identical to the inverse of

**P**provided that η is chosen such that

**P**is invertible. The reason for using

*N*intensity measurements where

*N*> 4 that yield an over-determined system for Eq. (1) is to reduce the system sensitivity to systematic errors as well as image noises. The drawback of such a strategy is an increase in acquisition, and data processing and handling time.

**I**incorporates all sources of possible noises.

**I**as the sum of two components, i.e. δ

**I**=δ

**I**

_{p}+ ε. The first noise term is proportional to the measured intensities while the second one includes other independent additive noises. ε is assumed to follow a Gaussian distribution and δ

**I**

*is attributed to photon noise that can be estimated by observing that, for each intensity measurement*

_{p}*I*, δ

_{i}**I**

*/*

_{i}*I*is proportional to the inverse of the square root of the number of detected photons (

_{i}*n*):

^{ph}_{i}*Q*) of the used photosensitive device (e.g. CCD camera) which in turn depends on the wavelength. This dependency can be expressed formally as:

_{E}*N*is the number of photons impinging onto the detector and λ is the considered wavelength.

^{ph}_{i}*N*intensity measurements where each one contains nearly the same amount of statistical noise, we can assume that the same number of photons

*n*(λ) contributes to each one. Finally, δ

_{ph}**I**

*can be written as:*

_{p}**I**

*= (1+ρ)*

_{d}**PS**

*.*

_{in}### 2.2 Underlying segmentation model

*S*(

_{in}*j*,

_{y}*j*,

_{x}*j*)=ϕ(

_{e}*j*,

_{e}*I*

_{0}(

*j*,

_{y}*j*)) and consequently, the k

_{x}^{th}column of the matrix

**Φ**represents the Stokes vector that corresponds to the k

^{th}class in the image. This hypothesis amounts to considering each Stokes parameter as being constant inside each class.

13. A. Gray, J. Kay, and D. Titterington, “An empirical study of the simulation of various models used for images,” IEEE Trans. PAMI ,**16**,507–513 (1994). [CrossRef]

14. A. Dunmur and D. Titterington, “Computational Bayesian analysis of hidden Markov mesh models,” IEEE PAMI. **19**,1296–1300 (1997). [CrossRef]

*β*

_{0}, which controls the size of the clusters: values of

*β*

_{0}close to 1 yield monochrome images, whereas

*β*

_{0}= 1/

*K*corresponds to an independent pixel wise prior. The proposed probabilistic model reads:

*k*∈[1,

_{i}*K*] are the four neighbors labels values. This third-order Markov mesh model, where a given pixel admits three causal neighbors, is equivalent to a second order MRF where a pixel admits eight neighbors.

*n*depend only on

_{i}**I**

_{0}.

### 2.3 Noise model

^{(jn)}, σ

^{(jn)}) are the noise parameters (mean and standard deviation) that affect the k

^{th}class in the jn channel image.

**Θ**as the set

**Θ**={

**Θ**

^{(jn)};

*j*∈ [0,

_{n}*N*−1]}, where

### 2.4 General probabilistic model

**I**

*represent fixed hyperparameters. Hyperparameters values are set so that the corresponding priors are weakly informative in the region of interest (see for example p. 735 of Ref. [15*

_{m}15. S. Richardson and P. Green, “On Bayesian analysis of mixtures with an unknown number of components (with discussion),” J. R. Stat. Soc. Ser. B. **59**,731–792 (1997). [CrossRef]

*p*(θ) is conjugate to the likelihood

*p*(

*y*∣θ) if the posterior distribution

*p*(

*y*∣θ) is in the same class of distributions as the prior

*p*(θ). Conjugate priors are often good approximations and they simplify computations. In the 10, 11] for extensive details on this topic. Let us emphasize that all probability laws derive directly from the Gaussian model used for the noise and from the MRF model used for

**I**

_{0}. Considering these assumptions and the conjugacy property, the other laws are fixed. This explains for example why the μ

^{(jn)}

*’s have a Gaussian distribution and why the σ*

_{k}^{−2(jn)}

*’s follow a Gamma distribution.*

_{k}15. S. Richardson and P. Green, “On Bayesian analysis of mixtures with an unknown number of components (with discussion),” J. R. Stat. Soc. Ser. B. **59**,731–792 (1997). [CrossRef]

*C*is the conjunction of conditions ϕ(0,

_{k}*k*) ≥ 0 and

*U*,

*Ga*, and

*Be*define respectively the uniform, Gamma, and Beta laws, as:

## 3. Estimation procedure

*p*(

**Θ**,

**Φ**,

**I**

_{0})

*ρ*

*α*

_{4}

*K*β

_{0}∣

**I**

*which verifies:*

_{m}**I**

_{0}is updated using a single raster scan (Iterated Conditional Modes, see Ref. [16]) and that a local optimization method is used to update ϕ (:,

*k*) if the solution to the unconstrained problem, e.g., not accounting for

*C*, violates the conditions

_{k}*C*. The number

_{k}*K*of classes is not updated during the procedure but optimizations considering all values of

*K*∈[

*K*,

_{min}*K*are achieved.

_{max}**S**̂

*obtained by Eq. (3). The initialization of*

_{in}**I**

_{0}relies on the

**S**̂

*clustering by a variant of the C-means algorithms family where the used distances were redefined in relation with the Stokes images specificities [8*

_{in}8. S. Ainouz, J. Zallat, A. de Martino, and C. Collet., “Physical interpretation of polarization-encoded images by color preview,” Opt. Express **14**,5916–5927 (2006). [CrossRef] [PubMed]

**S**̂

*corresponding to each class. Initializing the other variables follows straightforwardly as given in the Appendix.*

_{in}## 4. Discussion and analysis

### 4.1 Simulation results

_{1}=(1.0,0.5,0,0.866)

*was assigned to the black pixels, s*

^{t}_{2}=(0.8,0,0,0.8)

*to the heavy gray pixels, s*

^{t}_{3}=(0.9,0.39, −0.675,0.45)

*to the light gray pixels, and s*

^{t}_{4}=(1.2,0.85,0.85,0) to the white pixels in the label map. Fig. 2(b) shows these polarization state locations on the Poincaré sphere.

7. J. S. Tyo, “Design of optimal polarimeters: maximization of signal-to-noise ratio and minimization of systematic error,” Appl. Opt. **41**619–630 (2002). [CrossRef] [PubMed]

*θ*= ±15.12°,±51.7°}) where δ stands for the retardance of the birifringent wave-plate while

_{i}*θ*represent its angular positions.

_{i}*σ*

^{2}

*was added to each intensity image for this simulation providing the necessary inputs to our model, e.g. intensity images as well as the PMM P. The aim of this section is to compare the estimated Stokes image*

_{n}**S**̂

*resulting from Eq. (3) with the one resulting from our Bayesian model (*

_{LS}**S**̂

*).*

_{B}**S**̂

*for a noise variance*

_{LS}*σ*

^{2}

*= 0.01. Figures 4(a) and 4(b) show the polarization states locations corresponding to*

_{n}**S**̂

*over the Poincaré sphere for*

_{LS}*σ*

^{2}

*= 0.01(SNR ∼ 20dB) and*

_{n}*σ*

^{2}

*= 0.1(SNR ∼ 10dB). Points that lie outside the Poincaré sphere correspond to unphysical estimated states. The percentage of unphysical estimated pixel by the min-norm method are respectively 54 % for*

_{n}*σ*

^{2}

*= 0.01 and 57 % for*

_{n}*σ*

^{2}

*= 0.1.*

_{n}*σ*

^{2}

*= 0.1. Table 1 summarizes the estimated Stokes vectors for the different classes. The last two rows lists the RMSE defined as 100×∥s*

_{n}*−ŝ*

_{i}*/s*

_{i}*∥ as well as the ratio of misclassified pixels of each class for the two considered variances. As expected the estimated ρ factor is practically zero. Figure 5 shows the true polarization locations over the Poincaré sphere as well as the estimated ones for the two considered noise variances.*

_{i}*σ*

^{2}

*)*

_{n}_{i=1,4}= 0.1,0.01,0.05, and 0.2 that attain the four classes defined in Fig. 2(a).

*σ*̂

^{2}

*)= 0.1,0.01,0.049, and 0.203. Figure 7 shows the location of the estimated polarization states and the true ones over the Poincaré sphere. Again, we observe the excellent behaviour of our model for the case of different noises that reach each class.*

_{i=1,4}### 4.2 True measurement case

_{0}element and are given by:

**s**̂

_{1}corresponds to the smallest patch,

**s**̂

_{2}to the trapezoidal shape, and

**s**̂

_{3}to the background.

- ● The polarization state corresponding to
**s**̂_{1}is located near the center (Fig. 10) indicates that the transmission axis of this patch is oriented horizontally which extinguish a large amount of the vertically polarized incident beam. - ● The polarization state corresponding to
**s**̂_{2}lies in the (S_{1}-S_{2}) plan indicating a linear state as expected. - ● The incident beam undergoes isotropic depolarization induced by the diffusing glass slide. This is confirmed by observing that the degree of polarization (DOP) of
**s**̂_{2}and**s**̂_{3}are nearly equal*DOP*(**s**̂_{1}) ≈*DOP*(**s**̂_{2}).

## 5. Conclusion and future extensions

## Appendix: hyperparameters choice

15. S. Richardson and P. Green, “On Bayesian analysis of mixtures with an unknown number of components (with discussion),” J. R. Stat. Soc. Ser. B. **59**,731–792 (1997). [CrossRef]

**S**̂

*is the min-norm estimate of the Stokes image.*

_{LS}## References and Links

1. | D. Miyazaki, M. Saito, Y. Sato, and K. Ikeuchi, “Determining surface orientations of transparent objects based on polarization degrees in visible and infrared wavelengths,” J. Opt. Soc. Am. A |

2. | D. Miyazaki, M. Kagesawa, and K. Ikeuchi, “Transparent surface modeling from a pair of polarization images,” IEEE Trans. PAMI |

3. | J. M. Bueno and P. Artal, “Double-pass imaging polarimetry in the human eye,” Opt. Letters. |

4. | S. D. Giattina, et al., “Assessment of coronary plaque collagen with polarization sensitive optical coherence tomography (PS-OCT),” Int. J. Cardiol. |

5. | D. H. Goldstein and D. B. Chenault, and Society of Photo-optical Instrumentation Engineers, Polarization: measurement, analysis, and remote sensing II, 19–21 July, 1999, Denver, Colorado. 1999, Bellingham, Washington: SPIE. ix, 426 p. |

6. | M. H. Smith, “Optimizing a dual-rotating-retarder Mueller matrix polarimeter,” in |

7. | J. S. Tyo, “Design of optimal polarimeters: maximization of signal-to-noise ratio and minimization of systematic error,” Appl. Opt. |

8. | S. Ainouz, J. Zallat, A. de Martino, and C. Collet., “Physical interpretation of polarization-encoded images by color preview,” Opt. Express |

9. | S. N. Savenkov, “Optimization and structuring of the instrument matrix for polarimetric measurements,” Opt. Eng. |

10. | J. Bernardo and A. Smith, |

11. | A. Gelman, J. Carlin, H. Stern, and D. Rubin, |

12. | S. Z. Li, |

13. | A. Gray, J. Kay, and D. Titterington, “An empirical study of the simulation of various models used for images,” IEEE Trans. PAMI , |

14. | A. Dunmur and D. Titterington, “Computational Bayesian analysis of hidden Markov mesh models,” IEEE PAMI. |

15. | S. Richardson and P. Green, “On Bayesian analysis of mixtures with an unknown number of components (with discussion),” J. R. Stat. Soc. Ser. B. |

16. | J. Besag, “On the statistical analysis of dirty pictures (with discussion),” J. R. Stat. Soc. B |

**OCIS Codes**

(000.3860) General : Mathematical methods in physics

(100.3190) Image processing : Inverse problems

(110.2960) Imaging systems : Image analysis

(120.5410) Instrumentation, measurement, and metrology : Polarimetry

**ToC Category:**

Imaging Systems

**History**

Original Manuscript: September 7, 2006

Manuscript Accepted: October 16, 2006

Published: January 8, 2007

**Citation**

Jihad Zallat and Christian Heinrich, "Polarimetric data reduction: a Bayesian approach," Opt. Express **15**, 83-96 (2007)

http://www.opticsinfobase.org/oe/abstract.cfm?URI=oe-15-1-83

Sort: Year | Journal | Reset

### References

- D. Miyazaki, M. Saito, Y. Sato, and K. Ikeuchi, "Determining surface orientations of transparent objects based on polarization degrees in visible and infrared wavelengths," J. Opt. Soc. Am. A 19, 687-694 (2002). [CrossRef]
- D. Miyazaki, M. Kagesawa, and K. Ikeuchi, "Transparent surface modeling from a pair of polarization images," IEEE Trans. PAMI 26,920-932 (2004). [CrossRef]
- J. M. Bueno and P. Artal, "Double-pass imaging polarimetry in the human eye," Opt. Letters. 2464-66 (1999). [CrossRef]
- S. D. Giattina, et al., "Assessment of coronary plaque collagen with polarization sensitive optical coherence tomography (PS-OCT)," Int. J. Cardiol. 107, 400-409 (2006). [CrossRef] [PubMed]
- D. H. Goldstein, D. B. Chenault, and Society of Photo-optical Instrumentation Engineers, Polarization: measurement, analysis, and remote sensing II, 19-21 July, 1999, Denver, Colorado. 1999, Bellingham, Washington: SPIE. ix, 426 p.
- M. H. Smith, "Optimizing a dual-rotating-retarder Mueller matrix polarimeter," in Polarization Analysis and Measurements IV, SPIE (2001).
- J. S. Tyo, "Design of optimal polarimeters: maximization of signal-to-noise ratio and minimization of systematic error," Appl. Opt. 41619-630 (2002). [CrossRef] [PubMed]
- S. Ainouz, J. Zallat, A. de Martino, and C. Collet., "Physical interpretation of polarization-encoded images by color preview," Opt. Express 14, 5916-5927 (2006). [CrossRef] [PubMed]
- S. N. Savenkov, "Optimization and structuring of the instrument matrix for polarimetric measurements," Opt. Eng. 41, 965-972 (2002). [CrossRef]
- J. Bernardo and A. Smith, Bayesian Theory, (Wiley, 2000).
- A. Gelman, J. Carlin, H. Stern, and D. Rubin, Bayesian data analysis, Second ed., (CRC Press, 2003).
- S. Z. Li, Markov random field modeling in image analysis, Second ed., (Springer, 2001).
- A. Gray, J. Kay, and D. Titterington, "An empirical study of the simulation of various models used for images," IEEE Trans. PAMI, 16, 507-513 (1994). [CrossRef]
- A. Dunmur and D. Titterington, "Computational Bayesian analysis of hidden Markov mesh models," IEEE PAMI. 19, 1296-1300 (1997). [CrossRef]
- S. Richardson and P. Green, "On Bayesian analysis of mixtures with an unknown number of components (with discussion)," J. R. Stat. Soc. Ser. B. 59, 731-792 (1997). [CrossRef]
- J. Besag, "On the statistical analysis of dirty pictures (with discussion)," J. R. Stat. Soc. B 48, 259-302 (1986).

## Cited By |
Alert me when this paper is cited |

OSA is able to provide readers links to articles that cite this paper by participating in CrossRef's Cited-By Linking service. CrossRef includes content from more than 3000 publishers and societies. In addition to listing OSA journal articles that cite this paper, citing articles from other participating publishers will also be listed.

« Previous Article | Next Article »

OSA is a member of CrossRef.