A simple script to calculate average pore diameter.

avg_pore_dia.py

""" A simple script to calculate average pore diameter.

Written for Ayush on 13th Oct 2013.

Usage: python avg_pore_dia.py </dir/containing/pore-files>

This script assumes that the pore files are .txt files.  Point the
script to a directory containing all your pore files, and it should
print out the average pore diameter from all the files.

"""

from pylab import argmin, average, loadtxt

def diff(arr):
    """ Return difference of consecutive elements.

    Note: new array is 1 smaller in length.

    """

    return arr[:-1] - arr[1:]


def find_diameter(filename):
    """ Returns pore diameter, given a line profile file. """

    # Read the file, and get x, y co-ordinates.
    data = loadtxt(filename)
    x, y = data.transpose()
    left_peak = right_peak = min_index = argmin(y)

    # Calculate difference of consecutive elements in y
    diff_array = diff(y)
    left_array = diff_array[:min_index]
    right_array = diff_array[min_index:]

    # Iterate over the left array and get index
    for element in left_array[::-1]:
        if element > 0.0:
            left_peak -= 1
        else:
            break

    for element in right_array:
        if element < 0.0:
            right_peak += 1
        else:
            break

    return x[right_peak] - x[left_peak]

if __name__ == "__main__":

    import os
    from os.path import join
    from pprint import pprint
    import sys

    target_dir = sys.argv[1] if len(sys.argv) == 2 else '.'

    pore_files = [
        join(target_dir, path) for path in os.listdir(target_dir)
        if path.endswith('.txt')
    ]

    if len(pore_files) == 0:
        print 'Usage: python %s /dir/containing/pore-files' % sys.argv[0]
        sys.exit(1)

    diameters = [ find_diameter(pore) for pore in pore_files ]

    pprint(sorted(zip(pore_files, diameters)))

    print 'Average diameter: ', average(diameters)