function [Bmat, Amat] = gtCrystHKL2CartesianMatrix(lp)
%
%     [Bmat, Amat] = gtCrystHKL2CartesianMatrix(lp)
%     ---------------------------------------------
%     Returns a 3x3 matrix that can be used for coordinate transform from 
%     a signed HKL into a Cartesian vector in real space and for any lattice 
%     and unit cell.
%
%     The Cartesian (orthogonal) basis (X,Y,Z) is defined such that:
%       - lattice vectors of the unit cell are a1, a2, a3 (not orthogonal)
%       - X is parallel with a1
%       - Y lies in the a1-a2 plane
%       - Z is cross(X,Y)
% 
%     The definition does not correspond to the one in Poulsen's 3DXRD
%     book, Chapter 3, but it is an equivalent representation based on the
%     real space lattice.
%     This function essentially computes the B matrix in the equation 
%     Gc = B*Ghkl  (used with column vectors!)
%
%     INPUT:
%       lp = unit cell lattice parameters in real space 
%            [a b c alpha beta gamma] (angles are in degrees)
%
%     OUTPUT:
%       Bmat = the transformation matrix for column vectors (!)
%       Amat = the A matrix too (real space) for dealing with UVW lattice
%              directions.
%

    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    %%% Lattice vectors in real space, in a Cartesian basis
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    cos_alpha = cosd(lp(4));
    cos_beta  = cosd(lp(5));
    cos_gamma = cosd(lp(6));
    sin_gamma = sind(lp(6));

    % Calculate 'a', and 'b' lattice vectors in X, Y, Z coordinates:
    a = lp(1) * [1,         0,         0];
    b = lp(2) * [cos_gamma, sin_gamma, 0];

    % Calculate 'c1' and 'c2' components:
    c1 = lp(3) * cos_beta;
    c2 = lp(3) * (cos_alpha - cos_gamma * cos_beta) / sin_gamma;
    % Calculate 'c3' component as the result:
    c = [c1,  c2, sqrt(([lp(3), c1, c2] .^ 2) * [1; -1; -1])];

    % We keep it transposed for the computation of B, and then transpose it
    % back again in the end.
    Amat = [a; b; c];

    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    %%% Reciprocal lattice vectors in the Cartesian basis
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    % Reciprocal lattice vectors by Crystallographers' definition:
    %  a* = cross(b, c) / cell_vol
    %  b* = cross(c, a) / cell_vol
    %  c* = cross(a, b) / cell_vol
    cross_prods = gtMathsCross(Amat, Amat([2 3 1], :))';

    % Volume of the direct lattice unit cell: V = (a dot (b cross c))
    cell_vol = a * cross_prods(:, 2);

    % B matrix finally
    Bmat = cross_prods(:, [2 3 1]) / cell_vol;

    % A matrix (for output)
    Amat = Amat';

end % of function