Generalised March-Dollase preferred orientation corrections

This is a collection of macros that are used to model preferred orientation using the March-Dollase model, where the angle between the preferred orientation vector and the diffraction vector is able to change.

In these macros, you can deal explicitly with symmetric reflection, capillary transmission, and asymmetric reflection.

The maths for this was taken from [1].

Contributor: Matthew Rowles

[1] Ida, Takeshi. 2013. 'Effect of Preferred Orientation in Synchrotron X-Ray Powder Diffraction'. Annual Report of the Advanced Ceramics Research Center: Nagoya Institute of Technology..

The March-Dollase pole density profile is given by

P_p^*(r,rho) = (r² cos²rho + (sin²rho/r))^-3/2

where r is the preferred orientation factor, and rho is the polar angle between the preferred orientation direction, p^*, and the specimen direction, s^*. The specimen is assumed to have rotational symmetry about s^*.

To find the factor by which a the intensity of a reflection of any given direction should multiplied, we use the following equations:

f_d^*(r,alpha,Delta) = (1/(2*Pi)) * Int(g(r,alpha,Delta,phi), phi=0..2*Pi)

g(r,alpha,Delta,phi) = P_p^*(r,rho)

cos(rho) = cos(alpha)*cos(Delta) - sin(alpha)*sin(Delta)*sin(phi)

where alpha is the angle between the diffraction vector, d^*, of interest and the preferred orienation vector, and Delta is the angle between d^* and s^*.

Delta is 0 for symmetric reflection, Pi/2 for capillary transmission, and Abs(Th - omega) for asymmetric reflection, where the incident beam angle is given by omega.

The integral can be expressed as a summation:

f_d^*(r,alpha,Delta) = (1/N) * Sum(g(r,alpha,Delta,(j+(1/2))*Pi/N), j=0..(N-1))

however, convergence is slow if r » 1. N = 16 is usually sufficient.

The macros to use are the following:

macro PO_symmetric_reflection(r_c, r_v, hkl)
{
	generalised_PO_eqn(r_c, r_v, , hkl, , 0, 16)
}
 
macro PO_capillary_transmission(r_c, r_v, hkl)
{
	generalised_PO_eqn(r_c, r_v, , hkl, , (Pi/2), 16)
}
 
macro PO_asymmetric_reflection(r_c, r_v,omega_c,omega_v, hkl)
{
	#m_argu omega_c
	If_Prm_Eqn_Rpt(omega_c, omega_v, min 0.0001 max 90) 'incident angle in degrees'
	prm #m_unique Delta = Abs(Th - ((CeV(omega_c,omega_v))*Deg)) ;
	generalised_PO_eqn(r_c, r_v, , hkl, , Delta, 16)
}

An example of their use would be

PO_asymmetric_reflection(r, 0.6, !omega, 7.5, 0 0 1)

to refine the preferred orientation factor for the (001) plane with a fixed incident beam angle of 7.5°.

The required helper macros are as follows, and implement equation 10 (aka I) in [1]:

macro & cosrho(& alpha, & Delta, & phi)
{
	Cos(alpha)*Cos(Delta) - Sin(alpha)*Sin(Delta)*Sin(phi) 
}
 
macro & PO_P(& r, & alpha, & Delta, & phi)
{
	(r^2 * cosrho(alpha,Delta,phi)^2 + (1-cosrho(alpha,Delta,phi)^2)/r)^(-3/2)
}
 
macro & PO_f(& r, & alpha, & Delta, & N)
{
	'hard coded this for N==16 as the Sum formalism is not working '
	'(1/N)*Sum( PO_P(r,alpha,Delta,(j+(1/2))*Pi/N) , j=0,j<=N-1,j=j+1)'
	 (1/16)*
	 (PO_P(r,alpha,Delta,( 0+(1/2))*Pi/16) +
	  PO_P(r,alpha,Delta,( 1+(1/2))*Pi/16) +
	  PO_P(r,alpha,Delta,( 2+(1/2))*Pi/16) +
	  PO_P(r,alpha,Delta,( 3+(1/2))*Pi/16) +
	  PO_P(r,alpha,Delta,( 4+(1/2))*Pi/16) +
	  PO_P(r,alpha,Delta,( 5+(1/2))*Pi/16) +
	  PO_P(r,alpha,Delta,( 6+(1/2))*Pi/16) +
	  PO_P(r,alpha,Delta,( 7+(1/2))*Pi/16) +
	  PO_P(r,alpha,Delta,( 8+(1/2))*Pi/16) +
	  PO_P(r,alpha,Delta,( 9+(1/2))*Pi/16) +
	  PO_P(r,alpha,Delta,(10+(1/2))*Pi/16) +
	  PO_P(r,alpha,Delta,(11+(1/2))*Pi/16) +
	  PO_P(r,alpha,Delta,(12+(1/2))*Pi/16) +
	  PO_P(r,alpha,Delta,(13+(1/2))*Pi/16) +
	  PO_P(r,alpha,Delta,(14+(1/2))*Pi/16) +
	  PO_P(r,alpha,Delta,(15+(1/2))*Pi/16))	
}
 
 
'angles in radians'
macro generalised_PO_eqn(r_c, r_v, alpha, hkl, Delta_c, Delta_v, N) 
{
   #m_argu r_c
   #m_argu Delta_c
   If_Prm_Eqn_Rpt(r_c, r_v, min 0.0001 max = 2 Val + .5;) 'preferred orientation parameter'
   If_Prm_Eqn_Rpt(Delta_c, Delta_v, min 0.0001 max =Pi;) 'angle in radians:  ==0 symmetric reflection, ==Pi/2 capillary transmission, ==Abs(Th-Omega) asymmetric reflection'
 
   'The following code modified from S11.5.4 in the TechRef & PO macro in topas.inc'
   #m_ifarg alpha ""
      #m_unique_not_refine alpha
   #m_endif
 
   str_hkl_angle alpha hkl 'angle in radians between PO direction and current diffraction vector'
 
   scale_pks = Multiplicities_Sum(  PO_f(CeV(r_c,r_v), alpha, CeV(Delta_c,Delta_v), N)  );  
}