SCATMECH > Classes and Functions >
Polarization > MuellerMatrix
The class MuellerMatrix represents a Mueller matrix
in the Stokes-Mueller representation of polarization
states. Arithmetic operations between Mueller matrices and
Stokes vectors, various properties of a Mueller matrix, and
transformation operations are defined.
Include file:
#include "mueller.h"
See also:
SCATMECH Home,
Conventions, StokesVector, JonesMatrix, JonesVector
C.F. Bohren and D.R. Huffman, Absorption and Scattering of Light by Small Particles, (Wiley, New York, 1983).
H.C. van de Hulst, Light Scattering by Small Particles, (Dover, New York, 1981).
R.A. Chipman, "Polarimetry." in Handbook of Optics, (McGraw-Hill, New York, 1995).
C.R. Givens and A.B. Kostinski, "A simple necessary and sufficient condition on physically realizable Mueller matrices," J. Mod. Opt. 40, 471 (1993).
B.N. Simon, et al., "A complete characterization of pre-Mueller and Mueller matrices in polarization optics," J. Opt. Soc. Am. A 27(2), 188-199 (2010).
S.-Y. Lu and R.A. Chipman, "Interpretation of Mueller matrices based on polar decomposition," J. Opt. Soc. Am. A 13(5), 1106-1113 (1996).
Definition of public elements:
class MuellerMatrix {
public:
MuellerMatrix();
MuellerMatrix(const JonesMatrix& j);
MuellerMatrix(const MuellerMatrix& x);
MuellerMatrix& operator=(const MuellerMatrix& x);
double* operator[](int i);
const double* operator[](int i) const;
double Tmax() const;
double Tmin() const;
double diattenuation() const;
double linear_diattenuation() const;
double depolarization_index() const;
double polarization_dependent_loss() const;
double polarizance() const;
double extinction_ratio() const;
bool valid() const;
bool physically_valid() const;
MuellerMatrix rotate(double angle) const;
MuellerMatrix parity() const;
MuellerMatrix transpose() const;
MuellerMatrix inverse() const;
MuellerMatrix log() const;
MuellerMatrix exp() const;
MuellerMatrix normalized() const;
void Cloude_Decomposition(MuellerMatrix& M1,MuellerMatrix& M2,MuellerMatrix& M3, MuellerMatrix& M4) const;
MuellerMatrix Closest_NonDepolarizing() const;
void Lu_Chipman_Decomposition(MuellerMatrix& depolarizer,MuellerMatrix& retarder,MuellerMatrix& diattenuator) const;
MuellerMatrix operator*(const MuellerMatrix& matrix) const;
MuellerMatrix& operator*=(const MuellerMatrix& matrix);
MuellerMatrix operator*(double d) const;
friend MuellerMatrix operator*(double d,const MuellerMatrix &v);
StokesVector operator*(const StokesVector& s) const;
MuellerMatrix& operator*=(double d);
MuellerMatrix operator/(double d) const;
MuellerMatrix& operator/=(double d);
MuellerMatrix operator+(const MuellerMatrix& a) const;
MuellerMatrix operator-(const MuellerMatrix& a) const;
MuellerMatrix operator-() const;
MuellerMatrix& operator+=(MuellerMatrix& matrix);
MuellerMatrix& operator-=(MuellerMatrix& matrix);
};
ostream& operator<<(ostream& os,const MuellerMatrix& mm);
istream& operator>>(istream& is, MuellerMatrix& mm);
The default constructor for a MuellerMatrix does not
initialize any of the elements.
Top of Page
Constructor that converts a JonesMatrix to a
MuellerMatrix, using the algorithm described by Bohren and Huffman.
Example:
JonesMatrix j;
MuellerMatrix m(j); // Conversion from a Jones matrix to a Mueller matrix.
Top of Page
The copy constructor.
Example:
MuellerMatrix a;
MuellerMatrix b(a);
Top of Page
The assignment operator.
Example:
MuellerMatrix a;
MuellerMatrix b=a;
Top of Page
Operators to access specific elements of the
MuellerMatrix.
Example:
MuellerMatrix a;
double b=a[0][0]; // Access the 00 term of the Mueller matrix
Top of Page
Functions to return the maximum and minimum
transmissivity, the diattenuation, the linear
diattenuation, the polarization-dependent loss, the
polarizance, the extinction ratio, and the depolarization index, respectively, for the given
Mueller matrix. These terms are defined in Chipman.
Example:
MuellerMatrix a;
cout << "Tmax = " << a.Tmax() << endl
<< "Tmin = " << a.Tmin() << endl
<< "diattenuation = " << a.diattenuation() << endl
<< "lineardiattenuation = " << a.lineardiattenuation() << endl;
Top of Page
The function valid() uses the criterion described by
Givens and Kostinski to test if the
matrix is a mapping from the space of valid Stokes vectors
to itself. The function physically_valid() tests the
criterion that a Mueller matrix must be a convex sum (i.e.,
a sum of positive Mueller matrices) of Jones-Mueller
matrices, as described by Simon, et
al.. The function physically_valid() is a more
stringent requirement that a physical Mueller matrix must
satisfy.
Top of Page
Function that returns a MuellerMatrix rotated by
angle in the clockwise direction, looking into the
beam.
Example:
MuellerMatrix a,b;
b = a.rotate(45*PI/180.);
Top of Page
Function that returns a MuellerMatrix that has
undergone a mirror operation.
Example:
MuellerMatrix a,b;
b = a.parity();
Top of Page
Function that returns the transpose of a
MuellerMatrix.
Example:
MuellerMatrix a,b;
b = a.transpose();
Top of Page
Function that returns the inverse of a
MuellerMatrix.
Example:
MuellerMatrix a,b;
b = a.inverse();
Top of Page
Functions that return the matrix logarithm and exponent of a MuellerMatrix.
Top of Page
Function that returns the Mueller matrix normalized by the 00 element.
Top of Page
Function that decomposes the Mueller matrix into the sum of four non-depolarizing Mueller matrices, as described by Cloude.
Example:
MuellerMatrix m = MuellerMatrix(JonesMatrix(1,1,0,0))+MuellerMatrix(JonesMatrix(1,-2,0,0));
MuellerMatrix m1,m2,m3,m4;
m.Cloude_Decomposition(m1,m2,m3,m4);
cout << m - m1 - m2 - m3 - m4 << endl; // should be close to the zero matrix
Top of Page
Function that returns the largest of the matrices returned by Cloude_Decomposition.
MuellerMatrix m;
cout << m.Closest_NonDepolarizing() << endl;
Top of Page
Function that decomposes the Mueller matrix into the product of a depolarizer, a retarder, and a diattenuator,
according to Lu and Chipman.
Example:
MuellerMatrix a,depolarizer,retarder,diattenuator;
a.Lu_Chipman_Decomposition(depolarizer,retarder,diattenuator);
cout << a - depolarizer*retarder*diattenuator; // should be close to the zero matrix
Top of Page
MuellerMatrix operator*(const MuellerMatrix& matrix) const
MuellerMatrix& operator*=(const MuellerMatrix& matrix)
Binary multiplication between two MuellerMatrix
objects. If matrix2 is followed by matrix1,
then matrix1*matrix2 returns the net
MuellerMatrix. This operator is not commutative.
Example:
MuellerMatrix a,b,c;
c = a*b;
a *= b; // Same as a = a*b;
Top of Page
Multiplication and division of a MuellerMatrix by
a scalar.
Example:
double scalar;
MuellerMatrix matrix;
matrix = matrix * scalar;
matrix = scalar * matrix;
matrix *= scalar;
matrix = matrix / scalar;
matrix = scalar / matrix; // Not allowed.
matrix /= scalar;
Top of Page
Left multiplication of a StokesVector by a
MuellerMatrix.
Example:
MuellerMatrix matrix;
StokesVector vector;
vector = matrix*vector;
Top of Page
Addition and subtraction of MuellerMatrix objects.
Example:
MuellerMatrix matrix1,matrix2,matrix3;
matrix3=matrix1+matrix2;
matrix3=matrix1-matrix2;
matrix1+=matrix2;
matrix1-=matrix2;
Top of Page
Negation of a MuellerMatrix.
Example:
MuellerMatrix matrix;
matrix=-matrix;
Top of Page
Operator to send a MuellerMatrix to an
ostream. The result will contain four Stokes vectors, representing each of the four rows, separated by commas, and surrounded by
parentheses, e.g., ((1,0,0,0),(0,1,0,0),(0,0,1,0),(0,0,0,1)).
Example:
MuellerMatrix M;
cout << "M = " << M << endl;
Top of Page
Reads a MuellerMatrix from an input stream. The input string must contain
four Stokes vectors, representing each of the four rows, separated by commas, and surrounded by
parentheses, e.g., ((1,0,0,0),(0,1,0,0),(0,0,1,0),(0,0,0,1)).
Top of Page
For More Information
SCATMECH Technical Information and Questions
Sensor Science Division Home Page
Sensor Science Division Inquiries
Website Comments
Current SCATMECH version: 7.22 (April 2021)
|