atmos1D
index
/users/schrei_f/src/py4CAtS/art/atmos1D.py

atmos1D
 
Read atmospheric data file(s) in xy tabular ascii format:
extract profiles (columns), convert to cgs units, optionally truncate top or bottom levels, save in xy format.
If several files are given, the profiles from the second, third, ... file are interpolated to the grid of the first file.
 
Usage:
 
  atmos1D [options] file(s)
 
  -h              help
  -c    char      comment character(s) used in input,output file (default '#')
  -o    string    output file for saving of atmospheric data
  -i    string    interpolation method (default: numpy's linear interp)
  -r              replace duplicate profile from preceeding file(s) with new data
                  (default:  ignore duplicate profile(s) from second, third, ... file)
  -w    string    the leading word identifying the header row specifying the columns  (default: "what")
  -u    string    the leading word identifying the header row specifying the units    (default: "units")
  -x    string    eXtract profiles from file(s) (a string with comma seperated entries (default: "z,p,T,H2O"))
 --BoA   float    bottom-of-atmosphere altitude in km
 --ToA   float    top-of-atmosphere altitude in km
 --scale floats   multiply molecular concentrations with scaleFactors
                  (either a comma separated list of floats (no blanks) in the same order as for the line data files,
                  or just a single float to scale the profile of the molecule corresponding to the first lineFile.)
  -z     floats   regrid to (evtl. equidistant) altitude grid

 
Modules
       
numpy
re
sys

 
Functions
       
atmInfo(data)
Print some basic information about the atmospheric data.
atmMerge(atmData1, atmData2, replace=False, interpolate='linear', verbose=False)
Combine two atmospheric datasets (given as numpy structured array), interpolate second one if necessary.
 
ARGUMENTS:
----------
atmData1, atmData2:   the two atmospheres to merge
replace:              when a particular profile is given in both datasets,
                      use the first and ignore the second (default)
                      OR replace the first by the second
interpolate:          the interpolation method used for regridding (see the atmRegrid doc)
                      NOTE:  the final grid is given by the grid of atmData1
verbose:
 
RETURNS:              a unified atmospheric dataset (as a numpy structured array)
atmPlot(data, what='T', vertical='km', vmr='', **kwArgs)
Plot atmospheric profile(s) vs. altitude (or pressure).
 
ARGUMENTS:
----------
data:      either a numpy structured array of atmospheric data (z, p, T, gases)
           OR a list thereof
what:      the data field to plot, default T (=temperature)
vertical:  unit for vertical axis, default 'km' for altitude
           any of the altitude or pressure units defined in cgsUnits supported
vmr:       plot volume mixing ratio instead of number density (default number density)
           (ignored for non-gas profiles and when vmr is not a valid mixingRatioUnit)
kwArgs     passed directly to plot / semilogy and can be used to set colors, line styles and markers etc.
           ignored (cannot be used) in recursive calls with lists or dictionaries of atmos data.
atmRead(atmFile, commentChar='#', extract=None, zToA=0.0, zBoA=0.0, what='what', units='units', scaleFactors=None, returnUnits=False, verbose=False)
Read atmospheric data file, convert to cgs units, optionally truncate top/bottom or scale concentrations;
and return a structured numpy array.
 
Parameters:
-----------
atmFile        string   file(name) with data to be read
extract        strings  list of entries (essentially gases) to be returned (default None ==> read all)
zToA           float    top    of atmosphere altitude [km] (default: 0.0)
zBoA           float    bottom of atmosphere altitude [km] (default: 0.0)
what           string   first 'word' in file header line identifying column names (default: "what")
units          string   first 'word' in file header line identifying units        (default: "units")
scaleFactors   floats   change abundances, i.e. multiply molecular concentrations (default: 1.0 for all gases)
returnUnits    flag     optionally also return units (default: False)
 
Note: all data are returned in cgs units
      z[cm],  p[dyn/cm**2],  T[K],  n[molec/cm**3)
 
If you want to combine data from different files (e.g. 'standard' p&T and 'exotic' trace gas concentrations)
then read the files separately using this `atmRead` (or vmrRead) function and then use the `atmMerge` function
 
atmRead expects that pressure, temperature, and air density (at least two of these three!) are given;
If you want to read a file with only molecular VMRs or densities, use `vmrRead`.
 
Atmospheric data file(s):  see the data/atmos directory for examples.
atmRegrid(data, zNew, interpolate='linear')
Interpolate atmospheric data (given as structured array) to a new altitude grid.
 
ARGUMENTS:
----------
data          original atmospheric data (structured array)
zNew          the new altitude grid (list or numpy array)
interpolate   select a specific method,  default "" for linear Lagrange with numpy.interp
              0 | z   zero order spline
              1 | s   first order spline
              2 | q   second order spline
              3 | c   third order spline
              p | n   previous or next, i.e. use the left/right point
              h       Piecewise cubic Hermite (scipy.pchip)
              (see scipy.interpolate.interp1d)
 
HINT:
you can use the parseGridSpec function of the grid module to generate a new altitude grid,
e.g.:  zNew=parseGridSpec('0[2]16[4]60');  atmNew=atmRegrid(atmOld,zNew)
 
NOTE:
All atmospheric data incl. the altitude grid are in cgs units, hence z[cm].
However, for convenience zNew can be given in km:
if all new altitude values are 'small' (z<500), they are silently assumed to be kilometer and scaled to cm.
atmSave(atmos, outFile=None, units=None)
Write atmospheric data to ascii tabular file.
atmSave_namelist(atmos, outFile=None, vmr=False, zpT=False)
Write atmospheric data to namelist formatted file (used by GARLIC).
atmTruncate(atmos, zToA=None, zBoA=None, verbose=False)
Truncate atmosphere at top or bottom, i.e. delete levels for altitudes above zToA or below zBoA.
cmr(atmos, what=None, zMin=None, zMax=None, verbose=False)
Column mixing ratio:  quotient of molecular vertical column densities and air vertical column density.
densities(atmos, what=None)
Densities (1/cm**3) of molecules in atmospheric data, return a numpy array with shape (nLevels,nGases).
If a gas (or "air") is specified, return only this particular density.
dry_vmr(atmos, what=None)
Volume mixing ratios of molecules in atmospheric data with respect to dry air,
return a numpy array with shape (nLevels,nGases).
If a gas is specified, return only this particular VMR.
extract_profiles(data, units, extract='main')
Extract some profiles from the atmospheric data (given as structured array),
i.e., remove all other profiles.
 
Default is to return (extract) altitude, pressure, temperature and the 'main' gases H2O, CO2, O3, ...
gases(atmos, what=None)
Return names of molecules in atmospheric data.
If a gas is specified, only check if it is present.
pT(atmos)
Pressure and temperature:
given a numpy structured array of atmospheric data (z, p, T, gases)
return a numpy array with shape (nLevels,2).
scaleDensities(data, scaleFactors, gas=None, verbose=False)
Scale molecular densities (either a single gas or all).
If a single factor is given and no gas is specified, scale the very first.
scaleHeight(atmos)
Return the scale height [cm], i.e. the reciprocal of the slope of log(pressure) vs. altitude.
vcd(atmos, what=None, zMin=None, zMax=None)
Vertical column densities [molec/cm**2] for air and individual molecules.
If a gas (or "air") is specified, return only this particular vcd.
 
Returns integral_zMin^zMax density dz
 
ARGUMENTS:
what   default None, i.e. return all VCDs
zMin   lower integral limit (default BoA)
zMax   upper integral limit (default ToA)
 
Notes:
* small zMin, zMax values (z<250) are interpreted as km and automatically converted to cm;
* the integral(s) start/stop at the given altitude levels, i.e. zMin or zMax is "rounded" to the next grid point.
vmr(atmos, what=None)
Volume mixing ratios of molecules in atmospheric data, return a numpy array with shape (nLevels,nGases).
If a gas is specified, return only this particular VMR.
vmrRead(vmrFile, commentChar='#', extract=None, zToA=0.0, zBoA=0.0, what='what', units='units', scaleFactors=None, verbose=False)
Read volume mixing ratio (VMR) data file, convert to cgs units, optionally truncate top/bottom or scale concentrations;
and return a structured numpy array.
 
See the atmRead function for a full explanation of all arguments.
 
The main purpose of this function is to read trace gases concentrations from a file without pressure/temperature data;
then you can combine these VMRs with p, T, and the main gas number densities using the `atmMerge` function.
zpT(atmos)
Altitude, pressure and temperature:
given a numpy structured array of atmospheric data (z, p, T, gases)
return a numpy array with shape (nLevels,3).

 
Data
        ascii_uppercase = 'ABCDEFGHIJKLMNOPQRSTUVWXYZ'
goodNames = ['height', 'heights', 'altitude', 'altitudes', 'Height', 'Heights', 'Altitude', 'Altitudes', 'press', 'pressure', 'pressures', 'Press', 'Pressure', 'Pressures', 'temp', 'temperature', 'temperatures', 'Temp', 'Temperature', 'Temperatures', ...]
kBoltzmann = 1.3806504e-16
lengthUnits = {'A': 1e-08, 'au': 14959787070000, 'cm': 1.0, 'dm': 10.0, 'inch': 2.54, 'km': 100000.0, 'ly': 946073047258080000, 'm': 100.0, 'micrometer': 0.0001, 'mm': 0.1, ...}
mixingRatioUnits = {'%': 0.01, 'pp1': 1.0, 'ppV': 1.0, 'ppb': 1e-09, 'ppbV': 1e-09, 'ppbv': 1e-09, 'ppm': 1e-06, 'ppmV': 1e-06, 'ppmv': 1e-06, 'ppt': 1e-12, ...}
molecNames = ['H2O', 'HDO', 'CO2', 'O3', 'N2O', 'CO', 'CH4', 'CH3D', 'O2', 'NO', 'SO2', 'NO2', 'NH3', 'HNO3', 'OH', 'HF', 'HCl', 'HBr', 'HI', 'ClO', ...]
molecules = {'BrO': {'NumDeg': [1], 'TempExpGL': 0.5, 'TempExpQR': 1.0, 'VibFreq': [500.0], 'isotopes': ['69', '61'], 'mass': 95.0}, 'C2H2': {'NumDeg': [1, 1, 1, 2, 2], 'TempExpGL': 0.75, 'TempExpQR': 1.0, 'VibFreq': [3374.0, 1974.0, 3289.0, 629.0, 730.0], 'geisa': 24, 'hitran': 26, 'isotopes': ['1221', '1231', '1222'], 'mass': 26.03, 'sao': 26}, 'C2H4': {'NumDeg': [1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1], 'TempExpGL': 0.5, 'TempExpQR': 1.5, 'VibFreq': [3026, 1625, 1342, 1023, 3103, 1236, 949, 943, 3106, 826, 2989, 1444], 'geisa': 25, 'hitran': 38, 'isotopes': ['221', '311'], 'mass': 28.0}, 'C2H6': {'NumDeg': [1, 1, 1, 1, 1, 1, 2, 1, 2], 'TempExpGL': 0.75, 'TempExpQR': 1.9, 'VibFreq': [2899.0, 1375.0, 993.0, 275.0, 2954.0, 1379.0, 2994.0, 1486.0, 822.0], 'geisa': 22, 'hitran': 27, 'isotopes': ['1221', '1231'], 'mass': 30.07, 'sao': 27}, 'C2HD': {'NumDeg': [1, 1, 1, 2, 2], 'VibFreq': [3374.0, 1974.0, 3289.0, 629.0, 730.0], 'geisa': 48, 'isotopes': ['122'], 'mass': 17.0}, 'C2N2': {'NumDeg': [1, 1, 1, 1, 1], 'VibFreq': [2330, 846, 2158, 503, 234], 'geisa': 29, 'hitran': 48, 'isotopes': ['4224'], 'mass': 52.0}, 'C3H4': {'NumDeg': [1, 1, 1, 1, 1, 2, 2, 2, 2, 2], 'VibFreq': [3334, 2918, 2142, 1382, 931, 3008, 1452, 1053, 633, 328], 'geisa': 40, 'mass': 40.0}, 'C3H8': {'NumDeg': [1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, ...], 'VibFreq': [2977, 2962, 2887, 1476, 1462, 1392, 1158, 869, 369, 2967, 1451, 1278, 940, 216, 2968, 2887, 1464, 1378, 1338, 1054, ...], 'geisa': 28, 'isotopes': ['221'], 'mass': 44.0}, 'C4H2': {'NumDeg': [1, 1, 1, 1, 1, 1, 1, 1, 1], 'VibFreq': [3293, 2184, 874, 3329, 2020, 627, 482, 630, 231], 'geisa': 30, 'hitran': 43, 'isotopes': ['221'], 'mass': 50.0}, 'C6H6': {'NumDeg': [1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1], 'VibFreq': [3062, 992, 1326, 673, 3068, 1010, 995, 703, 1310, 1150, 849, 3063, 1486, 1038, 3047, 1596, 1178, 606, 975, 410], 'geisa': 47, 'isotopes': ['266'], 'mass': 78.0}, ...}
pressureUnits = {'N/m**2': 10.0, 'N/m^2': 10.0, 'Pa': 10.0, 'Torr': 1333.22, 'atm': 1013250.0, 'bar': 1000000.0, 'dyn/cm**2': 1.0, 'dyn/cm^2': 1.0, 'g/(cm*s**2)': 1.0, 'g/(cm*s^2)': 1.0, ...}
simpleNames = {'Altitude': 'z', 'Altitudes': 'z', 'Density': 'air', 'HGT': 'z', 'Height': 'z', 'Heights': 'z', 'PRE': 'p', 'Press': 'p', 'Pressure': 'p', 'Pressures': 'p', ...}