ReadMCA: Convert spectra gathered by MCA8000 into ROOT format (http://root.cern.ch/)

mca2root.py

#!/usr/bin/env python
# -*- coding: utf-8 -*-
'''Convert all alpha monitor spectrum into ROOT format.'''

from ROOT    import TFile, TH1F
from ReadMCA import readmca
import simplejson as sj
import os, pickle

def mca2root_file(fn,hn="",ht=""):
    '''Conver MCA file <fn> to TH1F object.'''
    print "Converting %s ..."%fn
    spc = readmca(fn)
    if not spc or not (spc.has_key('INFO') and spc.has_key('DATA')):
        return None
    binsX = len(spc['DATA'])
    if not hn: hn = "hmca"
    if not ht: ht = "Histogram converted from MCA file: %s"%fn
    h = TH1F(hn, ht, binsX, 0, float(binsX))
    for i in range(binsX): h.Fill(i, spc['DATA'][i])
    spc.pop('DATA')
    return [spc, h]

def mca2root_fdir(dn, frname="", isNewFile=False, isInDir=False, isOutJson=True):
    '''Convert all MCA spectra in directory <dn>.'''
    specialChars = [' ', '-', '.', '#', '(', ')', '+', '*']

    # Generate file list
    flist = []
    if os.path.isdir(dn):
        # Find valid files in dir <dn>
        for r,d,files in os.walk(dn):
            if not files:   continue
            for f in files:
                if f[-4:].lower() == ".mca":
                    flist.append(os.path.join(r,f))
    elif os.path.isfile(dn):
        flist.append(dn)
    else:
        print dn, "-- Not a valid file or directory!"
        return None

    if not flist:
        print "No valid MCA file found!"
        return None

    frbase = ""
    if frname == "":
        if dn == "." or dn == "./":
      frbase = os.path.abspath("./").split(os.path.sep)[-1]
        elif len(dn) > 4 and dn[-4:].lower() == ".mca":
            frbase = dn[:-4]
        else:
            frbase = dn
        frname = frbase+"_mca.root"
    elif len(frname)>5 and frname[-5:].lower() == ".root":
        frbase = frname[:-5]
    else:
        frbase = frname
        frname = frbase + ".root"

    # Open ROOT file
    opt = "RECREATE" if isNewFile else "UPDATE"
    froot = TFile(frname,opt)
    minfo = []

    # Convert files
    for fm in flist:
        # Create same directory hierarchy in ROOT file
        fdir = os.path.dirname(fm) if isInDir else \
            os.path.dirname(fm[len(os.path.join(dn,"")):])

        hdir = [] if not fdir else fdir.split(os.path.sep)
        rcdir = froot
        for d in hdir:
            rcdir = rcdir.mkdir(d) if not rcdir.Get(d) else rcdir.Get(d)
            rcdir.cd()

        # Convert from MCA to TH1F
        spc, h = mca2root_file(fm)
        fnbase = os.path.basename(fm)[:-4]
        for c in specialChars: fnbase = fnbase.replace(c,'_')
        h.SetName("h_"+fnbase)
        h.SetTitle("Spectrum from %s"%fm)
        h.Write()
        spc['file'] = fm
        minfo.append(spc)

    froot.Close()

    # output MCA info in JSON format
    if minfo:
        if isOutJson:
            print "Writing MCA information to file %s_mca.json ..."%frbase
            open(frbase+"_mca.json","w").write(sj.dumps(minfo,indent=4))
        else:
            print "Writing MCA information to file %s_mca.pkl ..."%frbase
            pickle.dump(minfo,open(frbase+"_mca.pkl","wb"))

    return minfo

import argparse as ap
def pmain(args):
    par = ap.ArgumentParser(
        prog='mca2root', prefix_chars = '-+',
        description = 'Convert MCA histograms into ROOT format')
    par.add_argument( '-o', '--output', dest = 'output', nargs='?',
                      help='output filename')
    par.add_argument( '-u', '--update', action='store_true',
                      help='update existing output (default: No)')
    par.add_argument( '-d', '--dirtree', action='store_true',
                      help='output histograms in dirs')
    par.add_argument( '-j', '--json', action='store_true',
                      help='output MCA info as JSON format (default: No)')
    par.add_argument( 'name', nargs='+', help="files or directories")

    vas = par.parse_args(args)
    frname  = "" if not vas.output else vas.output
    b_newf  = True
    b_indir = vas.dirtree
    b_json  = vas.json
    if not vas.update or ( len(vas.name)>1 and len(frname)>1 ):
        b_newf = False

    # loop name list
    for fn in vas.name:
        print "Processing %s ..."%fn
        mca2root_fdir(fn, frname, b_newf, b_indir, b_json)

if __name__=='__main__':
    import sys
    if len(sys.argv)<2: pmain(['-h',])
    else: pmain(sys.argv[1:])

ReadMCA.py

#!/usr/bin/env python
# -*- coding: utf-8 -*-
'''Read data acquired by MCA8000A through software "pcma".
File extension: .mca'''
import re

def readmca(fm):
    sraw=[l.strip() for l in open(fm,'r').readlines()]

    kidx=[]
    for i in range(len(sraw)):
        keys=re.findall('<<(.*)>>',sraw[i])
        if keys: kidx.append([i,keys[0]])
    knum = len(kidx)
    if knum<3:  print "WARNING: not a proper MCA spectrum!"

    spc = dict()
    for k in range(knum):
        istart = kidx[k][0]+1
        if k+1<knum: iend = kidx[k+1][0]
        else: iend = -1
        # PMCA SPECTRUM
        if kidx[k][1] == 'PMCA SPECTRUM':
            spc['INFO'] = dict()
            for l in sraw[istart:iend]:
                vs = [i.strip() for i in l.split('-')]
                if len(vs)>1:
                    vvs="".join(vs[1:])
                    if vvs.isdigit(): v = eval(vvs)
                    else: v = vvs
                    spc['INFO'][vs[0]] = v
                else:
                    print "WARNING: empty info key --", vs[0]
        # Calibration
        elif kidx[k][1] == 'CALIBRATION':
            spc['CALI'] = dict()
            spc['CALI']['info'] = [i.strip() for i in sraw[istart].split('-')]
            if iend != -1 and iend <= istart + 1: continue
            spc['CALI']['data'] = []
            for l in sraw[istart+1:iend]:
                spc['CALI']['data'].append([eval(i) for i in l.split()])
        # ROI
        elif kidx[k][1] == 'ROI':
            spc['ROI'] = []
            for l in sraw[istart:iend]:
                spc['ROI'].append([eval(i) for i in l.split()])
        # Data
        elif kidx[k][1] == 'DATA':
            spc['DATA'] = [int(i) for i in sraw[istart:iend]]
        else: pass

    return spc

if __name__=="__main__":
    import sys
    if len(sys.argv)<2:
        print "Usage: %s <file.mca>"%sys.argv[0]
        exit()
    spc = readmca(sys.argv[1])
    if spc.has_key('INFO'):
        for k in spc['INFO']: print k, '::', spc['INFO'][k]