Kikuchi bands from Silicon.

A rotating cube...

I have been learning Python. I want to work towards writing some scripts that together take EBSD-based orientation data which are relative to a sample surface reference frame and rotate the data for visualization relative to a kinematic or geographic reference frame.
I decided to start simple – learn something about quaternions and rotate 2D plots of a cube.

Below is some code to create an animation (.avi file)  by rotating cube axes.
Practicing. Learning. Far from the goal.

NOTE: IF you want to use this code, you'll need mencoder for the script to successfully create the animation file.

A glimpse of a video that this script produced: cube movie


__________________________________________________________________________

#!/usr/bin/python

##########################################################################

## Cube animation program using python, numpy, os, and matplotlib.

## This script demonstrates how rotate plots using quaternion conversion.

## 

## Written by Zach Michels. Last Edited 5/12/2012.

## Cheers!

##########################################################################

from pylab import *

import numpy as np

import subprocess                 # For issuing commands to the OS.

import os, time, stat, Image

import sys                        # For determining the Python version.

##########################################################################

## Produces a rotation matrix from the 3 last components of a quaternion.

def quaternion_to_matrix(myx):

    # '''converts from a quaternion representation (last 3 values) to rotation matrix.'''

    xb,xc,xd = myx

    xnormsq = xb*xb+xc*xc+xd*xd

    if xnormsq < 1:

        ## Inside the unit sphere, these are the components,

        ## and we only have to calculate 'a' to normalize

        ## the quaternion.

        b,c,d = xb,xc,xd

        a = np.sqrt(1-xnormsq)

    else:

        ## If we have invalid inputs, we reflect the vectors

        ## outside the unit sphere to the other side, and

        ## use the inverse norm and negative values of 'a'.

        b,c,d = -xb/xnormsq,-xc/xnormsq,-xd/xnormsq

        a = -np.sqrt(1-  1.0/xnormsq  )

        ## This generally defeats the whole concept that

        ## small (b,c,d) vector norms have small rotation angles. It's

        ## really just to work near the borders of the sphere. 

        ## Optimization algorithms should work with initalizations 

        ## just inside the spere, and avoid wandering outside of it.

    assert a >= -1

    assert a <= 1

    ## Notice we return a transpose matrix, because we work with line-vectors

    return np.array([ [(a*a+b*b-c*c-d*d), (2*b*c-2*a*d),     (2*b*d+2*a*c)      ],

                         [(2*b*c+2*a*d),     (a*a-b*b+c*c-d*d), (2*c*d-2*a*b)      ],

                         [(2*b*d-2*a*c),     (2*c*d+2*a*b),     (a*a-b*b-c*c+d*d)] ]  ).T  \

                         / (a*a+b*b+c*c+d*d)

                         

###########################################################################

## Calculate quaternion values using sinusoids over the frame index

## number, and get rotation matrix. Notice that the second parameter

## is the number of frames. That allows us to easily create an

## animation that can be played as a loop.

def crazy_rotation(ind,Nind):

    return quaternion_to_matrix(0.5*sin(pi*2*ind*Nind**-1.*array([1,2,3])))

## This function calculates the projected image coordinates form the 3D

## coordinates of the object vertices.

def project(D, vecs):

    return vvs[:,:2]/(vvs[:,[2,2]]-D)

## The cube vertices

vs = reshape(mgrid[-1:2:2,-1:2:2,-1:2:2].T, (8,3))

## Generate the list of connected vertices

ed=[(j,k)

    for j in range(8)

    for k in range(j,8)

    if sum(abs(vs[j]-vs[k]))==2 ]

## 'Camera' position

D=-5

## Create the figure. figsize needed to set the image sizes at the end...

figure(1, figsize=(6.4,4.8))

## create/get the axes object.

ax=subplot(1,1,1)

## Set image title.

suptitle("Zach's rotating cube")

## Number of frames. Lets user set number.

