%*****************************************************************************************
% Useful plot as it shows relative contributions positive or negative
%*****************************************************************************************
% The leading diagonal has the highest values and always positive.
% Typically off-diagonal values are much smaller or even negative.
%
% Stiffness Voigt Matrix in 3D Histogram or 2D Histogram
% in colour or black and white  - single crystal or aggregate stiffness
% tensors
%
% Tested with mtex-5.3beta3 and MATLAB 2016a
%
% Remember that the function
%**************************************************************************
% Stiffness_Voigt_Matrix_in_3D_Histo(Stiffness,Plot_title,col_map,iview)
%**************************************************************************
% should be in your current directory
%
% David Mainprice 2/05/2019
%
% Options for plots
% 1) plot    (col_map = 0) white2black      and (iview = 1) tilted 3D
% 2) plot    (col_map = 1) WhiteJetColorMap and (iview = 1) tilted 3D
% 3) plot    (col_map = 0) white2black      and (iview = 0) normal plane 2D
% 4) plot    (col_map = 1) WhiteJetColorMap and (iview = 0) normal plane 2D
%
clear all
close all
clc
%%
%************************************************
% Define elastic stiffness tensor (GPa)
%************************************************
%
% Last checked 28/1/2018 correct symbol '-3m'
%
% Reference Elastic constants
% Reference
% Ogi, H., Ohmori, T. Nakamura, N. and Hirao M., 2006. Elastic, anelastic and piezoelectric
% coefficients of alpha-quartz determined by resonance ultrasound spectroscopy.
% Journal of Applied Physics, 100, 053511.
%
% Reference Crystal Structure
% alpha quartz Trigonal SG P3221 25 C d=2.648 g/cm3
%
% Define density (g/cm3)
rho= 2.65;
% Space groups: the enantiomorphic space groups, P3121 (No. 152) and P3221 (No. 154)
%
% Define Cartesian Tensor crystal symmetry & frame
cs_Tensor = crystalSymmetry('-3m',[  4.9134  4.9134  5.4052],...
 [  90.00  90.00 120.00]*degree,'x||a','z||c',...
 'mineral','Alpha-Quartz 2006');
%
% elastic Cij stiffness tensor (GPa) as matrix M
 M =.... 
 [[   86.76    6.868   11.85  -18.02    0.00    0.00];...
 [     6.868   86.76   11.85  +18.02    0.00    0.00];...
 [    11.85    11.85  105.46    0.00    0.00    0.00];...
 [   -18.02   +18.02    0.00   58.14    0.00    0.00];...
 [     0.00     0.00    0.00    0.00   58.14  -18.02];...
 [     0.00     0.00    0.00    0.00  -18.02   39.946]];
% M as stiffness tensor C with MTEX tensor command
C = stiffnessTensor(M,cs_Tensor,'density',rho);
%% ************************************************************************
%  Plot 1 : Quartz Plot 3D Voigt Stiffness Matrix
%  with options (col_map = 0)=> white2black and (iview = 1)=> tilted view
%**************************************************************************
% Tensor = C = MTEX elastic stiffness tensor (GPa)
% Plot_title = example 'Alpha-quartz stiffness tensor in Voigt matix C(i,j) in GPa'
% col_map = 1 WhiteJetColorMap 0= white2black
% iview =   1 tilted 3D view   0= normal to plane 2D
% Title of plot
Plot_title='Alpha-quartz stiffness tensor in Voigt matix Cij in GPa';
% Option WhiteJetColorMap
col_map = 0;
% Option tilted 3D view
iview = 1;
close all
[Tmax_Pos,Tmin_Pos,Tmax_Neg,Tmin_Neg]=Stiffness_Voigt_Matrix_in_3D_Histo(C,Plot_title,col_map,iview);
%% ************************************************************************
%  Plot 2 : Olivine Plot 3D Voigt Stiffness Matrix
%  with options (col_map = 1)=> WhiteJetColorMap and (iview = 1)=> tilted view
%**************************************************************************
% Tensor = C = MTEX elastic stiffness tensor (GPa)
% Plot_title = example 'Olivine stiffness tensor in Voigt matix C(i,j) in GPa'
% col_map = 1 WhiteJetColorMap 0= white2black
% iview =   1 tilted 3D view   0= normal to plane 2D
% Title of plot
Plot_title='Alpha-quartz stiffness tensor in Voigt matix Cij in GPa';
% Option WhiteJetColorMap
col_map = 1;
% Option tilted 3D view
iview = 1;
close all
[Tmax_Pos,Tmin_Pos,Tmax_Neg,Tmin_Neg]=Stiffness_Voigt_Matrix_in_3D_Histo(C,Plot_title,col_map,iview);
%% ************************************************************************
%  Plot 3 : Olivine Plot 3D Voigt Stiffness Matrix
%  with options (col_map = 0)=> white2black and (iview = 0)=> 2D view
%**************************************************************************
% Tensor = C = MTEX elastic stiffness tensor (GPa)
% Plot_title = example 'Olivine stiffness tensor in Voigt matix C(i,j) in GPa'
% col_map = 1 WhiteJetColorMap 0= white2black
% iview =   1 tilted 3D view   0= normal to plane 2D
% Title of plot
Plot_title='Alpha-quartz stiffness tensor in Voigt matix Cij in GPa';
% Option white2black
col_map = 0;
% Option normal to plane 2D
iview = 0;
close all
[Tmax_Pos,Tmin_Pos,Tmax_Neg,Tmin_Neg]=Stiffness_Voigt_Matrix_in_3D_Histo(C,Plot_title,col_map,iview);
%% ************************************************************************
%  Olivine Plot 3D Voigt Stiffness Matrix
%  with options (col_map = 1)=> WhiteJetColorMap and (iview = 0)=> 2D view
%**************************************************************************
% Tensor = C = MTEX elastic stiffness tensor (GPa)
% Plot_title = example 'Olivine stiffness tensor in Voigt matix C(i,j) in GPa'
% col_map = 1 WhiteJetColorMap 0= white2black
% iview =   1 tilted 3D view   0= normal to plane 2D
% Title of plot
Plot_title='Alpha-quartz stiffness tensor in Voigt matix Cij in GPa';
% Option WhiteJetColorMap
col_map = 1;
% Option normal to plane 2D
iview = 0;
close all
[Tmax_Pos,Tmin_Pos,Tmax_Neg,Tmin_Neg]=Stiffness_Voigt_Matrix_in_3D_Histo(C,Plot_title,col_map,iview);
