Analytical calculation of the axis angle 2V from extinction measurements on the spindle stage
Florian Dufey



Analytical calculation of the axis angle 2V from extinction measurements on the spindle stage
    F. Dufey
    Received: date / Accepted: date
=============================================================================================


	A concise derivation of the "Joel equations", which allow for the determination of the axis angle 2V from measurements of extinction directions on a spindle stage, is provided starting from the wave-equation. 
	Only analytic methods and no geometric arguments referring to stereographic projections are invoked. 
	For error free data, the resulting equations allow for a closed form solution.  
	If angle data with measurement error are to be used, a maximum likelihood estimation methodology is proposed, which may be solved using e.g. iterative reweighting. 
	The method is tested with and compared to  published data. 



§ INTRODUCTION


Since its introduction in the beginning 1950s,  the one-circle or spindle stage goniometer for the determination of the angle 2V between the optical axis has become one of the standard methods for this purpose. 
While the determination of the axis angle relied in the first years on rather complicated geometric constructions in the stereographic projection, the technique gained further acceptance when it became possible to solve the underlying equation on a computer <cit.>. 
The derivation of these equations is not easy to comprehend today, as it involves a mix of geometric and analytic constructs, like the equivibration curve,  and is scattered several articles <cit.>. 
It is one goal of this notice to give a concise derivation of these equations in a purely analytical fashion starting from the Fresnel equation. 
The second result to be presented is that these equations principally allow for a closed form solution for the components of the inverse dielectric tensor Φ, and hence also for 2V, when error free extinction positions are available and to derive a maximum likelihood estimator for the components of the inverse dielectric tensor, when data prone to statistical error are used. 



§ DERIVATION OF A FORMULA FOR THE ANGLE 2V


The tensor Φ is the inverse of the optical dielectric tensor,

    Φ=ϵ^-1.

It can be decomposed as


    Φ=A_2 1 -1/2 (A_1-A_3) φ,

Here, A_i=1/ϵ_i=1/n_i^2 (i ∈{1,2,3}) where ϵ_i are the principal components of the dielectric  tensor and n_i the main indices of refraction. 
The A_i are ordered as A_1>A_2>A_3 (uniaxial crystals are not considered here).
The tensor φ=a_1a_2^T+a_2a_1^T is made up from the two vectors a_1 and a_2 of unit length which are perpendicular to the two circular sections of the indicatrix and parallel to the optical axes. 
We are interested in calculating  the tensor Φ and the angle 2V between the vectors a_1 and a_2, but will not use these vectors any further in the course of the following calculations. 


According to <cit.>, eq. 21, the directions of extinction coincide with the dielectric displacement vector D, which fulfills the wave equation

    (1-n̂n̂^T)Φ D =1/n^2D,

where n̂=k/|k| is the unit vector parallel to the wavevector k  of the wave. 
The resulting equation for the refractive index n is Fresnel's equation, but we are more interested in the measured direction of 
D, which is always perpendicular to n̂. 


So if we hold n̂ and D fixed and choose the coordinate axes as e_y=n̂ and D || e_z, (the direction e_x is then parallel to the second possible polarisation.) the equation for the z-component is trivially fulfilled and the equation for the y-component can always be fulfilled with a suitable choice of n, leaving the only restriction, that the left hand side does not introduce a y-component for D, or,

    Φ_xy=e_y^TΦ e_x=0.