Nind = input("How many frames would you like to make in your animation? Each frame corresponds to a percent of the overall intended rotations, such that a greater number of frames ultimately produces a slower animation. This script will write a single .PNG image for each frame into the same directory as the .py file. These .PNG files are used to make the .AVI file, but are not automatically deleted. You will be prompted later for the directory PATH if you would like to delete them. NOTE: a greater number of frames will require longer time for script processing. NUMBER OF FRAMES (200+ recommended): ")

##

## Now start a loop over the animation frames. The index is used both

## to name the images and to calculate the rotation matrix. You could

## use this index to make some physics stuff...

for ind in range(Nind):

    print ind

    ## This is crucial, clears the figure for the new plot. 

    ax.clear() 

    ## Get the rotation matrix for the current frame.

    rotM = crazy_rotation(ind, Nind)

    ## Calculate the 3D coordinates of the vertices of the rotated

    ## cube. Just multiply the vectors by the rotation matrix...

    vvs=dot(vs,rotM)

    ## Now calculate the image coordinates of the points.

    pt = project(D,vvs)

    ## Plot the edges.

    for j,k in ed:

        ax.plot(pt[[j,k],0], pt[[j,k],1], 'm-', lw=3)

    ## Plot the vertices.

    ax.plot(pt[:,0], pt[:,1], 'bo')

    ## Set axes limits

    ax.axis('equal')

    ax.axis([-0.5,0.5,-0.5,0.5])

    ## Save the current frame. We need the dpi (along with the figure

    ## size up there) to set the image size.

    savefig('anim%03d.png'%ind, dpi=100)

# Now stitch graphed images together using Mencoder to create a movie.

# Each image will become a single frame in the movie.

#

# Use Python to make what would normally be a command line

# call to Mencoder.  Specifically, the command line call we want to

# emulate is (without the initial '#'):

# mencoder mf://*.png -mf type=png:w=800:h=600:fps=25 -ovc lavc -lavcopts 

# vcodec=mpeg4 -oac copy -o output.avi

# See the MPlayer and Mencoder documentation for details.

#

command = ('mencoder',

           'mf://*.png',

           '-mf',

           'type=png:w=800:h=600:fps=25',

           '-ovc',

           'lavc',

           '-lavcopts',

           'vcodec=mpeg4',

           '-oac',

           'copy',

           '-o',

           'cube_zm.avi')

# os.spawnvp(os.P_WAIT, 'mencoder', command)

print "\n\nabout to execute:\n%s\n\n" % ' '.join(command)

subprocess.check_call(command)

# explain compiled video output file and format

print "\n\n The movie was written to 'cube_zm.avi'"

print "\n\n You may want to delete *.png now.\n\n"

# opportunity to delete .PNG files used to create .AVI animation. 

# Note: If the directory PATH supplied is incorrect, the script 

# will return an error and the files will not be dleted. 

# However, the file cube.avi should still be created.

if raw_input("Delete all .PNG files in the directory? (y/n): ") == "y":

path = raw_input("Path to directory? (ex: '/Users/Usernname/folder'):  ")

if raw_input("are you sure you want to delete these files? Ths cannot be undone. (y/n): ") == "y":

dircount=0

pngcount=0

# deletes .PNG files by stepping through working directory 

# (the directory that contains the .py file).

for file in os.listdir(path):

    dircount=dircount+1

    basename=os.path.basename(file)

    NoExtName=os.path.splitext(file)[0]

   

    if basename.endswith('.png'):

pngcount=pngcount+1

print 'Deleting '+basename

os.remove(file)

Go Micro.

Serendipity and science...

As I get to know more folks who are well established in their research niches, I am struck by the trend that many of them recount influential serendipitous events in their academic/professional careers which guided them down a new path, perhaps to ultimately do lots of work in a field they never intended or imagined. Indeed this should not surprise me –– it seems commonplace that strong minds in sciences embrace interconnectedness and feedback between disciplines and tend to follow their interests such that they get their fingers into all kinds of stuff. And, I have no doubt that whatever drives such research diversification and/or the desire to descend new rabbit holes is different for everyone. However, I think that this trend likely reflects the heart of the scientific method. That is, science is an evolving dialogue that (ideally) does not have a specific trajectory and rewards scientists with adaptable attitudes.

