%==========================================================================
%==========================================================================
%
%              /$$$$$$$$ /$$$$$$$   /$$$$$$  /$$$$$$$                       
%             | $$_____/| $$__  $$ /$$__  $$| $$__  $$                      
%             | $$      | $$  \ $$| $$  \__/| $$  \ $$                      
%             | $$$$$   | $$$$$$$ |  $$$$$$ | $$  | $$                      
%             | $$__/   | $$__  $$ \____  $$| $$  | $$                      
%             | $$      | $$  \ $$ /$$  \ $$| $$  | $$                      
%             | $$$$$$$$| $$$$$$$/|  $$$$$$/| $$$$$$$/                      
%             |________/|_______/  \______/ |_______/                       
%                                                                           
%                                                                           
%                                                                           
%                                                    /$$                    
%                                                   |__/                    
%   /$$$$$$  /$$   /$$  /$$$$$$   /$$$$$$   /$$$$$$$ /$$  /$$$$$$$  /$$$$$$ 
%  /$$__  $$|  $$ /$$/ /$$__  $$ /$$__  $$ /$$_____/| $$ /$$_____/ /$$__  $$
% | $$$$$$$$ \  $$$$/ | $$$$$$$$| $$  \__/| $$      | $$|  $$$$$$ | $$$$$$$$
% | $$_____/  >$$  $$ | $$_____/| $$      | $$      | $$ \____  $$| $$_____/
% |  $$$$$$$ /$$/\  $$|  $$$$$$$| $$      |  $$$$$$$| $$ /$$$$$$$/|  $$$$$$$
%  \_______/|__/  \__/ \_______/|__/       \_______/|__/|_______/  \_______/
%                                                                           
%                                                                           
%                                                                           
%                     /$$$$$$   /$$$$$$   /$$$$$$   /$$$$$$                 
%                    /$$__  $$ /$$$_  $$ /$$__  $$ /$$$_  $$                
%                   |__/  \ $$| $$$$\ $$|__/  \ $$| $$$$\ $$                
%                     /$$$$$$/| $$ $$ $$  /$$$$$$/| $$ $$ $$                
%                    /$$____/ | $$\ $$$$ /$$____/ | $$\ $$$$                
%                   | $$      | $$ \ $$$| $$      | $$ \ $$$                
%                   | $$$$$$$$|  $$$$$$/| $$$$$$$$|  $$$$$$/                
%                   |________/ \______/ |________/ \______/                 
%                                                                           
%                                                                   
%                                                                           
%==========================================================================
%==========================================================================
 
% The aim of this exercise is to get familiar with mtex and look at some ebsd
% data