Now for any two vectors q and q' which fulfill (q+q')|| e_x and (q-q')||e_y, we get

    (q+q')Φ(q-q')=0,

which are the equations given by Joel. 

Before solving eq. <ref>, we note that the tensor Φ is a symmetric tensor which can be parametrized by 6 constants. 
It is clear that the direction of D cannot depend on the isotropic part of Φ, 

    Φ_iso= 1 (Φ_xx+Φ_yy+Φ_zz)/3

 and is also invariant under a rescaling Φ→ cΦ with an arbitrary constant c.
Therefore, we can hope at best to recover from measurements on the spindle stage the anisotropic part 

    Φ_aniso=Φ -Φ_iso
 
which depends on 5 parameters of which one may be chosen at will to fix the scaling. 
An explicit expression for Φ_aniso is 

    Φ_aniso=	([  d_x^2-y^2-d_z^2             d_xy             d_xz;             d_xy -d_x^2-y^2-d_z^2             d_yz;             d_xz             d_yz           2d_z^2 ]),

with 


    d_z^2   =   (2Φ_zz-Φ_xx-Φ_yy)/6,
    
    	d_x^2-y^2   =   (Φ_xx-Φ_yy)/2,
    
    	d_xy   =   Φ_xy,
    
    	d_xz   =   Φ_xz,
    
    	d_yz   =   Φ_yz.



Various orientations parameterized by the spindle angle S have to be analysed 
which results in the substitution. 

    Φ→Φ̃(S, E_S)=R_y(E_S)R_z(S) Φ R^T_z(S)R^T_y(E_S).

Here, R_y and R_z are the usual rotation matrices around the y- and z-axes, respectively,

    R_y=( [  cos E_S        0 -sin E_S;        0        1        0;  sin E_S        0  cos E_S ]),
 
and 

    R_z=( [  cos S  sin S      0; -sin S  cos S      0;      0      0      1 ]).
 
Specifically, z is chosen to be the spindle axis and y the direction perpendicular to the rotation desk of the microscope with S and E_S being the respective rotation angles. 


Substitution into eq. <ref> yields an implicit equation for the equivibration curve E_S(S), 


    -3d_z^2tan 2E_S+2d_xzcos S+ 2d_yzsin S+
    + d_x^2-y^2cos 2 Stan 2E_S+d_xysin 2Stan 2E_S=0.


Hence one of the  d-terms may be fixed advantageously, e.g. d_z^2=1. 
Then, if measurements are made for 4 angles E, this linear system for the other d_i can be solved exactly. 
However, if more than 4 angles are measured, then some regression technique has to be employed.
The programs available <cit.> minimize an ad hoc sum of squared Joel equations. 
This arbitrarily assigns weights (=1) to all observations. 
To put the regression on a statistical basis, it seems more appropriate to start from a reasonable error model for the dependent variables E_S. 
With the expectation value of the angle E_S being 1/2 arctanV/W where V=2d_xzcos S+ 2d_yzsin S and W=3d_z^2-d_x^2-y^2cos 2 S-d_xysin 2S, 
we assume the measurement error of the angles to follow a von Mises distribution <cit.>, 

    f(2E_S)= exp(±κcos(2E_S-arctan(V/W)))/2π I_0(κ)

where κ is a dispersion parameter and I_0 is a modified Bessel function.
The sign of the exponent depends on which of the two solutions E_S is chosen with points on the same extinction curve giving rise to the same sign. 
The parameter κ may be estimated as the empirical variance of the observed angles from the calculated values. 
The likelihood to be maximized with respect to the parameters d_j, j∈{x^2-y^2,xy, xz, yz} given measured values E_Si from one of the two extinction curves with i∈ [1,n]  is then

    L({d_j})=±κ∑_i=1^N W_icos 2E_Si +V_i sin 2E_Si/√(V_i^2+W_i^2)

The optimum is attained when the four equations 

    0=∂ L/∂ d_j=±κ∑_i (W_i sin 2E_Si-V_i cos2E_Si) (V_i^2+W_i^2)^-3/2∂ V_iW_i/∂ d_j

are fulfilled, which can be seen to be weighed sums of the original Joel equations. 
Hence, iterative reweighting seems to be a promising and numerically stable option to obtain maximum likelihood estimates of the parameters d_j. 
The second derivatives of L with respect to the d_j yield the information matrix from which estimates of the variances of the d_j may be obtained. 

The eigenvectors of Φ_aniso, constructed from the estimates d_j,  are the directions of the main axes of the indicatrix, 
corresponding to the eigenvalues μ_h, μ_m and μ_l, ordered by decreasing value. 
From the latter, the axis angle 2V can be calculated <cit.>,

    tan^2V=μ_h-μ_m/μ_m-μ_l.

It is obvious that this result is not influenced by neither a  constant scaling nor a shift of all eigenvalues. 



§ COMPARISON WITH PUBLISHED DATA


To exemplify the feasibility of the method, angles 2V calculated from the measured angles E_S taken from  the  articles <cit.>
are being compared to values 2V reported there: 

<cit.>, Fig. 5: 64.09^∘ (63.66^∘)

<cit.>, Fig. 2: 76.88^∘ (77.57^∘)

<cit.>, Fig. 4:   48.98^∘ (49.36^∘)

The numerical optimization of the likelihood eq. <ref> was performed using the optmodel procedure of  <cit.>.  

Expected values E_S calculated with these parameters compare very well with the reported ones with the mean squared deviations from the measured E_S being even slightly smaller. 



apalike