"So WHAT?" you say... "You are rambling about the obvious."

Well, I'm fixated on this because I feel that many grad students –– myself included –– set out in their academic careers with their eye on some prize (everyone's proverbial carrot is a bit different)... and this may not be the most useful approach. That said, it may get you in and out of school faster and on to bigger and better things for yourself, and in that case it is wonderful for some folks. But for the scientists that inspire me, it is clear it has behooved them to be sensitive and malleable to the research climate and needs of the community. Mostly, such an attitude rewards them with more work and more resources (intellectual and financial). 

I arrived at Madison with a Masters from SJSU and no interest in EBSD, although I have kept a soft spot for microstructures and microspace for a quite while. I imagined a relatively straightforward path for myself, felt like I had a clear plan. I consider/ed myself primarily a field geologist... which may be just a cooler way of saying that I aspire to be a generalist, or maybe that I don't know enough to further specialize my skill set. 

I am a late-comer to geoscience (more about that another time), and I sometimes remark at the disconnect between my passion for geology and my overall lack of preparation for it. I often feel like I took a stab in the dark and struck gold. Yet, it is clear to me that among my fellow grad students (most of whom are rather brilliant at something or other), I am not even a strong "generalist" –– I have a lot to learn. 

And I am learning. I'm learning that I love EBSD –– that it augments my love for microstructure. And that may change. But it is exciting to be jazzed about something new –– it is suddenly all I want to do. That likely *will* change. But it is already opening my vision and opportunities, and I know that it is going to become a big part of what I do over the next few years.

Of course this is what school is for. And you probably find my blathering to be more predictable than remarkable. But now that I have written it down, I think what I am trying to say is that it it is powerful to remind ourselves that our roles in science are not predetermined. And that the community often benefits from the work of those who follow their intrigue and can flow.

MAC

It seems that the assignment of mass absorption coefficients (MAC) for elements involves some serious hocus-pocus and much assumption.

All of my officemates are out this week – conference in Italy; drilling on a ship in Japan.

Maybe I'll get more done this week, but I doubt it.

Filament burned out..

An update on my experience as an EBSD noob:

So... yesterday one of the postdocs was kind enough to let me sit in on some EBSP (Electron Backscatter Pattern) analysis he was conducting on some dunite from the Twin Sisters. He set it up to map orientation contrast (OC) of the olivine lattice planes across a 2x3 matrix with a 5um step-size. He started the maping at about 2:30pm and checked in on it at about 11pm only to find that the tungsten filament in the SEM burned out. A couple of the map "blocks" (in the matrix) completed prior to the burnout and are useful, but he will likely be back in to finish the mapping next week. I am going to join him for that, too.

Also, my academic advisor just mentioned trying to put together a short-course this sumer that Seth will lead... a crash-course type of thing on "all you ever needed to know about EBSD."

Every Tuesday...

... is SEM day!

Location:W Dayton St,Madison,United States

Postgres regardless of progress...

Today I continued some work on a django interface for my Postgres database. It's actually quite easy, but I'm a complete noob... so only recently have I got it doing what I want. Basically, I now have a way to enter observations into forms that save all the info into its proper place in a relational database schema. I think it will be worthwhile in the future when I accumulate larger volumes of data (more mapping, sampling, thin sections, and EBSD data), and I look forward to the GIS integration with QGIS and GRASS. Anywhooo... Where's that paper I should be reading...?

EBSD

Welcome. This first post is just an excitement.
I signed up to attend the MAS workshop on EBSD at Carnegie Mellon in June.
I am excited about it.
A hefty portion of my (near-future) research will involve analyzing crystallographic axes of mantle materials to infer dominant slip-systems during deformation.