[yt-svn] commit/yt: 10 new changesets
commits-noreply at bitbucket.org
commits-noreply at bitbucket.org
Mon Jun 1 08:38:17 PDT 2015
10 new commits in yt:
https://bitbucket.org/yt_analysis/yt/commits/5e319f9fa9d6/
Changeset: 5e319f9fa9d6
Branch: yt
User: brittonsmith
Date: 2015-05-20 21:51:30+00:00
Summary: Removing hard-coded unit conversions from absorption spectrum generator.
Affected #: 2 files
diff -r 4cba8216c42868a394e99e3b10cedd31bec9e704 -r 5e319f9fa9d67997fd9c1ec0fe2ecbc293d071b4 yt/analysis_modules/absorption_spectrum/absorption_line.py
--- a/yt/analysis_modules/absorption_spectrum/absorption_line.py
+++ b/yt/analysis_modules/absorption_spectrum/absorption_line.py
@@ -171,26 +171,28 @@
## shift lam0 by deltav
if deltav is not None:
- lam1 = lam0 * (1 + deltav / c)
+ lam1 = lam0 * (1 + deltav / ccgs)
elif delta_lambda is not None:
lam1 = lam0 + delta_lambda
else:
lam1 = lam0
## conversions
- vdop = vkms * cm_per_km # in cm/s
- lam0cgs = lam0 / 1.e8 # rest wavelength in cm
- lam1cgs = lam1 / 1.e8 # line wavelength in cm
- nu1 = ccgs / lam1cgs # line freq in Hz
- nudop = vdop / ccgs * nu1 # doppler width in Hz
- lamdop = vdop / ccgs * lam1 # doppler width in Ang
+ vdop = vkms.in_units("cm/s") # in cm/s
+ lam0cgs = lam0.in_units("cm") # rest wavelength in cm
+ lam1cgs = lam1.in_units("cm") # line wavelength in cm
+ nu1 = (ccgs / lam1cgs).in_units("1/s") # line freq in Hz
+ nudop = (vdop / ccgs * nu1).in_units("1/s") # doppler width in Hz
+ lamdop = (vdop / ccgs * lam1).in_units("angstrom") # doppler width in Ang
+ #import pdb ; pdb.set_trace()
+
## create wavelength
if lambda_bins is None:
lambda_bins = lam1 + \
np.arange(n_lambda, dtype=np.float) * dlambda - \
n_lambda * dlambda / 2 # wavelength vector (angstroms)
- nua = ccgs / (lambda_bins / 1.e8) # frequency vector (Hz)
+ nua = (ccgs / lambda_bins).in_units("1/s") # frequency vector (Hz)
## tau_0
tau_X = np.sqrt(np.pi) * e**2 / (me * ccgs) * \
diff -r 4cba8216c42868a394e99e3b10cedd31bec9e704 -r 5e319f9fa9d67997fd9c1ec0fe2ecbc293d071b4 yt/analysis_modules/absorption_spectrum/absorption_spectrum.py
--- a/yt/analysis_modules/absorption_spectrum/absorption_spectrum.py
+++ b/yt/analysis_modules/absorption_spectrum/absorption_spectrum.py
@@ -4,7 +4,7 @@
"""
-from __future__ import print_function
+
from __future__ import absolute_import
#-----------------------------------------------------------------------------
@@ -21,13 +21,12 @@
from .absorption_line import tau_profile
from yt.funcs import get_pbar, mylog
-from yt.units.yt_array import YTArray
+from yt.units.yt_array import YTArray, YTQuantity
from yt.utilities.physical_constants import \
amu_cgs, boltzmann_constant_cgs, \
- speed_of_light_cgs, km_per_cm
+ speed_of_light_cgs
from yt.utilities.on_demand_imports import _astropy
-speed_of_light_kms = speed_of_light_cgs * km_per_cm
pyfits = _astropy.pyfits
class AbsorptionSpectrum(object):
@@ -49,8 +48,10 @@
self.tau_field = None
self.flux_field = None
self.spectrum_line_list = None
- self.lambda_bins = np.linspace(lambda_min, lambda_max, n_lambda)
- self.bin_width = (lambda_max - lambda_min) / float(n_lambda - 1)
+ self.lambda_bins = YTArray(np.linspace(lambda_min, lambda_max, n_lambda),
+ "angstrom")
+ self.bin_width = YTQuantity((lambda_max - lambda_min) /
+ float(n_lambda - 1), "angstrom")
self.line_list = []
self.continuum_list = []
@@ -76,7 +77,8 @@
mass of atom in amu.
"""
self.line_list.append({'label': label, 'field_name': field_name,
- 'wavelength': wavelength, 'f_value': f_value,
+ 'wavelength': YTQuantity(wavelength, "angstrom"),
+ 'f_value': f_value,
'gamma': gamma, 'atomic_mass': atomic_mass,
'label_threshold': label_threshold})
@@ -209,16 +211,15 @@
# include factor of (1 + z) because our velocity is in proper frame.
delta_lambda += line['wavelength'] * (1 + field_data['redshift']) * \
field_data['velocity_los'] / speed_of_light_cgs
- thermal_b = km_per_cm * np.sqrt((2 * boltzmann_constant_cgs *
- field_data['temperature']) /
- (amu_cgs * line['atomic_mass']))
- thermal_b.convert_to_cgs()
+ thermal_b = np.sqrt((2 * boltzmann_constant_cgs *
+ field_data['temperature']) /
+ (amu_cgs * line['atomic_mass']))
center_bins = np.digitize((delta_lambda + line['wavelength']),
self.lambda_bins)
# ratio of line width to bin width
width_ratio = ((line['wavelength'] + delta_lambda) * \
- thermal_b / speed_of_light_kms / self.bin_width).value
+ thermal_b / speed_of_light_cgs / self.bin_width).in_units("").d
if (width_ratio < 1.0).any():
mylog.warn(("%d out of %d line components are unresolved, " +
@@ -240,11 +241,13 @@
my_bin_ratio = spectrum_bin_ratio
while True:
lambda_bins, line_tau = \
- tau_profile(line['wavelength'], line['f_value'],
- line['gamma'], thermal_b[lixel],
- column_density[lixel],
- delta_lambda=delta_lambda[lixel],
- lambda_bins=self.lambda_bins[left_index[lixel]:right_index[lixel]])
+ tau_profile(
+ line['wavelength'], line['f_value'],
+ line['gamma'], thermal_b[lixel].in_units("km/s"),
+ column_density[lixel],
+ delta_lambda=delta_lambda[lixel],
+ lambda_bins=self.lambda_bins[left_index[lixel]:right_index[lixel]])
+
# Widen wavelength window until optical depth reaches a max value at the ends.
if (line_tau[0] < max_tau and line_tau[-1] < max_tau) or \
(left_index[lixel] <= 0 and right_index[lixel] >= self.n_lambda):
@@ -260,7 +263,7 @@
if line['label_threshold'] is not None and \
column_density[lixel] >= line['label_threshold']:
if use_peculiar_velocity:
- peculiar_velocity = km_per_cm * field_data['velocity_los'][lixel]
+ peculiar_velocity = field_data['velocity_los'][lixel].in_units("km/s")
else:
peculiar_velocity = 0.0
self.spectrum_line_list.append({'label': line['label'],
@@ -271,7 +274,7 @@
'redshift': field_data['redshift'][lixel],
'v_pec': peculiar_velocity})
pbar.update(i)
- pbar.finish()
+ pbar.finish()
del column_density, delta_lambda, thermal_b, \
center_bins, width_ratio, left_index, right_index
@@ -280,7 +283,7 @@
"""
Write out list of spectral lines.
"""
- print("Writing spectral line list: %s." % filename)
+ mylog.info("Writing spectral line list: %s." % filename)
self.spectrum_line_list.sort(key=lambda obj: obj['wavelength'])
f = open(filename, 'w')
f.write('#%-14s %-14s %-12s %-12s %-12s %-12s\n' %
@@ -295,7 +298,7 @@
"""
Write spectrum to an ascii file.
"""
- print("Writing spectrum to ascii file: %s." % filename)
+ mylog.info("Writing spectrum to ascii file: %s." % filename)
f = open(filename, 'w')
f.write("# wavelength[A] tau flux\n")
for i in range(self.lambda_bins.size):
@@ -307,7 +310,7 @@
"""
Write spectrum to a fits file.
"""
- print("Writing spectrum to fits file: %s." % filename)
+ mylog.info("Writing spectrum to fits file: %s." % filename)
col1 = pyfits.Column(name='wavelength', format='E', array=self.lambda_bins)
col2 = pyfits.Column(name='flux', format='E', array=self.flux_field)
cols = pyfits.ColDefs([col1, col2])
@@ -319,7 +322,7 @@
Write spectrum to an hdf5 file.
"""
- print("Writing spectrum to hdf5 file: %s." % filename)
+ mylog.info("Writing spectrum to hdf5 file: %s." % filename)
output = h5py.File(filename, 'w')
output.create_dataset('wavelength', data=self.lambda_bins)
output.create_dataset('tau', data=self.tau_field)
https://bitbucket.org/yt_analysis/yt/commits/287d5d45255b/
Changeset: 287d5d45255b
Branch: yt
User: brittonsmith
Date: 2015-05-21 17:48:17+00:00
Summary: Cleaning up voigt profile generator and making it use units more properly.
Affected #: 1 file
diff -r 5e319f9fa9d67997fd9c1ec0fe2ecbc293d071b4 -r 287d5d45255bc8fd3a2ea707e5cdb8dd89daeb2f yt/analysis_modules/absorption_spectrum/absorption_line.py
--- a/yt/analysis_modules/absorption_spectrum/absorption_line.py
+++ b/yt/analysis_modules/absorption_spectrum/absorption_line.py
@@ -16,12 +16,9 @@
import numpy as np
from yt.utilities.physical_constants import \
charge_proton_cgs, \
- cm_per_km, \
- km_per_cm, \
mass_electron_cgs, \
speed_of_light_cgs
-
def voigt(a,u):
"""
NAME:
@@ -140,69 +137,75 @@
k1 = k1.astype(np.float64).clip(0)
return k1
-def tau_profile(lam0, fval, gamma, vkms, column_density,
- deltav=None, delta_lambda=None,
+def tau_profile(lamba_0, f_value, gamma, v_doppler, column_density,
+ delta_v=None, delta_lambda=None,
lambda_bins=None, n_lambda=12000, dlambda=0.01):
- """
+ r"""
Create an optical depth vs. wavelength profile for an
absorption line using a voigt profile.
- :param lam0 (float): central wavelength (angstroms).
- :param fval (float): f-value.
- :param gamma (float): gamma value.
- :param vkms (float): doppler b-parameter.
- :param column_density (float): column density (cm^-2).
- :param deltav (float): velocity offset from lam0 (km/s).
- Default: None (no shift).
- :param delta_lambda (float): wavelength offset in angstroms.
- Default: None (no shift).
- :param lambda_bins (array): array of wavelengths in angstroms.
- Default: None
- :param n_lambda (float): size of lambda bins to create
- array if lambda_bins is None. Default: 12000
- :param dlambda (float): lambda bin width if lambda_bins is
- None. Default: 0.01
+
+ Parameters
+ ----------
+
+ lamba_0 : float YTQuantity in length units
+ central wavelength.
+ f_value : float
+ absorption line f-value.
+ gamma : float
+ absorption line gamma value.
+ v_doppler : float YTQuantity in velocity units
+ doppler b-parameter.
+ column_density : float YTQuantity in (length units)^-2
+ column density.
+ delta_v : float YTQuantity in velocity units
+ velocity offset from lamba_0.
+ Default: None (no shift).
+ delta_lambda : float YTQuantity in length units
+ wavelength offset.
+ Default: None (no shift).
+ lambda_bins : YTArray in length units
+ wavelength array for line deposition. If None, one will be
+ created using n_lambda and dlambda.
+ Default: None.
+ n_lambda : int
+ size of lambda bins to create if lambda_bins is None.
+ Default: 12000.
+ dlambda : float
+ lambda bin width in angstroms if lambda_bins is None.
+ Default: 0.01.
+
"""
- ## constants
- me = mass_electron_cgs # grams mass electron
- e = charge_proton_cgs # esu
- c = speed_of_light_cgs * km_per_cm # km/s
- ccgs = speed_of_light_cgs # cm/s
-
- ## shift lam0 by deltav
- if deltav is not None:
- lam1 = lam0 * (1 + deltav / ccgs)
+ ## shift lamba_0 by delta_v
+ if delta_v is not None:
+ lam1 = lamba_0 * (1 + delta_v / speed_of_light_cgs)
elif delta_lambda is not None:
- lam1 = lam0 + delta_lambda
+ lam1 = lamba_0 + delta_lambda
else:
- lam1 = lam0
+ lam1 = lamba_0
## conversions
- vdop = vkms.in_units("cm/s") # in cm/s
- lam0cgs = lam0.in_units("cm") # rest wavelength in cm
- lam1cgs = lam1.in_units("cm") # line wavelength in cm
- nu1 = (ccgs / lam1cgs).in_units("1/s") # line freq in Hz
- nudop = (vdop / ccgs * nu1).in_units("1/s") # doppler width in Hz
- lamdop = (vdop / ccgs * lam1).in_units("angstrom") # doppler width in Ang
+ nu1 = speed_of_light_cgs / lam1 # line freq in Hz
+ nudop = v_doppler / speed_of_light_cgs * nu1 # doppler width in Hz
+ lamdop = v_doppler / speed_of_light_cgs * lam1 # doppler width in Ang
- #import pdb ; pdb.set_trace()
-
## create wavelength
if lambda_bins is None:
lambda_bins = lam1 + \
np.arange(n_lambda, dtype=np.float) * dlambda - \
- n_lambda * dlambda / 2 # wavelength vector (angstroms)
- nua = (ccgs / lambda_bins).in_units("1/s") # frequency vector (Hz)
+ n_lambda * dlambda / 2 # wavelength vector (angstroms)
+ nua = (speed_of_light_cgs / lambda_bins) # frequency vector (Hz)
## tau_0
- tau_X = np.sqrt(np.pi) * e**2 / (me * ccgs) * \
- column_density * fval / vdop
- tau0 = tau_X * lam0cgs
+ tau_X = np.sqrt(np.pi) * charge_proton_cgs**2 / \
+ (mass_electron_cgs * speed_of_light_cgs) * \
+ column_density * f_value / v_doppler
+ tau0 = tau_X * lamba_0
# dimensionless frequency offset in units of doppler freq
- x = (nua - nu1) / nudop
- a = gamma / (4 * np.pi * nudop) # damping parameter
- phi = voigt(a, x) # profile
- tauphi = tau0 * phi # profile scaled with tau0
+ x = ((nua - nu1) / nudop).in_units("")
+ a = (gamma / (4 * np.pi * nudop)).in_units("s") # damping parameter
+ phi = voigt(a, x) # line profile
+ tauphi = (tau0 * phi).in_units("") # profile scaled with tau0
return (lambda_bins, tauphi)
https://bitbucket.org/yt_analysis/yt/commits/400801f9a26b/
Changeset: 400801f9a26b
Branch: yt
User: brittonsmith
Date: 2015-05-21 18:13:04+00:00
Summary: Fixing a tab.
Affected #: 1 file
diff -r 287d5d45255bc8fd3a2ea707e5cdb8dd89daeb2f -r 400801f9a26b821b668e99c355503527eaabeb8f yt/analysis_modules/absorption_spectrum/absorption_spectrum.py
--- a/yt/analysis_modules/absorption_spectrum/absorption_spectrum.py
+++ b/yt/analysis_modules/absorption_spectrum/absorption_spectrum.py
@@ -193,7 +193,7 @@
column_density[lixel] / continuum['normalization']
self.tau_field[left_index[lixel]:right_index[lixel]] += line_tau
pbar.update(i)
- pbar.finish()
+ pbar.finish()
def _add_lines_to_spectrum(self, field_data, use_peculiar_velocity):
"""
https://bitbucket.org/yt_analysis/yt/commits/515a01ad24f1/
Changeset: 515a01ad24f1
Branch: yt
User: brittonsmith
Date: 2015-05-21 19:20:37+00:00
Summary: It was right all along.
Affected #: 1 file
diff -r 400801f9a26b821b668e99c355503527eaabeb8f -r 515a01ad24f13f22619c3226c8f7d72e2d1a945e yt/analysis_modules/absorption_spectrum/absorption_spectrum.py
--- a/yt/analysis_modules/absorption_spectrum/absorption_spectrum.py
+++ b/yt/analysis_modules/absorption_spectrum/absorption_spectrum.py
@@ -193,7 +193,7 @@
column_density[lixel] / continuum['normalization']
self.tau_field[left_index[lixel]:right_index[lixel]] += line_tau
pbar.update(i)
- pbar.finish()
+ pbar.finish()
def _add_lines_to_spectrum(self, field_data, use_peculiar_velocity):
"""
https://bitbucket.org/yt_analysis/yt/commits/2276a782446a/
Changeset: 2276a782446a
Branch: yt
User: brittonsmith
Date: 2015-05-24 15:56:09+00:00
Summary: Atomic masses now are YTQuantities.
Affected #: 1 file
diff -r 515a01ad24f13f22619c3226c8f7d72e2d1a945e -r 2276a782446a539fc9f925f9ef97cd8a87f58156 yt/analysis_modules/absorption_spectrum/absorption_spectrum.py
--- a/yt/analysis_modules/absorption_spectrum/absorption_spectrum.py
+++ b/yt/analysis_modules/absorption_spectrum/absorption_spectrum.py
@@ -23,7 +23,7 @@
from yt.funcs import get_pbar, mylog
from yt.units.yt_array import YTArray, YTQuantity
from yt.utilities.physical_constants import \
- amu_cgs, boltzmann_constant_cgs, \
+ boltzmann_constant_cgs, \
speed_of_light_cgs
from yt.utilities.on_demand_imports import _astropy
@@ -79,7 +79,8 @@
self.line_list.append({'label': label, 'field_name': field_name,
'wavelength': YTQuantity(wavelength, "angstrom"),
'f_value': f_value,
- 'gamma': gamma, 'atomic_mass': atomic_mass,
+ 'gamma': gamma,
+ 'atomic_mass': YTQuantity(atomic_mass, "amu"),
'label_threshold': label_threshold})
def add_continuum(self, label, field_name, wavelength,
@@ -213,7 +214,7 @@
field_data['velocity_los'] / speed_of_light_cgs
thermal_b = np.sqrt((2 * boltzmann_constant_cgs *
field_data['temperature']) /
- (amu_cgs * line['atomic_mass']))
+ line['atomic_mass'])
center_bins = np.digitize((delta_lambda + line['wavelength']),
self.lambda_bins)
https://bitbucket.org/yt_analysis/yt/commits/d26b1b67dee3/
Changeset: d26b1b67dee3
Branch: yt
User: brittonsmith
Date: 2015-05-24 19:55:26+00:00
Summary: Updating cookbook recipe.
Affected #: 1 file
diff -r 2276a782446a539fc9f925f9ef97cd8a87f58156 -r d26b1b67dee311dfe1b6b4f4384d1eb036c03ea3 doc/source/cookbook/fit_spectrum.py
--- a/doc/source/cookbook/fit_spectrum.py
+++ b/doc/source/cookbook/fit_spectrum.py
@@ -100,6 +100,6 @@
order_fits = ['OVI', 'HI']
# Fit spectrum and save fit
-fitted_lines, fitted_flux = generate_total_fit(wavelength,
+fitted_lines, fitted_flux = generate_total_fit(wavelength.d,
flux, order_fits, species_dicts,
output_file='spectrum_fit.h5')
https://bitbucket.org/yt_analysis/yt/commits/5e426c5aaf38/
Changeset: 5e426c5aaf38
Branch: yt
User: brittonsmith
Date: 2015-05-29 19:46:17+00:00
Summary: Undoing change to cookbook.
Affected #: 1 file
diff -r d26b1b67dee311dfe1b6b4f4384d1eb036c03ea3 -r 5e426c5aaf389407192b07f07f0a7cb0d9cbabba doc/source/cookbook/fit_spectrum.py
--- a/doc/source/cookbook/fit_spectrum.py
+++ b/doc/source/cookbook/fit_spectrum.py
@@ -100,6 +100,6 @@
order_fits = ['OVI', 'HI']
# Fit spectrum and save fit
-fitted_lines, fitted_flux = generate_total_fit(wavelength.d,
+fitted_lines, fitted_flux = generate_total_fit(wavelength,
flux, order_fits, species_dicts,
output_file='spectrum_fit.h5')
https://bitbucket.org/yt_analysis/yt/commits/ee1cb9fa80bd/
Changeset: ee1cb9fa80bd
Branch: yt
User: brittonsmith
Date: 2015-05-29 19:50:48+00:00
Summary: Stripping units on wavelength array if they're there.
Affected #: 1 file
diff -r 5e426c5aaf389407192b07f07f0a7cb0d9cbabba -r ee1cb9fa80bd207896d77470a1365469521a7b56 yt/analysis_modules/absorption_spectrum/absorption_spectrum_fit.py
--- a/yt/analysis_modules/absorption_spectrum/absorption_spectrum_fit.py
+++ b/yt/analysis_modules/absorption_spectrum/absorption_spectrum_fit.py
@@ -79,6 +79,10 @@
absorption profiles. Same size as x.
"""
+ # convert to NumPy array if we have a YTArray
+ if isinstance(x, YTArray):
+ x = x.d
+
#Empty dictionary for fitted lines
allSpeciesLines = {}
https://bitbucket.org/yt_analysis/yt/commits/691da9acc759/
Changeset: 691da9acc759
Branch: yt
User: brittonsmith
Date: 2015-05-29 20:49:05+00:00
Summary: Importing YTArray and cleaning up imports.
Affected #: 1 file
diff -r ee1cb9fa80bd207896d77470a1365469521a7b56 -r 691da9acc75937e1ba62aa3620c7858130d82d7c yt/analysis_modules/absorption_spectrum/absorption_spectrum_fit.py
--- a/yt/analysis_modules/absorption_spectrum/absorption_spectrum_fit.py
+++ b/yt/analysis_modules/absorption_spectrum/absorption_spectrum_fit.py
@@ -1,9 +1,14 @@
-from __future__ import print_function
+import h5py
import numpy as np
-import h5py
-from yt.analysis_modules.absorption_spectrum.absorption_line \
- import voigt
-from yt.utilities.on_demand_imports import _scipy
+
+from yt.analysis_modules.absorption_spectrum.absorption_line import \
+ voigt
+from yt.funcs import \
+ mylog
+from yt.units.yt_array import \
+ YTArray
+from yt.utilities.on_demand_imports import \
+ _scipy
optimize = _scipy.optimize
@@ -1011,6 +1016,5 @@
f.create_dataset("{0}/b".format(ion),data=params['b'])
f.create_dataset("{0}/z".format(ion),data=params['z'])
f.create_dataset("{0}/complex".format(ion),data=params['group#'])
- print('Writing spectrum fit to {0}'.format(file_name))
+ mylog.info('Writing spectrum fit to {0}'.format(file_name))
f.close()
-
https://bitbucket.org/yt_analysis/yt/commits/1dee4eb06343/
Changeset: 1dee4eb06343
Branch: yt
User: xarthisius
Date: 2015-06-01 15:38:09+00:00
Summary: Merged in brittonsmith/yt (pull request #1596)
Fixing unit support in AbsorptionSpectrum
Affected #: 4 files
diff -r d6983e51e584eb0e996dfd2a4475ff2a8e0a6014 -r 1dee4eb0634354c90ace1ba14a6d2b4ed4237882 yt/analysis_modules/absorption_spectrum/absorption_line.py
--- a/yt/analysis_modules/absorption_spectrum/absorption_line.py
+++ b/yt/analysis_modules/absorption_spectrum/absorption_line.py
@@ -16,12 +16,9 @@
import numpy as np
from yt.utilities.physical_constants import \
charge_proton_cgs, \
- cm_per_km, \
- km_per_cm, \
mass_electron_cgs, \
speed_of_light_cgs
-
def voigt(a,u):
"""
NAME:
@@ -140,67 +137,75 @@
k1 = k1.astype(np.float64).clip(0)
return k1
-def tau_profile(lam0, fval, gamma, vkms, column_density,
- deltav=None, delta_lambda=None,
+def tau_profile(lamba_0, f_value, gamma, v_doppler, column_density,
+ delta_v=None, delta_lambda=None,
lambda_bins=None, n_lambda=12000, dlambda=0.01):
- """
+ r"""
Create an optical depth vs. wavelength profile for an
absorption line using a voigt profile.
- :param lam0 (float): central wavelength (angstroms).
- :param fval (float): f-value.
- :param gamma (float): gamma value.
- :param vkms (float): doppler b-parameter.
- :param column_density (float): column density (cm^-2).
- :param deltav (float): velocity offset from lam0 (km/s).
- Default: None (no shift).
- :param delta_lambda (float): wavelength offset in angstroms.
- Default: None (no shift).
- :param lambda_bins (array): array of wavelengths in angstroms.
- Default: None
- :param n_lambda (float): size of lambda bins to create
- array if lambda_bins is None. Default: 12000
- :param dlambda (float): lambda bin width if lambda_bins is
- None. Default: 0.01
+
+ Parameters
+ ----------
+
+ lamba_0 : float YTQuantity in length units
+ central wavelength.
+ f_value : float
+ absorption line f-value.
+ gamma : float
+ absorption line gamma value.
+ v_doppler : float YTQuantity in velocity units
+ doppler b-parameter.
+ column_density : float YTQuantity in (length units)^-2
+ column density.
+ delta_v : float YTQuantity in velocity units
+ velocity offset from lamba_0.
+ Default: None (no shift).
+ delta_lambda : float YTQuantity in length units
+ wavelength offset.
+ Default: None (no shift).
+ lambda_bins : YTArray in length units
+ wavelength array for line deposition. If None, one will be
+ created using n_lambda and dlambda.
+ Default: None.
+ n_lambda : int
+ size of lambda bins to create if lambda_bins is None.
+ Default: 12000.
+ dlambda : float
+ lambda bin width in angstroms if lambda_bins is None.
+ Default: 0.01.
+
"""
- ## constants
- me = mass_electron_cgs # grams mass electron
- e = charge_proton_cgs # esu
- c = speed_of_light_cgs * km_per_cm # km/s
- ccgs = speed_of_light_cgs # cm/s
-
- ## shift lam0 by deltav
- if deltav is not None:
- lam1 = lam0 * (1 + deltav / c)
+ ## shift lamba_0 by delta_v
+ if delta_v is not None:
+ lam1 = lamba_0 * (1 + delta_v / speed_of_light_cgs)
elif delta_lambda is not None:
- lam1 = lam0 + delta_lambda
+ lam1 = lamba_0 + delta_lambda
else:
- lam1 = lam0
+ lam1 = lamba_0
## conversions
- vdop = vkms * cm_per_km # in cm/s
- lam0cgs = lam0 / 1.e8 # rest wavelength in cm
- lam1cgs = lam1 / 1.e8 # line wavelength in cm
- nu1 = ccgs / lam1cgs # line freq in Hz
- nudop = vdop / ccgs * nu1 # doppler width in Hz
- lamdop = vdop / ccgs * lam1 # doppler width in Ang
+ nu1 = speed_of_light_cgs / lam1 # line freq in Hz
+ nudop = v_doppler / speed_of_light_cgs * nu1 # doppler width in Hz
+ lamdop = v_doppler / speed_of_light_cgs * lam1 # doppler width in Ang
## create wavelength
if lambda_bins is None:
lambda_bins = lam1 + \
np.arange(n_lambda, dtype=np.float) * dlambda - \
- n_lambda * dlambda / 2 # wavelength vector (angstroms)
- nua = ccgs / (lambda_bins / 1.e8) # frequency vector (Hz)
+ n_lambda * dlambda / 2 # wavelength vector (angstroms)
+ nua = (speed_of_light_cgs / lambda_bins) # frequency vector (Hz)
## tau_0
- tau_X = np.sqrt(np.pi) * e**2 / (me * ccgs) * \
- column_density * fval / vdop
- tau0 = tau_X * lam0cgs
+ tau_X = np.sqrt(np.pi) * charge_proton_cgs**2 / \
+ (mass_electron_cgs * speed_of_light_cgs) * \
+ column_density * f_value / v_doppler
+ tau0 = tau_X * lamba_0
# dimensionless frequency offset in units of doppler freq
- x = (nua - nu1) / nudop
- a = gamma / (4 * np.pi * nudop) # damping parameter
- phi = voigt(a, x) # profile
- tauphi = tau0 * phi # profile scaled with tau0
+ x = ((nua - nu1) / nudop).in_units("")
+ a = (gamma / (4 * np.pi * nudop)).in_units("s") # damping parameter
+ phi = voigt(a, x) # line profile
+ tauphi = (tau0 * phi).in_units("") # profile scaled with tau0
return (lambda_bins, tauphi)
diff -r d6983e51e584eb0e996dfd2a4475ff2a8e0a6014 -r 1dee4eb0634354c90ace1ba14a6d2b4ed4237882 yt/analysis_modules/absorption_spectrum/absorption_spectrum.py
--- a/yt/analysis_modules/absorption_spectrum/absorption_spectrum.py
+++ b/yt/analysis_modules/absorption_spectrum/absorption_spectrum.py
@@ -4,7 +4,7 @@
"""
-from __future__ import print_function
+
from __future__ import absolute_import
#-----------------------------------------------------------------------------
@@ -21,13 +21,12 @@
from .absorption_line import tau_profile
from yt.funcs import get_pbar, mylog
-from yt.units.yt_array import YTArray
+from yt.units.yt_array import YTArray, YTQuantity
from yt.utilities.physical_constants import \
- amu_cgs, boltzmann_constant_cgs, \
- speed_of_light_cgs, km_per_cm
+ boltzmann_constant_cgs, \
+ speed_of_light_cgs
from yt.utilities.on_demand_imports import _astropy
-speed_of_light_kms = speed_of_light_cgs * km_per_cm
pyfits = _astropy.pyfits
class AbsorptionSpectrum(object):
@@ -49,8 +48,10 @@
self.tau_field = None
self.flux_field = None
self.spectrum_line_list = None
- self.lambda_bins = np.linspace(lambda_min, lambda_max, n_lambda)
- self.bin_width = (lambda_max - lambda_min) / float(n_lambda - 1)
+ self.lambda_bins = YTArray(np.linspace(lambda_min, lambda_max, n_lambda),
+ "angstrom")
+ self.bin_width = YTQuantity((lambda_max - lambda_min) /
+ float(n_lambda - 1), "angstrom")
self.line_list = []
self.continuum_list = []
@@ -76,8 +77,10 @@
mass of atom in amu.
"""
self.line_list.append({'label': label, 'field_name': field_name,
- 'wavelength': wavelength, 'f_value': f_value,
- 'gamma': gamma, 'atomic_mass': atomic_mass,
+ 'wavelength': YTQuantity(wavelength, "angstrom"),
+ 'f_value': f_value,
+ 'gamma': gamma,
+ 'atomic_mass': YTQuantity(atomic_mass, "amu"),
'label_threshold': label_threshold})
def add_continuum(self, label, field_name, wavelength,
@@ -209,16 +212,15 @@
# include factor of (1 + z) because our velocity is in proper frame.
delta_lambda += line['wavelength'] * (1 + field_data['redshift']) * \
field_data['velocity_los'] / speed_of_light_cgs
- thermal_b = km_per_cm * np.sqrt((2 * boltzmann_constant_cgs *
- field_data['temperature']) /
- (amu_cgs * line['atomic_mass']))
- thermal_b.convert_to_cgs()
+ thermal_b = np.sqrt((2 * boltzmann_constant_cgs *
+ field_data['temperature']) /
+ line['atomic_mass'])
center_bins = np.digitize((delta_lambda + line['wavelength']),
self.lambda_bins)
# ratio of line width to bin width
width_ratio = ((line['wavelength'] + delta_lambda) * \
- thermal_b / speed_of_light_kms / self.bin_width).value
+ thermal_b / speed_of_light_cgs / self.bin_width).in_units("").d
if (width_ratio < 1.0).any():
mylog.warn(("%d out of %d line components are unresolved, " +
@@ -240,11 +242,13 @@
my_bin_ratio = spectrum_bin_ratio
while True:
lambda_bins, line_tau = \
- tau_profile(line['wavelength'], line['f_value'],
- line['gamma'], thermal_b[lixel],
- column_density[lixel],
- delta_lambda=delta_lambda[lixel],
- lambda_bins=self.lambda_bins[left_index[lixel]:right_index[lixel]])
+ tau_profile(
+ line['wavelength'], line['f_value'],
+ line['gamma'], thermal_b[lixel].in_units("km/s"),
+ column_density[lixel],
+ delta_lambda=delta_lambda[lixel],
+ lambda_bins=self.lambda_bins[left_index[lixel]:right_index[lixel]])
+
# Widen wavelength window until optical depth reaches a max value at the ends.
if (line_tau[0] < max_tau and line_tau[-1] < max_tau) or \
(left_index[lixel] <= 0 and right_index[lixel] >= self.n_lambda):
@@ -260,7 +264,7 @@
if line['label_threshold'] is not None and \
column_density[lixel] >= line['label_threshold']:
if use_peculiar_velocity:
- peculiar_velocity = km_per_cm * field_data['velocity_los'][lixel]
+ peculiar_velocity = field_data['velocity_los'][lixel].in_units("km/s")
else:
peculiar_velocity = 0.0
self.spectrum_line_list.append({'label': line['label'],
@@ -271,7 +275,7 @@
'redshift': field_data['redshift'][lixel],
'v_pec': peculiar_velocity})
pbar.update(i)
- pbar.finish()
+ pbar.finish()
del column_density, delta_lambda, thermal_b, \
center_bins, width_ratio, left_index, right_index
@@ -280,7 +284,7 @@
"""
Write out list of spectral lines.
"""
- print("Writing spectral line list: %s." % filename)
+ mylog.info("Writing spectral line list: %s." % filename)
self.spectrum_line_list.sort(key=lambda obj: obj['wavelength'])
f = open(filename, 'w')
f.write('#%-14s %-14s %-12s %-12s %-12s %-12s\n' %
@@ -295,7 +299,7 @@
"""
Write spectrum to an ascii file.
"""
- print("Writing spectrum to ascii file: %s." % filename)
+ mylog.info("Writing spectrum to ascii file: %s." % filename)
f = open(filename, 'w')
f.write("# wavelength[A] tau flux\n")
for i in range(self.lambda_bins.size):
@@ -307,7 +311,7 @@
"""
Write spectrum to a fits file.
"""
- print("Writing spectrum to fits file: %s." % filename)
+ mylog.info("Writing spectrum to fits file: %s." % filename)
col1 = pyfits.Column(name='wavelength', format='E', array=self.lambda_bins)
col2 = pyfits.Column(name='flux', format='E', array=self.flux_field)
cols = pyfits.ColDefs([col1, col2])
@@ -319,7 +323,7 @@
Write spectrum to an hdf5 file.
"""
- print("Writing spectrum to hdf5 file: %s." % filename)
+ mylog.info("Writing spectrum to hdf5 file: %s." % filename)
output = h5py.File(filename, 'w')
output.create_dataset('wavelength', data=self.lambda_bins)
output.create_dataset('tau', data=self.tau_field)
diff -r d6983e51e584eb0e996dfd2a4475ff2a8e0a6014 -r 1dee4eb0634354c90ace1ba14a6d2b4ed4237882 yt/analysis_modules/absorption_spectrum/absorption_spectrum_fit.py
--- a/yt/analysis_modules/absorption_spectrum/absorption_spectrum_fit.py
+++ b/yt/analysis_modules/absorption_spectrum/absorption_spectrum_fit.py
@@ -1,9 +1,14 @@
-from __future__ import print_function
+import h5py
import numpy as np
-import h5py
-from yt.analysis_modules.absorption_spectrum.absorption_line \
- import voigt
-from yt.utilities.on_demand_imports import _scipy
+
+from yt.analysis_modules.absorption_spectrum.absorption_line import \
+ voigt
+from yt.funcs import \
+ mylog
+from yt.units.yt_array import \
+ YTArray
+from yt.utilities.on_demand_imports import \
+ _scipy
optimize = _scipy.optimize
@@ -79,6 +84,10 @@
absorption profiles. Same size as x.
"""
+ # convert to NumPy array if we have a YTArray
+ if isinstance(x, YTArray):
+ x = x.d
+
#Empty dictionary for fitted lines
allSpeciesLines = {}
@@ -1007,6 +1016,5 @@
f.create_dataset("{0}/b".format(ion),data=params['b'])
f.create_dataset("{0}/z".format(ion),data=params['z'])
f.create_dataset("{0}/complex".format(ion),data=params['group#'])
- print('Writing spectrum fit to {0}'.format(file_name))
+ mylog.info('Writing spectrum fit to {0}'.format(file_name))
f.close()
-
Repository URL: https://bitbucket.org/yt_analysis/yt/
--
This is a commit notification from bitbucket.org. You are receiving
this because you have the service enabled, addressing the recipient of
this email.
More information about the yt-svn
mailing list