Questions about seismic anisotropy plots

[sauer.k...@gmail.com]

Hi all,

I am working on plotting seismic anisotropy from qtz-feld-mica mylonite samples. With the update to MTEX 5.0.3, I found the Multi_seismic_plot.m no longer working, but the SeismicsMultiplot.m works well.

However I am getting results that don't make sense to me. For example, the Vp plots show the fast direction along the primitive, where I would expect this (especially in strongly foliated and mica-rich mylonites) to be along the XY-plane with the slow direction to the north.
This may be too similar to David Mainprice's question from may 3!

This is an example of the script I'm using (just the part for Vp), after defining the stiffness tensors, ODFs, and Voight/Reuss/Hill averages:

%% Aggregate elastic properties

C_Aggregate_Hill
 = C_Hill_qtz*Quartz_Area_fraction +
C_Hill_olig*Oligoclase_Area_fraction + C_Hill_musc*Musc_Area_fraction +
C_Hill_bio*Bio_Area_fraction + C_Hill_calcite*Calcite_Area_fraction +
C_Hill_clino*Clino_Area_fraction + C_Hill_pyr*Pyrope_Area_fraction


%%
[vp,vs1,vs2,pp,ps1,ps2] = C_Aggregate_Hill.velocity('harmonic');

%%

% extrema
[maxVp, maxVpPos] = max(vp)
[minVp, minVpPos] = min(vp)

% percentage anisotropy
AVp = 200*(maxVp-minVp) / (maxVp+minVp)


%% Plotting section


% plotting convention - plot a-axis to east
plota2east;


Label_X='Z'; Label_Y='X';Label_Z='Y';


% set colour map to seismic color map : blue2redColorMap
setMTEXpref('defaultColorMap',blue2redColorMap)

% some options
blackMarker = {'Marker','s','MarkerSize',10,'antipodal',...
  'MarkerEdgeColor','white','MarkerFaceColor','black','doNotDraw'};
whiteMarker = {'Marker','o','MarkerSize',10,'antipodal',...
  'MarkerEdgeColor','black','MarkerFaceColor','white','doNotDraw'};

% some global options for the titles
%titleOpt = {'FontSize',getMTEXpref('FontSize'),'visible','on'}; %{'FontSize',15};
titleOpt = {'visible','on','color','k'};

% Setup multiplot
% define plot size [origin X,Y,Width,Height]
mtexFig = mtexFigure('position',[0 0 1000 1000]);

% set up spacing between subplots default is 10 pixel
%mtexFig.innerPlotSpacing = 20;


%% Vp : Plot P-wave velocity (km/s)


% Plot P-wave velocity (km/s)
plot(vp,'contourf','complete','upper')

mtexTitle('Vp (km/s)',titleOpt{:})

% extrema
[maxVp, maxVpPos] = max(vp);
[minVp, minVpPos] = min(vp);

% percentage anisotropy
AVp = 200*(maxVp-minVp) / (maxVp+minVp);

% mark axes labels
hold on
text([xvector,yvector,zvector],{Label_X,Label_Y,Label_Z},...
  'backgroundcolor','w');

% mark maximum with black square and minimum with white circle
hold on
plot(maxVpPos,blackMarker{:}) % cannot get it to plot only one black square
plot(minVpPos,whiteMarker{:})
hold off
colorbar
% subTitle
xlabel(['Vp Anisotropy = ',num2str(AVp,'%6.1f')],titleOpt{:})

Additionally, I was using the bit below, but the .symmetrise no longer works

plot(maxVpPos.symmetrise,blackMarker{:})
plot(minVpPos.symmetrise,whiteMarker{:})

Not enough input arguments.

Error in vector3d/symmetrise (line 59)
  v = reshape(S * v,length(S),length(v));

Error in vector3d/subsref (line 23)
    [varargout{1:nargout}] = builtin('subsref',v,s);

So, I get this:
















