#!/usr/bin/env python
from optparse import OptionParser
 
from Dacapo import Dacapo
from ASE.Utilities.DOS import DOS
 
cmd = OptionParser(usage = '%prog [-s] input_nc_file output_text_file')
cmd.add_option('-r', '--reverse', action = "store_true", default = False,
               help = 'reverse out put for minor spin states')
(opt, argv) = cmd.parse_args()
 
if len(argv) != 2:
    cmd.print_help()
    raise SystemExit
 
ncfile = argv[0]
textfile = argv[1]
sfmt0 = '%-15s %-15s\n'
sfmt1 = '%-15s %-15s %-15s %-15s\n'
efmt0 = '%15.7E %15.7E\n'
efmt1 = '%15.7E %15.7E %15.7E %15.7E\n'
 
atoms = Dacapo.ReadAtoms(ncfile)
calc = atoms.GetCalculator()
isSpinPolarized = calc.GetSpinPolarized()
dos = DOS(calc)
x0 = dos.GetEnergies()
if isSpinPolarized: 
  y0 = dos.GetDOS(0)
  y1 = dos.GetDOS(1)
  fp = open(textfile, 'w')
  if opt.reverse:
    fp.write(sfmt1 % ('Energy', 'DOS0', 'DOS1', 'DIFF'))
  else:
    fp.write(sfmt1 % ('Energy', 'DOS0', 'DOS1', 'TOTAL'))
  for i in range(len(x0)):
    if opt.reverse:
      fp.write(efmt1 % (x0[i], y0[i], -y1[i], y0[i] - y1[i]))        
    else:
      fp.write(efmt1 % (x0[i], y0[i],  y1[i], y0[i] + y1[i]))        
  fp.close()
else:
  y0 = dos.GetDOS(0)
  fp = open(textfile, 'w')
  fp.write(sfmt0 % ('Energy', 'DOS'))
  for i in range(len(x0)):
    fp.write(efmt0 % (x0[i], y0[i]))
  fp.close()

$ getdos
usage: getdos [-r] input_nc_file output_text_file
 
options:
  -h, --help     show this help message and exit
  -r, --reverse  reverse output for minor spin states
$