% Some general Notes:
% - usually you can find help by typing `help` command
%   if a command with an identical name might be available for different
%   data types try e.g. `help orientation/calcMDF` or `help odf/calcMDF` or
%   help orientation/mean` (`help mean` will be Matlabs number-valued mean)
% - using TAB-completion makes many tasks faster, prevents typos and helps
%   to explore functions and properties available to an object
%   e.g. ebsd.[hit the tabulator key] will present a list of features
%  available with the ebsd object
%% install/startup mtex (in case it isn't already)
addpath /Path/to/mtex 
startup_mtex

%% 0) load the data -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
% we can use the import_wizard and generate a script which imports the ebsd
% data. Use the file ebsd1.ctf


import_wizard
% selcte EBSD and click Next ...
% when asked for the reference frame, x should be east and y south when
% using convertEuler2SpatialReferenceFrame


%% 1) explore the ebsd data *-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
%  What can you tell about the ebsd data?
% hint: if you use ebsd.[TAB] you will see some properties and functions
% available for this object


% %########################################################################
% a) How many indexed pixels does the map have?
% %########################################################################


% %########################################################################
% b) Which properties does the ebsd map have?
% %########################################################################

ebsd.prop

%       mad = mean angular deviation:misfit between detected bands and
%       simulation (in degree, high mad = 'bad' solution)
%       bc = band contrast: ratio of actual intensity in hough space vs
%       average intensity of kikuchi bands (high value = high contrast)
%       bs =band slope : higest intensity gradient at the morgin of kikuchi
%       bands (high value = sharp band)


% %########################################################################
% c) How many percent of the map have pixels with a mad < 0.8 degree? 
% %########################################################################


% %########################################################################
% d) How many pixels are indexed within the first 100 points?
% %########################################################################


% %########################################################################
% e) How large is the map? How many pixels does it have in x,y direction?
% %########################################################################

ebsd.gridify
ebsd.extend

% %########################################################################
% f) What is the stepsize of the map?
% %########################################################################

ebsd.unitCell

%% % 2) look at the map -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*

% use the plot command to plot almost anything in mtex
plot(ebsd)
% maybe not so spectacular at first sight, it will just plot the phase

% mtex will always try to honor the convention
% plot(where,what)
% in case of EBSD data:
% "where" are x,y coordiantes of each measurement point
% "what" can be anything as long as it has the same size/length as the first argument.

% %########################################################################
% Question:
% a) Colorcode the band contrast in a grey-scale color-coding
% hint1: when plotting (where,what) values, both need to have the same
% length
% hint2: to overlay two plots, use 'hold on' and 'hold off'
% %########################################################################

plot(ebsd,ebsd.prop.bc)
mtexColorMap white2black
mtexColorbar(gca,'Title','band contrast')

% %########################################################################
% b) now add to the plot above the non-indexed pixels in red.
% %########################################################################

hold on  % this keeps the present figure axis
plot(ebsd('notIndexed'),'FaceColor',[1 0.2 0],'FaceAlpha',0.6)
hold off 



%% 2) inspect orientations -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*

% %########################################################################
% a) access the first 100 Euler angles
% Hint: orientations are stored in ebsd(phase).orientations
% %########################################################################

ebsd('q').orientations
ebsd('q').orientations(1:100).Euler

% %########################################################################
% b) plot a polefigure (point plot)
%    - define crystal directions (@Miller)
%    - either use plotPDF or remember that orientiation*Miller -> specimen direction
% %########################################################################


cs = ebsd('q').CS;
hs= {... % 
 Miller(0,0,0,1,cs,'UVTW'),... %c
 Miller(1,1,-2,0,cs,'UVTW'),... %a
 Miller(1,0,-1,0,cs,'hkil'),... %m
 Miller(0,1,-1,1,cs,'hkil'),... %z
 Miller(1,0,-1,1,cs,'hkil')... %r
 };
plotPDF(ebsd('q').orientations,hs)

% %########################################################################
% c) plot an inverse pole figure (point plot) with respect to (0,0,1), (-5,1,0)
%     and (2,-3,0)
%    - define a specimen direction (@vector3d)
%    - either use plotIPDF or remember that orientation*vector3d -> crystal direction
% %########################################################################

refv = [zvector vector3d(-5,1,0) vector3d(2,-3,0)]

plotIPDF(ebsd('q').orientations,refv,'MarkerSize',2)
%%
close all
%% 3) orientation mapping  -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
%  Colorcode the orientations in the map

% %########################################################################
% a) plot three maps, colorcoded by the three Euler angles.  Do you
% understand what you see?
% %########################################################################

eu = ebsd('q').orientations.Euler;
plot(ebsd('q'),eu(:,1)/degree); nextAxis
plot(ebsd('q'),eu(:,2)/degree); nextAxis
plot(ebsd('q'),eu(:,3)/degree)

% %########################################################################
% b) make an inverse pole figure colorcoding with respect to one of the
% directions in 2c) and plot it as a map 
% useful commands: ipfHSVKey the property .inversePoleFigureDirection
% %########################################################################


oM = ipfHSVKey(ebsd('q').CS);
oM.inversePoleFigureDirection = refv(2);
col=oM.orientation2color(ebsd('q').orientations);
plot(ebsd('q'),col)

% %########################################################################
% c) colorcode points in a polefigure and in an inverse polefigure indentical
%    to the colorcoding in the map from b)
% %########################################################################

for i=1:length(refv)
    oM.inversePoleFigureDirection = refv(i);
    col=oM.orientation2color(ebsd('q').orientations);
    figure
    plotPDF(ebsd('q').orientations,col, hs([1 2 4]),'MarkerSize',4)
    nextAxis
    plotIPDF(ebsd('q').orientations, col, refv(i),'MarkerSize',4)
    nextAxis 
    plot(ebsd('q'),col)
end

% %########################################################################
% d) plot a color mapping for the c-axis directions
% useful commands: HSVDirectionKey
% note: the c-axis is just a direction, use .direction2color
% %########################################################################


dM = HSVDirectionKey(specimenSymmetry('-1'))
cdir = ebsd('q').orientations.*hs(1);
plot(ebsd('q'),dM.direction2color(cdir))
nextAxis
plot(dM)

% %########################################################################
% e) plot a Euler angle mapping, an axis-angle mapping and compare with the
%    previous ones and compare. Which is best?
% %########################################################################


eM = BungeColorKey(cs)
plot(ebsd('q'),eM.orientation2color(ebsd('q').orientations))
figure
aM = axisAngleColorKey(cs)
plot(ebsd('q'),aM.orientation2color(ebsd('q').orientations))



%% 4) Selecting data  -----------------------------------------------------
% %########################################################################
% a) Question:
% Do pixels indexed with just few bands (e.g. <=6) show a different 
% polefigure compared to those with more bands (e.g. >9) ? Where are they
% locatd in the map?
% %########################################################################

ebsdQ=ebsd('q')
ind6=ebsdQ.prop.bands <= 6;
ind9=ebsdQ.prop.bands > 9;

figure
plotPDF(ebsdQ(ind6).orientations,ebsdQ(ind6).prop.mad,hs,'MarkerSize',3)
mtexColorbar('Title','mad')
CLim(gcm,[min(ebsd('q').prop.mad) max(ebsd('q').prop.mad)])
figure
plotPDF(ebsdQ(ind9).orientations,ebsdQ(ind9).prop.mad,hs,'MarkerSize',3)
mtexColorbar('Title','mad')
CLim(gcm,[min(ebsd('q').prop.mad) max(ebsd('q').prop.mad)])

figure
plotPDF(ebsdQ(ind6).orientations,hs,'MarkerSize',3)
f=gcm;
f.currentId = 1
hold on
plotPDF(ebsdQ(ind9).orientations,hs,'MarkerSize',3)
hold off
%%
 close all
 
%% manually select data - ebsd area
% %########################################################################
% b) Select a triangular area in the upper left of the map.
% useful commands: EBSD/inpolygon
% %########################################################################

plot(ebsd('q'),ebsd('q').orientations)
ebsd_xy = ginput; % click around a polygon
ebsd_sel = ebsd(inpolygon(ebsd,ebsd_xy))
plot(ebsd_sel('q'),ebsd_sel('q').orientations)

% %########################################################################
% c) Plot the c-axis inclination along a N-S profile in the ebsd map.
% useful commands: EBSD/spatialProfile
% %########################################################################

plot(ebsd('q'),ebsd('q').orientations)
xy = ginput(2)
[ebsdLine,dist] = spatialProfile(ebsd,xy)
% find the indexed points
ind = ebsdLine.phase==1;
cinc = angle(ebsdLine(ind).orientations*Miller(0,0,0,1,cs),zvector,'antipodal')
plot(dist(ind),cinc/degree)
 
%% -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*
% 5) More questions:
% -------------------------------------------------------------------------
% a) The map is taken with respect to a kinematic referce system:
%  - the flow/shear plane is E-W (1,0,0)
%  - the instantaneous shortening direction is (2,-3,0)
%  - the plane of maximum shear stress is (-5,1,0)
%  - the long axis of the finite strain ellipsoid is at 7 degree with in
%    the map plane (rho =7*degree)




% %########################################################################
% a) Question: What can you say about the a-axis alignement in the sample
%  with respect to the kinemtic reference frame? 
%  Colorcode the a-axis direction as a function of the distance from the
%  shear plane, the plane of maximum shear stress and the X direction of
%  finite strain ellipsoid.
%
%  Hint: rotate the vector into crystal coordinates and compute the angle 
%  with the reference crystal direction. This way, you don't have to worry
%  about symmetrically different solutions.
% %########################################################################



% %########################################################################
% b) Plot a spot color mapping , color all orientations with <a> 
%  (i)roughly and (ii) exactly up to 25 degree away from the a-maximum  
% %########################################################################
% Hint: this will give you the direction of the maximum in the <a> axes
% polefigure
cs = ebsd('q').CS
a= Miller(1,1,-2,0,cs)
aF  =calcDensity(ebsd('q').orientations*a,'halfwidth',10*degree)
[~,maxpos] = max(aF);
% %########################################################################

%
% option 1
om =ipfSpotKey(cs),
om.center=a.unique;
om.inversePoleFigureDirection=maxpos;
om.psi=S2DeLaValleePoussin('halfwidth',15*degree);
col = om.orientation2color(ebsdQ.orientations);
plot(ebsdQ,col)

% option 2
dist_a = angle(inv(ebsdQ.orientations)*maxpos,Miller(1,1,-2,0,cs));
cond = dist_a <=25*degree;
plot(ebsdQ(cond),dist_a(cond)/degree)

figure
plotPDF(ebsdQ.orientations,dist_a,hs)
annotate(maxpos,'add2all')



% %########################################################################
% c) Compare the <1,1,-2,0> distributions of those orientation which contain 
% c-axes at the periphery of the pole figure to those with c-axes in the center
% Maybe we need to wait till Friday to answer this question.
% %########################################################################

cs =   ebsdQ.CS;
cdir = ebsdQ.orientations.*Miller(0,0,0,1,cs);
cond1 = angle(cdir,zvector,'antipodal')<45*degree;
a=symmetrise(Miller(1,1,-2,0,cs));
aF1  =calcDensity(ebsdQ(cond1).orientations*a,'halfwidth',10*degree)
aF2  =calcDensity(ebsdQ(~cond1).orientations*a,'halfwidth',10*degree)

plot(aF1,'minmax','TR',{'pfJ:', num2str(round(10*mean(aF1.^2))/10)})
nextAxis
plot(aF2,'minmax','TR',{'pfJ:', num2str(round(10*mean(aF2.^2))/10)})
%%
 close all
 

% 