But I expect something closer to this (from Demsey et al., 2011, on very similar samples to what we're working with):


















Thanks for any help!

cheers,
Kat

[David Mainprice]

Dear All,

I am attaching to this e-mail the latest version of Multi_seismic_plot_MTEX5, specifically optimised for MTEX 5.0.3. Many of the new feature have been implemented by Ralf.

This is includes several options, see below for the details

all the best David

% Tested with MTEX 5.0.3 and MATLAB 2016a
%
% David Mainprice 10/05/2018
%
% Plot matrix of anisotropic seismic properties
%        1 Vp        2 AVs      3 S1 polarizations  4 Vs1
%        5 Vs2       6 dVs      7 Vp/Vs1            8 Vp/Vs2
%
% INPUT
%
% MTEX tensor defined with density
% C = stiffnessTensor(M,cs_tensor,'density',rho)
% User defined labels at MTEX X,Y, and Z positions reference frame
% on Vp plot (top right plot)
% You can define the labels to suit your application single crystal or
% polycrystalline stiffness tensor with triclinic symmetry('1')
%
% X_Label = for example 'a', 'X1' or 'X'
% Y_Label = for example 'b', 'X2' or 'Y'
% Z_Label = for example 'c', 'X3' or 'Z'
%
% hemisphere = 'upper' or 'lower'
%
% New feature is generating symmetically equivalent max and min Vp,Vs1,vs2
% etc using the harmonic method for symmetry of single crystals
%
% symm = 1 only one max and min values on the plot
% symm = 2 for generating symmetrise max and min values on the plot
%
% plot font_size = set your fontsize  25,20,15  etc
%
%**************************************************************************
% Remember to set the position of x-axis in specimen coordinates
% before calling this function
% plot2xeast,  plot2xnorth, plot2xwest or plot2xsouth
% to suit your single crystal or polycrystalline elastic stiffness tensor
%**************************************************************************
%
%
% Compute seismic velocities as functions
% using option velocity(C,[]) or C.velocity 
% BOTH generate velocities on sphere using harmonic method
% OR using the traditional velocity(C,x) where x is the propagation
% directions , many directions or a grid.
% such as XY_grid = equispacedS2Grid('upper','resolution',1*degree)
% [vp,vs1,vs2,pp,ps1,ps2] = velocity(C,XY_grid)
% Using a user defined grid use classical method not harmonic method.
 
% Seismic velocities as functions on the sphere with harmonic method
% is recommended for a smooth representation that actuately respects
% the symmetry of single crystals and can be plotted in 3D.
%
[vp,vs1,vs2,pp,ps1,ps2] = C.velocity('harmonic');
%
% you can sample the velocity values in any x,y,z direction by using
% vp_xyz = eval(vp,vector3d(x,y,z)) or
% vp_rho_theta = eval(vp,vector3d('rho',45*degree,'theta',90*degree))
% etc
%

[sauer.k...@gmail.com]

Hi David,

This is great, thanks! Now the .symmetrise works perfectly.

However, I'm still concerned about the results that I get. For strongly foliated quartzofeldspathic mylonites, I expect the fast direction to be in the foliation plane, not around the primitive. (see attached)

Any ideas on why I'm getting the plots that I am?

cheers,
Kat

[Luiz F. G. Morales]

Hi Kat

Is it possible that somewhere, somehow, in your in your script you have rotated the orientation data from, lets say, the XZ plane (your smaller figure with numbers) to the XY (your plot)? Or that you have sectioned this sample parallel to both foliation and lineation and measured the orientations in the XY place? Because to have a figure like that you would need to have the mica (001) poles well orientated parallel to the centre of your Vp plot

does this make sense?

best

Luiz

-- 
If you want to reduce the number of emails you get through this forum login to https://groups.google.com/forum/?fromgroups=#!forum/mtexmail, click "My membership" and select "Don't send me email updates". You can still get emails on selected topics by staring them.
--- 
You received this message because you are subscribed to the Google Groups "MTEX" group.
To unsubscribe from this group and stop receiving emails from it, send an email to mtexmail+u...@googlegroups.com.
Visit this group at https://groups.google.com/group/mtexmail.
For more options, visit https://groups.google.com/d/optout.
<actualResult.png><expectedResult.png>

ETH Zürich
Luiz F. G. Morales
Scientific Center for Optical and Electron Microscopy (ScopeM)
HPT D 9
Auguste-Piccard-Hof 1
8093 Zurich, Switzerland
Phone +41 44 633 37 46
E-mail: luiz.m...@scopem.ethz.ch
URL: www.scopem.ethz.ch
URL: https://sites.google.com/site/luizfgmorales/
Skype: lfgmorales

[Katrina Sauer]

Hi Luiz,

 

That was the only thing I could come up with as well, but I don’t see anywhere where I could have rotated anything in the code. Our thin section is cut in the XZ plane:

In the original import, I use :
setMTEXpref('xAxisDirection','east');

setMTEXpref('zAxisDirection','intoPlane');

 

and when setting the preferences for the anisotropy plots I use:

plota2east;

I agree that it appears it has been rotated – just not really sure how or where that would occur??

 

 

cheers,

Kat

You received this message because you are subscribed to a topic in the Google Groups "MTEX" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/mtexmail/S-PduPEWz8E/unsubscribe.
To unsubscribe from this group and all its topics, send an email to mtexmail+u...@googlegroups.com.

[Rüdiger Kilian]

Hi Kat,
how do some pole figures for the mica and for the blue phase (too small to read, most of the map) looks like?

Cheers,
Rüdiger

[Katrina Sauer]

Hey Rüdiger,

 

Quartz pole figs (the blue phase):

 

Muscovite – I filled in the non-indexed points with muscovite, to all be aligned, hence the very strong fabric:

 

Whereas Biotite is:

 

 

So looking at this it appears I’ve likely made a mistake in setting the muscovite orientation when using fill:

 

% n = notindexed

ebsd('n').phaseId=5 

 

% set some rotation e.g. a random one, but it could also be any

ebsd('n').rotations=rotation.rand(length(ebsd('n').rotations))

 

I noticed that setting the phaseId first gave all the muscovites the same orientation, whereas if I wanted random then setting the rotation first worked. I want them to all be the same orientation – then we get a maximum anisotropy. What would the correct way to do this be?

 

Cheers,

Kat

 

From: <mtex...@googlegroups.com> on behalf of Rüdiger Kilian <ruedige...@unibas.ch>
Reply-To:

<mtex...@googlegroups.com>

--

[ruediger Kilian]

Hi Kat,

so to understand correctly: you replace all unindexed pixels with muscovite (that would explain your polefigure which looks like all the orientations would be a bout Euler angles [0,0,0] )?

And you actually want them to be oriented in a way that (001) is parallel X?

You have a muscovite symmetry in your CSList? and it is phaseId 5, respectively phase 4?


Another strange thing is that I do not entirely understand the X,Y,Z annotation in your muscovite pole figures. Do you set the labels manually and this is not the plotting reference frame?

Cheers,
Rüdiger

[sauer.k...@gmail.com]

Hey Rüdiger,

I must have linked an older version of the muscovite PFs where I was messing around with the labeling to sort myself out. Just re-ran it to be sure, but the labels should be as you expect (like the quartz and biotite PFs).

Yes, I have muscovite in the list, CS{5} or the the 5th phase:
crystalSymmetry('12/m1', [5.189 9.004 20.256], [90,95.7,90]*degree, 'X||a*', 'Y||b', 'Z||c', 'mineral', 'Muscovite', 'color', 'magenta'),...

And yes, I would like (001) to be parallel to x.

Sorry for the confusion!

cheers,
Kat

[ruediger Kilian]

Hi Kat,
here's an example, replacing in the forsterite set all nonIndexed to diopside and inventing some common orientation on them.

mtexdata forsterite
% set notIndexed to e.g. diopside
ebsd('n').phase=3
% make up an orientation for diopside
ebsd('d').rotations=rotation('Euler',[0 90 10]*degree)

However, I'm not sure if this is really useful since looking at the amount of unindexed pixels in your maps and having them set to all the same orientation, you could also use an unimodal ODF of muscovite to begin with.


Cheers,
Rüdiger

[Katrina Sauer]

Hi Rüdiger,

This is definitely where the problem was, and it now gives exactly what I expect – Thank you SO SO much for your help!

Cheers,
Kat