#!/usr/bin/python

# Convert a spectrum data file with vacuum wavelengths to air wavelengths
# File format is xy pairs with x in Angstroms
# Uses Edlen's formula as given by  Morton (Ap. J. Suppl. 77, 119)


import os
import sys
import numpy as np

if len(sys.argv) == 1:
  print " "
  print "Usage: vactoair.py infile.dat outfile.dat  "
  print " "
  sys.exit("Convert spectrum to air wavelenghts from vacuum wavelengths \n")
elif len(sys.argv) == 3:
  infile = sys.argv[1]
  outfile = sys.argv[2]
else:
  print " "
  print "Usage: vactoair.py infile.dat outfile.dat  "
  print " "
  sys.exit("Convert spectrum to air wavelengths from vacuum wavelengths \n")


# Read the spectrum

lambda_vac, signal = np.loadtxt(infile, dtype = 'float', unpack=True)

# Compute conversion factors for each wavelength
# Neglect second order corrections by using air wavelengths to find sigma

sigma2 = (1.e4/lambda_vac)**2.  

factor = 1. + 6.4328e-5 + 2.94981e-2/(146.-sigma2) + 2.5540e-4/(41.-sigma2)

factor = factor*(lambda_vac >= 2000.) + 1.*(lambda_vac < 2000.)

lambda_air = lambda_vac/factor  

# Save spectrum

dataout = np.column_stack((lambda_air, signal))
np.savetxt(outfile, dataout)
exit()