[yt-svn] commit/yt: 4 new changesets
commits-noreply at bitbucket.org
commits-noreply at bitbucket.org
Wed Apr 20 11:20:37 PDT 2016
4 new commits in yt:
https://bitbucket.org/yt_analysis/yt/commits/5595d7cdcb46/
Changeset: 5595d7cdcb46
Branch: yt
User: migueldvb
Date: 2016-03-14 20:21:10+00:00
Summary: remove trailing space
Affected #: 1 file
diff -r e1b2e0a5997f0cce7ff52547347e9e5214ad0d5c -r 5595d7cdcb461ff2f15a4b931c7065e6d11d9fbd 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
@@ -12,14 +12,14 @@
optimize = _scipy.optimize
-def generate_total_fit(x, fluxData, orderFits, speciesDicts,
+def generate_total_fit(x, fluxData, orderFits, speciesDicts,
minError=1E-4, complexLim=.995,
- fitLim=.97, minLength=3,
+ fitLim=.97, minLength=3,
maxLength=1000, splitLim=.99,
output_file=None):
"""
- This function is designed to fit an absorption spectrum by breaking
+ This function is designed to fit an absorption spectrum by breaking
the spectrum up into absorption complexes, and iteratively adding
and optimizing voigt profiles to each complex.
@@ -31,46 +31,46 @@
array of flux corresponding to the wavelengths given
in x. (needs to be the same size as x)
orderFits : list
- list of the names of the species in the order that they
+ list of the names of the species in the order that they
should be fit. Names should correspond to the names of the species
given in speciesDicts. (ex: ['lya','OVI'])
speciesDicts : dictionary
Dictionary of dictionaries (I'm addicted to dictionaries, I
confess). Top level keys should be the names of all the species given
- in orderFits. The entries should be dictionaries containing all
- relevant parameters needed to create an absorption line of a given
+ in orderFits. The entries should be dictionaries containing all
+ relevant parameters needed to create an absorption line of a given
species (f,Gamma,lambda0) as well as max and min values for parameters
to be fit
complexLim : float, optional
- Maximum flux to start the edge of an absorption complex. Different
- from fitLim because it decides extent of a complex rather than
- whether or not a complex is accepted.
+ Maximum flux to start the edge of an absorption complex. Different
+ from fitLim because it decides extent of a complex rather than
+ whether or not a complex is accepted.
fitLim : float,optional
- Maximum flux where the level of absorption will trigger
+ Maximum flux where the level of absorption will trigger
identification of the region as an absorption complex. Default = .98.
(ex: for a minSize=.98, a region where all the flux is between 1.0 and
.99 will not be separated out to be fit as an absorbing complex, but
a region that contains a point where the flux is .97 will be fit
as an absorbing complex.)
minLength : int, optional
- number of cells required for a complex to be included.
+ number of cells required for a complex to be included.
default is 3 cells.
maxLength : int, optional
number of cells required for a complex to be split up. Default
is 1000 cells.
splitLim : float, optional
if attempting to split a region for being larger than maxlength
- the point of the split must have a flux greater than splitLim
+ the point of the split must have a flux greater than splitLim
(ie: absorption greater than splitLim). Default= .99.
output_file : string, optional
- location to save the results of the fit.
+ location to save the results of the fit.
Returns
-------
allSpeciesLines : dictionary
- Dictionary of dictionaries representing the fit lines.
+ Dictionary of dictionaries representing the fit lines.
Top level keys are the species given in orderFits and the corresponding
- entries are dictionaries with the keys 'N','b','z', and 'group#'.
+ entries are dictionaries with the keys 'N','b','z', and 'group#'.
Each of these corresponds to a list of the parameters for every
accepted fitted line. (ie: N[0],b[0],z[0] will create a line that
fits some part of the absorption spectrum). 'group#' is a similar list
@@ -87,7 +87,7 @@
# convert to NumPy array if we have a YTArray
if isinstance(x, YTArray):
x = x.d
-
+
#Empty dictionary for fitted lines
allSpeciesLines = {}
@@ -129,8 +129,8 @@
#Check if any flux at partner sites
if not _line_exists(speciesDict['wavelength'],
- fluxData,z,x0,xRes,fitLim):
- continue
+ fluxData,z,x0,xRes,fitLim):
+ continue
#Fit Using complex tools
newLinesP,flag=_complex_fit(xBounded,yDatBounded,yFitBounded,
@@ -139,7 +139,7 @@
#If flagged as a bad fit, species is lyman alpha,
# and it may be a saturated line, use special tools
if flag and species=='lya' and min(yDatBounded)<.1:
- newLinesP=_large_flag_fit(xBounded,yDatBounded,
+ newLinesP=_large_flag_fit(xBounded,yDatBounded,
yFitBounded,z,speciesDict,
fitLim,minError)
@@ -151,7 +151,7 @@
#Check existence of partner lines if applicable
if len(speciesDict['wavelength']) != 1:
- newLinesP = _remove_unaccepted_partners(newLinesP, x, fluxData,
+ newLinesP = _remove_unaccepted_partners(newLinesP, x, fluxData,
b, minError, x0, xRes, speciesDict)
@@ -177,15 +177,15 @@
return (allSpeciesLines,yFit)
-def _complex_fit(x, yDat, yFit, initz, minSize, errBound, speciesDict,
+def _complex_fit(x, yDat, yFit, initz, minSize, errBound, speciesDict,
initP=None):
""" Fit an absorption complex by iteratively adding and optimizing
voigt profiles.
-
+
A complex is defined as a region where some number of lines may be present,
or a region of non zero of absorption. Lines are iteratively added
and optimized until the difference between the flux generated using
- the optimized parameters has a least squares difference between the
+ the optimized parameters has a least squares difference between the
desired flux profile less than the error bound.
Parameters
@@ -199,7 +199,7 @@
array of flux profile fitted for the wavelength
space given by x already. Same size as x.
initz : float
- redshift to try putting first line at
+ redshift to try putting first line at
(maximum absorption for region)
minsize : float
minimum absorption allowed for a line to still count as a line
@@ -219,7 +219,7 @@
Returns
-------
linesP : (3,) ndarray
- Array of best parameters if a good enough fit is found in
+ Array of best parameters if a good enough fit is found in
the form [[N1,b1,z1], [N2,b2,z2],...]
flag : bool
boolean value indicating the success of the fit (True if unsuccessful)
@@ -227,8 +227,8 @@
#Setup initial line guesses
if initP is None: #Regular fit
- initP = [0,0,0]
- if min(yDat)<.01: #Large lines get larger initial guess
+ initP = [0,0,0]
+ if min(yDat)<.01: #Large lines get larger initial guess
initP[0] = speciesDict['init_N']*10**2
elif min(yDat)<.5:
initP[0] = speciesDict['init_N']*10**1
@@ -246,9 +246,9 @@
wl0 = speciesDict['wavelength'][0]
#Check if first line exists still
- if min(yDat-yFit+1)>minSize:
+ if min(yDat-yFit+1)>minSize:
return [],False
-
+
#Values to proceed through first run
errSq,prevErrSq,prevLinesP=1,10*len(x),[]
@@ -283,9 +283,9 @@
if any(linesP[:,1]==speciesDict['init_b']):
flag = True
break
-
+
#If good enough, break
- if errSq < errBound:
+ if errSq < errBound:
break
#If last fit was worse, reject the last line and revert to last fit
@@ -297,7 +297,7 @@
linesP = prevLinesP
break
- #If too many lines
+ #If too many lines
if np.shape(linesP)[0]>8 or np.size(linesP)+3>=len(x):
#If its fitable by flag tools and still bad, use flag tools
if errSq >1E2*errBound and speciesDict['name']=='HI lya':
@@ -311,11 +311,11 @@
prevLinesP = linesP
#Set up initial condition for new line
- newP = [0,0,0]
+ newP = [0,0,0]
yAdjusted = 1+yFit*yNewFit-yDat
-
- if min(yAdjusted)<.01: #Large lines get larger initial guess
+
+ if min(yAdjusted)<.01: #Large lines get larger initial guess
newP[0] = speciesDict['init_N']*10**2
elif min(yAdjusted)<.5:
newP[0] = speciesDict['init_N']*10**1
@@ -333,7 +333,7 @@
remove=[]
for i,p in enumerate(linesP):
check=_check_params(np.array([p]),speciesDict,x)
- if check:
+ if check:
remove.append(i)
linesP = np.delete(linesP,remove,axis=0)
@@ -358,7 +358,7 @@
array of flux profile fitted for the wavelength
space given by x already. Same size as x.
initz : float
- redshift to try putting first line at
+ redshift to try putting first line at
(maximum absorption for region)
speciesDict : dictionary
dictionary containing all relevant parameters needed
@@ -374,15 +374,15 @@
Returns
-------
bestP : (3,) ndarray
- array of best parameters if a good enough fit is found in
- the form [[N1,b1,z1], [N2,b2,z2],...]
+ array of best parameters if a good enough fit is found in
+ the form [[N1,b1,z1], [N2,b2,z2],...]
"""
#Set up some initial line guesses
lineTests = _get_test_lines(initz)
#Keep track of the lowest achieved error
- bestError = 1000
+ bestError = 1000
#Iterate through test line guesses
for initLines in lineTests:
@@ -403,7 +403,7 @@
bestError = errSq
bestP = linesP
- if bestError>10*errBound*len(x):
+ if bestError>10*errBound*len(x):
return []
else:
return bestP
@@ -411,11 +411,11 @@
def _get_test_lines(initz):
"""
Returns a 3d numpy array of lines to test as initial guesses for difficult
- to fit lyman alpha absorbers that are saturated.
-
+ to fit lyman alpha absorbers that are saturated.
+
The array is 3d because
the first dimension gives separate initial guesses, the second dimension
- has multiple lines for the same guess (trying a broad line plus a
+ has multiple lines for the same guess (trying a broad line plus a
saturated line) and the 3d dimension contains the 3 fit parameters (N,b,z)
Parameters
@@ -426,7 +426,7 @@
Returns
-------
testP : (,3,) ndarray
- numpy array of the form
+ numpy array of the form
[[[N1a,b1a,z1a], [N1b,b1b,z1b]], [[N2a,b2,z2a],...] ...]
"""
@@ -454,8 +454,8 @@
return testP
def _get_bounds(z, b, wl, x0, xRes):
- """
- Gets the indices of range of wavelength that the wavelength wl is in
+ """
+ Gets the indices of range of wavelength that the wavelength wl is in
with the size of some initial wavelength range.
Used for checking if species with multiple lines (as in the OVI doublet)
@@ -466,21 +466,21 @@
z : float
redshift
b : (3) ndarray/list
- initial bounds in form [i0,i1,i2] where i0 is the index of the
- minimum flux for the complex, i1 is index of the lower wavelength
+ initial bounds in form [i0,i1,i2] where i0 is the index of the
+ minimum flux for the complex, i1 is index of the lower wavelength
edge of the complex, and i2 is the index of the higher wavelength
edge of the complex.
wl : float
- unredshifted wavelength of the peak of the new region
+ unredshifted wavelength of the peak of the new region
x0 : float
wavelength of the index 0
xRes : float
difference in wavelength for two consecutive indices
-
+
Returns
-------
indices : (2) tuple
- Tuple (i1,i2) where i1 is the index of the lower wavelength bound of
+ Tuple (i1,i2) where i1 is the index of the lower wavelength bound of
the new region and i2 is the index of the higher wavelength bound of
the new region
"""
@@ -492,7 +492,7 @@
return indices
-def _remove_unaccepted_partners(linesP, x, y, b, errBound,
+def _remove_unaccepted_partners(linesP, x, y, b, errBound,
x0, xRes, speciesDict):
"""
Given a set of parameters [N,b,z] that form multiple lines for a given
@@ -507,18 +507,18 @@
Parameters
----------
linesP : (3,) ndarray
- array giving sets of line parameters in
+ array giving sets of line parameters in
form [[N1, b1, z1], ...]
x : (N) ndarray
wavelength array [nm]
y : (N) ndarray
normalized flux array of original data
b : (3) tuple/list/ndarray
- indices that give the bounds of the original region so that another
+ indices that give the bounds of the original region so that another
region of similar size can be used to determine the goodness
of fit of the other wavelengths
errBound : float
- size of the error that is appropriate for a given region,
+ size of the error that is appropriate for a given region,
adjusted to account for the size of the region.
Returns
@@ -569,7 +569,7 @@
#Remove all bad line fits
linesP = np.delete(linesP,removeLines,axis=0)
- return linesP
+ return linesP
@@ -588,12 +588,12 @@
xRes : float
difference in wavelength between consecutive cells in flux array
fluxMin : float
- maximum flux to count as a line existing.
+ maximum flux to count as a line existing.
Returns
-------
- flag : boolean
+ flag : boolean
value indicating whether all lines exist. True if all lines exist
"""
@@ -618,7 +618,7 @@
def _find_complexes(x, yDat, complexLim=.999, fitLim=.99,
minLength =3, maxLength=1000, splitLim=.99):
"""Breaks up the wavelength space into groups
- where there is some absorption.
+ where there is some absorption.
Parameters
----------
@@ -628,33 +628,33 @@
array of flux corresponding to the wavelengths given
in x. (needs to be the same size as x)
complexLim : float, optional
- Maximum flux to start the edge of an absorption complex. Different
- from fitLim because it decides extent of a complex rather than
- whether or not a complex is accepted.
+ Maximum flux to start the edge of an absorption complex. Different
+ from fitLim because it decides extent of a complex rather than
+ whether or not a complex is accepted.
fitLim : float,optional
- Maximum flux where the level of absorption will trigger
+ Maximum flux where the level of absorption will trigger
identification of the region as an absorption complex. Default = .98.
(ex: for a minSize=.98, a region where all the flux is between 1.0 and
.99 will not be separated out to be fit as an absorbing complex, but
a region that contains a point where the flux is .97 will be fit
as an absorbing complex.)
minLength : int, optional
- number of cells required for a complex to be included.
+ number of cells required for a complex to be included.
default is 3 cells.
maxLength : int, optional
number of cells required for a complex to be split up. Default
is 1000 cells.
splitLim : float, optional
if attempting to split a region for being larger than maxlength
- the point of the split must have a flux greater than splitLim
+ the point of the split must have a flux greater than splitLim
(ie: absorption greater than splitLim). Default= .99.
Returns
-------
- cBounds : (3,)
- list of bounds in the form [[i0,i1,i2],...] where i0 is the
+ cBounds : (3,)
+ list of bounds in the form [[i0,i1,i2],...] where i0 is the
index of the maximum flux for a complex, i1 is the index of the
- beginning of the complex, and i2 is the index of the end of the
+ beginning of the complex, and i2 is the index of the end of the
complex. Indexes refer to the indices of x and yDat.
"""
@@ -743,7 +743,7 @@
x : (N) ndarray
Array of wavelength
linesP: (3,) ndarray
- Array giving sets of line parameters in
+ Array giving sets of line parameters in
form [[N1, b1, z1], ...]
speciesDict : dictionary
Dictionary containing all relevant parameters needed
@@ -762,7 +762,7 @@
g=speciesDict['Gamma'][i]
wl=speciesDict['wavelength'][i]
y = y+ _gen_tau(x,p,f,g,wl)
- if firstLine:
+ if firstLine:
break
flux = np.exp(-y)
@@ -772,15 +772,15 @@
"""This calculates a flux distribution for given parameters using the yt
voigt profile generator"""
N,b,z= p
-
+
#Calculating quantities
tau_o = 1.4973614E-15*N*f*lambda_unshifted/b
a=7.95774715459E-15*Gamma*lambda_unshifted/b
x=299792.458/b*(lambda_unshifted*(1+z)/t-1)
-
+
H = np.zeros(len(x))
H = voigt(a,x)
-
+
tau = tau_o*H
return tau
@@ -797,8 +797,8 @@
Parameters
----------
- pTotal : (3,) ndarray
- Array with form [[N1, b1, z1], ...]
+ pTotal : (3,) ndarray
+ Array with form [[N1, b1, z1], ...]
x : (N) ndarray
array of wavelengths [nm]
yDat : (N) ndarray
@@ -816,7 +816,7 @@
error : (N) ndarray
the difference between the fit generated by the parameters
given in pTotal multiplied by the previous fit and the desired
- flux profile, w/ first index modified appropriately for bad
+ flux profile, w/ first index modified appropriately for bad
parameter choices and additional penalty for fitting with a lower
flux than observed.
"""
@@ -834,13 +834,13 @@
def _check_params(p, speciesDict,xb):
"""
- Check to see if any of the parameters in p fall outside the range
+ Check to see if any of the parameters in p fall outside the range
given in speciesDict or on the boundaries
Parameters
----------
p : (3,) ndarray
- array with form [[N1, b1, z1], ...]
+ array with form [[N1, b1, z1], ...]
speciesDict : dictionary
dictionary with properties giving the max and min
values appropriate for each parameter N,b, and z.
@@ -865,20 +865,20 @@
any(p[:,2] >= maxz) or\
any(p[:,2] <= minz):
check = 999
-
+
return check
def _check_optimization_init(p,speciesDict,initz,xb,yDat,yFit,minSize,errorBound):
"""
Check to see if any of the parameters in p are the
- same as initial paramters and if so, attempt to
+ same as initial paramters and if so, attempt to
split the region and refit it.
Parameters
----------
p : (3,) ndarray
- array with form [[N1, b1, z1], ...]
+ array with form [[N1, b1, z1], ...]
speciesDict : dictionary
dictionary with properties giving the max and min
values appropriate for each parameter N,b, and z.
@@ -902,7 +902,7 @@
# Attempt to split region
split = _split_region(yDat,init_bounds,newSplitLim)
-
+
# If we can't split it, just reject it. Its unphysical
# to just keep the default parameters and we're out of
# options at this point
@@ -939,7 +939,7 @@
Parameters
----------
p : (3,) ndarray
- array with form [[N1, b1, z1], ...]
+ array with form [[N1, b1, z1], ...]
speciesDict : dictionary
dictionary with properties giving the max and min
values appropriate for each parameter N,b, and z.
@@ -979,7 +979,7 @@
# we just set it to no fit because we've tried
# everything else at this point. this region just sucks :(
remove_lines.append(i)
-
+
if remove_lines:
p = np.delete(p, remove_lines, axis=0)
@@ -988,7 +988,7 @@
def _output_fit(lineDic, file_name = 'spectrum_fit.h5'):
"""
This function is designed to output the parameters of the series
- of lines used to fit an absorption spectrum.
+ of lines used to fit an absorption spectrum.
The dataset contains entries in the form species/N, species/b
species/z, and species/complex. The ith entry in each of the datasets
@@ -999,11 +999,11 @@
Parameters
----------
lineDic : dictionary
- Dictionary of dictionaries representing the fit lines.
+ Dictionary of dictionaries representing the fit lines.
Top level keys are the species given in orderFits and the corresponding
- entries are dictionaries with the keys 'N','b','z', and 'group#'.
+ entries are dictionaries with the keys 'N','b','z', and 'group#'.
Each of these corresponds to a list of the parameters for every
- accepted fitted line.
+ accepted fitted line.
fileName : string, optional
Name of the file to output fit to. Default = 'spectrum_fit.h5'
https://bitbucket.org/yt_analysis/yt/commits/0e09458b79c9/
Changeset: 0e09458b79c9
Branch: yt
User: migueldvb
Date: 2016-03-14 20:25:11+00:00
Summary: Fix indentation
Affected #: 1 file
diff -r 5595d7cdcb461ff2f15a4b931c7065e6d11d9fbd -r 0e09458b79c963950d3c2ff5e68dfd9d32035382 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
@@ -864,7 +864,7 @@
any(p[:,1] <= speciesDict['minb']) or\
any(p[:,2] >= maxz) or\
any(p[:,2] <= minz):
- check = 999
+ check = 999
return check
@@ -894,38 +894,38 @@
any(p[:,1] == speciesDict['init_b']) or\
any(p[:,1] == speciesDict['maxb']):
- # These are the initial bounds
- init_bounds = [yDat.argmin(),0,len(xb)-1]
+ # These are the initial bounds
+ init_bounds = [yDat.argmin(),0,len(xb)-1]
- # Gratitutous limit for splitting region
- newSplitLim = 1 - (1-min(yDat))*.5
+ # Gratitutous limit for splitting region
+ newSplitLim = 1 - (1-min(yDat))*.5
- # Attempt to split region
- split = _split_region(yDat,init_bounds,newSplitLim)
+ # Attempt to split region
+ split = _split_region(yDat,init_bounds,newSplitLim)
- # If we can't split it, just reject it. Its unphysical
- # to just keep the default parameters and we're out of
- # options at this point
- if not split:
- return []
+ # If we can't split it, just reject it. Its unphysical
+ # to just keep the default parameters and we're out of
+ # options at this point
+ if not split:
+ return []
- # Else set up the bounds for each region and fit separately
- b1,b2 = split[0][2], split[1][1]
+ # Else set up the bounds for each region and fit separately
+ b1,b2 = split[0][2], split[1][1]
- p1,flag = _complex_fit(xb[:b1], yDat[:b1], yFit[:b1],
- initz, minSize, errorBound, speciesDict)
+ p1,flag = _complex_fit(xb[:b1], yDat[:b1], yFit[:b1],
+ initz, minSize, errorBound, speciesDict)
- p2,flag = _complex_fit(xb[b2:], yDat[b2:], yFit[b2:],
- initz, minSize, errorBound, speciesDict)
+ p2,flag = _complex_fit(xb[b2:], yDat[b2:], yFit[b2:],
+ initz, minSize, errorBound, speciesDict)
- # Make the final line parameters. Its annoying because
- # one or both regions may have fit to nothing
- if np.size(p1)> 0 and np.size(p2)>0:
- p = np.r_[p1,p2]
- elif np.size(p1) > 0:
- p = p1
- else:
- p = p2
+ # Make the final line parameters. Its annoying because
+ # one or both regions may have fit to nothing
+ if np.size(p1)> 0 and np.size(p2)>0:
+ p = np.r_[p1,p2]
+ elif np.size(p1) > 0:
+ p = p1
+ else:
+ p = p2
return p
https://bitbucket.org/yt_analysis/yt/commits/e14f3ad3b7f2/
Changeset: e14f3ad3b7f2
Branch: yt
User: migueldvb
Date: 2016-04-07 02:00:08+00:00
Summary: Write logging functions as format string + args
This should have a little better performance than string formatting
Affected #: 1 file
diff -r 0e09458b79c963950d3c2ff5e68dfd9d32035382 -r e14f3ad3b7f2e27dac2c2141f7e76b0922fd8a5f 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
@@ -54,7 +54,7 @@
# lambda, flux, and tau are wavelength, flux, and optical depth
self.lambda_min = lambda_min
self.lambda_max = lambda_max
- self.lambda_field = YTArray(np.linspace(lambda_min, lambda_max,
+ self.lambda_field = YTArray(np.linspace(lambda_min, lambda_max,
n_lambda), "angstrom")
self.tau_field = None
self.flux_field = None
@@ -237,7 +237,7 @@
def _apply_observing_redshift(self, field_data, use_peculiar_velocity,
observing_redshift):
"""
- Change the redshifts of individual absorbers to account for the
+ Change the redshifts of individual absorbers to account for the
redshift at which the observer sits.
The intermediate redshift that is seen by an observer
@@ -275,9 +275,9 @@
"""
Add continuum features to the spectrum.
"""
- # Change the redshifts of continuum sources to account for the
+ # Change the redshifts of continuum sources to account for the
# redshift at which the observer sits
- redshift, redshift_eff = self._apply_observing_redshift(field_data,
+ redshift, redshift_eff = self._apply_observing_redshift(field_data,
use_peculiar_velocity, observing_redshift)
# Only add continuum features down to tau of 1.e-4.
@@ -319,9 +319,9 @@
Add the absorption lines to the spectrum.
"""
- # Change the redshifts of individual absorbers to account for the
+ # Change the redshifts of individual absorbers to account for the
# redshift at which the observer sits
- redshift, redshift_eff = self._apply_observing_redshift(field_data,
+ redshift, redshift_eff = self._apply_observing_redshift(field_data,
use_peculiar_velocity, observing_redshift)
# Widen wavelength window until optical depth falls below this tau
@@ -349,13 +349,13 @@
# we want to know the bin index in the lambda_field array
# where each line has its central wavelength after being
# redshifted. however, because we don't know a priori how wide
- # a line will be (ie DLAs), we have to include bin indices
- # *outside* the spectral range of the AbsorptionSpectrum
+ # a line will be (ie DLAs), we have to include bin indices
+ # *outside* the spectral range of the AbsorptionSpectrum
# object. Thus, we find the "equivalent" bin index, which
# may be <0 or >the size of the array. In the end, we deposit
# the bins that actually overlap with the AbsorptionSpectrum's
# range in lambda.
-
+
# this equation gives us the "equivalent" bin index for each line
# if it were placed into the self.lambda_field array
center_index = (lambda_obs.in_units('Angstrom').d - self.lambda_min) \
@@ -399,10 +399,10 @@
# a note to the user about which lines components are unresolved
if (thermal_width < self.bin_width).any():
- mylog.info(("%d out of %d line components will be " + \
- "deposited as unresolved lines.") %
- ((thermal_width < self.bin_width).sum(),
- n_absorbers))
+ mylog.info("%d out of %d line components will be " +
+ "deposited as unresolved lines.",
+ (thermal_width < self.bin_width).sum(),
+ n_absorbers)
# provide a progress bar with information about lines processsed
pbar = get_pbar("Adding line - %s [%f A]: " % \
@@ -411,13 +411,13 @@
# for a given transition, step through each location in the
# observed spectrum where it occurs and deposit a voigt profile
for i in parallel_objects(np.arange(n_absorbers), njobs=-1):
-
- # the virtual window into which the line is deposited initially
- # spans a region of 2 coarse spectral bins
+
+ # the virtual window into which the line is deposited initially
+ # spans a region of 2 coarse spectral bins
# (one on each side of the center_index) but the window
# can expand as necessary.
# it will continue to expand until the tau value in the far
- # edge of the wings is less than the min_tau value or it
+ # edge of the wings is less than the min_tau value or it
# reaches the edge of the spectrum
window_width_in_bins = 2
@@ -425,18 +425,18 @@
left_index = (center_index[i] - window_width_in_bins/2)
right_index = (center_index[i] + window_width_in_bins/2)
n_vbins = (right_index - left_index) * n_vbins_per_bin[i]
-
+
# the array of virtual bins in lambda space
vbins = \
- np.linspace(self.lambda_min + self.bin_width.d * left_index,
- self.lambda_min + self.bin_width.d * right_index,
+ np.linspace(self.lambda_min + self.bin_width.d * left_index,
+ self.lambda_min + self.bin_width.d * right_index,
n_vbins, endpoint=False)
# the virtual bins and their corresponding opacities
vbins, vtau = \
tau_profile(
- lambda_0, line['f_value'], line['gamma'],
- thermb[i], cdens[i],
+ lambda_0, line['f_value'], line['gamma'],
+ thermb[i], cdens[i],
delta_lambda=dlambda[i], lambda_bins=vbins)
# If tau has not dropped below min tau threshold by the
@@ -468,7 +468,7 @@
continue
# otherwise, determine how much of the original spectrum
- # is intersected by the expanded line window to be deposited,
+ # is intersected by the expanded line window to be deposited,
# and deposit the Equivalent Width data into that intersecting
# window in the original spectrum's tau
else:
@@ -517,7 +517,7 @@
"""
if filename is None:
return
- mylog.info("Writing absorber list: %s." % filename)
+ mylog.info("Writing absorber list: %s.", filename)
self.absorbers_list.sort(key=lambda obj: obj['wavelength'])
f = open(filename, 'w')
f.write('#%-14s %-14s %-12s %-14s %-15s %-9s %-10s\n' %
@@ -534,7 +534,7 @@
"""
Write spectrum to an ascii file.
"""
- mylog.info("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_field.size):
@@ -547,7 +547,7 @@
"""
Write spectrum to a fits file.
"""
- mylog.info("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_field)
col2 = pyfits.Column(name='flux', format='E', array=self.flux_field)
cols = pyfits.ColDefs([col1, col2])
@@ -560,7 +560,7 @@
Write spectrum to an hdf5 file.
"""
- mylog.info("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_field)
output.create_dataset('tau', data=self.tau_field)
https://bitbucket.org/yt_analysis/yt/commits/5d118c06fa7e/
Changeset: 5d118c06fa7e
Branch: yt
User: ngoldbaum
Date: 2016-04-20 18:20:28+00:00
Summary: Merged in migueldvb/yt (pull request #2114)
absorption spectrum module minor changes
Affected #: 2 files
diff -r 8f91fd7329595d87b927cec203b1246b1de963f9 -r 5d118c06fa7eecc494eff5d231ebf3c47d9d4ab9 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
@@ -54,7 +54,7 @@
# lambda, flux, and tau are wavelength, flux, and optical depth
self.lambda_min = lambda_min
self.lambda_max = lambda_max
- self.lambda_field = YTArray(np.linspace(lambda_min, lambda_max,
+ self.lambda_field = YTArray(np.linspace(lambda_min, lambda_max,
n_lambda), "angstrom")
self.tau_field = None
self.flux_field = None
@@ -237,7 +237,7 @@
def _apply_observing_redshift(self, field_data, use_peculiar_velocity,
observing_redshift):
"""
- Change the redshifts of individual absorbers to account for the
+ Change the redshifts of individual absorbers to account for the
redshift at which the observer sits.
The intermediate redshift that is seen by an observer
@@ -275,9 +275,9 @@
"""
Add continuum features to the spectrum.
"""
- # Change the redshifts of continuum sources to account for the
+ # Change the redshifts of continuum sources to account for the
# redshift at which the observer sits
- redshift, redshift_eff = self._apply_observing_redshift(field_data,
+ redshift, redshift_eff = self._apply_observing_redshift(field_data,
use_peculiar_velocity, observing_redshift)
# Only add continuum features down to tau of 1.e-4.
@@ -319,9 +319,9 @@
Add the absorption lines to the spectrum.
"""
- # Change the redshifts of individual absorbers to account for the
+ # Change the redshifts of individual absorbers to account for the
# redshift at which the observer sits
- redshift, redshift_eff = self._apply_observing_redshift(field_data,
+ redshift, redshift_eff = self._apply_observing_redshift(field_data,
use_peculiar_velocity, observing_redshift)
# Widen wavelength window until optical depth falls below this tau
@@ -349,13 +349,13 @@
# we want to know the bin index in the lambda_field array
# where each line has its central wavelength after being
# redshifted. however, because we don't know a priori how wide
- # a line will be (ie DLAs), we have to include bin indices
- # *outside* the spectral range of the AbsorptionSpectrum
+ # a line will be (ie DLAs), we have to include bin indices
+ # *outside* the spectral range of the AbsorptionSpectrum
# object. Thus, we find the "equivalent" bin index, which
# may be <0 or >the size of the array. In the end, we deposit
# the bins that actually overlap with the AbsorptionSpectrum's
# range in lambda.
-
+
# this equation gives us the "equivalent" bin index for each line
# if it were placed into the self.lambda_field array
center_index = (lambda_obs.in_units('Angstrom').d - self.lambda_min) \
@@ -399,10 +399,10 @@
# a note to the user about which lines components are unresolved
if (thermal_width < self.bin_width).any():
- mylog.info(("%d out of %d line components will be " + \
- "deposited as unresolved lines.") %
- ((thermal_width < self.bin_width).sum(),
- n_absorbers))
+ mylog.info("%d out of %d line components will be " +
+ "deposited as unresolved lines.",
+ (thermal_width < self.bin_width).sum(),
+ n_absorbers)
# provide a progress bar with information about lines processsed
pbar = get_pbar("Adding line - %s [%f A]: " % \
@@ -411,13 +411,13 @@
# for a given transition, step through each location in the
# observed spectrum where it occurs and deposit a voigt profile
for i in parallel_objects(np.arange(n_absorbers), njobs=-1):
-
- # the virtual window into which the line is deposited initially
- # spans a region of 2 coarse spectral bins
+
+ # the virtual window into which the line is deposited initially
+ # spans a region of 2 coarse spectral bins
# (one on each side of the center_index) but the window
# can expand as necessary.
# it will continue to expand until the tau value in the far
- # edge of the wings is less than the min_tau value or it
+ # edge of the wings is less than the min_tau value or it
# reaches the edge of the spectrum
window_width_in_bins = 2
@@ -425,18 +425,18 @@
left_index = (center_index[i] - window_width_in_bins/2)
right_index = (center_index[i] + window_width_in_bins/2)
n_vbins = (right_index - left_index) * n_vbins_per_bin[i]
-
+
# the array of virtual bins in lambda space
vbins = \
- np.linspace(self.lambda_min + self.bin_width.d * left_index,
- self.lambda_min + self.bin_width.d * right_index,
+ np.linspace(self.lambda_min + self.bin_width.d * left_index,
+ self.lambda_min + self.bin_width.d * right_index,
n_vbins, endpoint=False)
# the virtual bins and their corresponding opacities
vbins, vtau = \
tau_profile(
- lambda_0, line['f_value'], line['gamma'],
- thermb[i], cdens[i],
+ lambda_0, line['f_value'], line['gamma'],
+ thermb[i], cdens[i],
delta_lambda=dlambda[i], lambda_bins=vbins)
# If tau has not dropped below min tau threshold by the
@@ -468,7 +468,7 @@
continue
# otherwise, determine how much of the original spectrum
- # is intersected by the expanded line window to be deposited,
+ # is intersected by the expanded line window to be deposited,
# and deposit the Equivalent Width data into that intersecting
# window in the original spectrum's tau
else:
@@ -517,7 +517,7 @@
"""
if filename is None:
return
- mylog.info("Writing absorber list: %s." % filename)
+ mylog.info("Writing absorber list: %s.", filename)
self.absorbers_list.sort(key=lambda obj: obj['wavelength'])
f = open(filename, 'w')
f.write('#%-14s %-14s %-12s %-14s %-15s %-9s %-10s\n' %
@@ -534,7 +534,7 @@
"""
Write spectrum to an ascii file.
"""
- mylog.info("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_field.size):
@@ -547,7 +547,7 @@
"""
Write spectrum to a fits file.
"""
- mylog.info("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_field)
col2 = pyfits.Column(name='flux', format='E', array=self.flux_field)
cols = pyfits.ColDefs([col1, col2])
@@ -560,7 +560,7 @@
Write spectrum to an hdf5 file.
"""
- mylog.info("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_field)
output.create_dataset('tau', data=self.tau_field)
diff -r 8f91fd7329595d87b927cec203b1246b1de963f9 -r 5d118c06fa7eecc494eff5d231ebf3c47d9d4ab9 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
@@ -12,14 +12,14 @@
optimize = _scipy.optimize
-def generate_total_fit(x, fluxData, orderFits, speciesDicts,
+def generate_total_fit(x, fluxData, orderFits, speciesDicts,
minError=1E-4, complexLim=.995,
- fitLim=.97, minLength=3,
+ fitLim=.97, minLength=3,
maxLength=1000, splitLim=.99,
output_file=None):
"""
- This function is designed to fit an absorption spectrum by breaking
+ This function is designed to fit an absorption spectrum by breaking
the spectrum up into absorption complexes, and iteratively adding
and optimizing voigt profiles to each complex.
@@ -31,46 +31,46 @@
array of flux corresponding to the wavelengths given
in x. (needs to be the same size as x)
orderFits : list
- list of the names of the species in the order that they
+ list of the names of the species in the order that they
should be fit. Names should correspond to the names of the species
given in speciesDicts. (ex: ['lya','OVI'])
speciesDicts : dictionary
Dictionary of dictionaries (I'm addicted to dictionaries, I
confess). Top level keys should be the names of all the species given
- in orderFits. The entries should be dictionaries containing all
- relevant parameters needed to create an absorption line of a given
+ in orderFits. The entries should be dictionaries containing all
+ relevant parameters needed to create an absorption line of a given
species (f,Gamma,lambda0) as well as max and min values for parameters
to be fit
complexLim : float, optional
- Maximum flux to start the edge of an absorption complex. Different
- from fitLim because it decides extent of a complex rather than
- whether or not a complex is accepted.
+ Maximum flux to start the edge of an absorption complex. Different
+ from fitLim because it decides extent of a complex rather than
+ whether or not a complex is accepted.
fitLim : float,optional
- Maximum flux where the level of absorption will trigger
+ Maximum flux where the level of absorption will trigger
identification of the region as an absorption complex. Default = .98.
(ex: for a minSize=.98, a region where all the flux is between 1.0 and
.99 will not be separated out to be fit as an absorbing complex, but
a region that contains a point where the flux is .97 will be fit
as an absorbing complex.)
minLength : int, optional
- number of cells required for a complex to be included.
+ number of cells required for a complex to be included.
default is 3 cells.
maxLength : int, optional
number of cells required for a complex to be split up. Default
is 1000 cells.
splitLim : float, optional
if attempting to split a region for being larger than maxlength
- the point of the split must have a flux greater than splitLim
+ the point of the split must have a flux greater than splitLim
(ie: absorption greater than splitLim). Default= .99.
output_file : string, optional
- location to save the results of the fit.
+ location to save the results of the fit.
Returns
-------
allSpeciesLines : dictionary
- Dictionary of dictionaries representing the fit lines.
+ Dictionary of dictionaries representing the fit lines.
Top level keys are the species given in orderFits and the corresponding
- entries are dictionaries with the keys 'N','b','z', and 'group#'.
+ entries are dictionaries with the keys 'N','b','z', and 'group#'.
Each of these corresponds to a list of the parameters for every
accepted fitted line. (ie: N[0],b[0],z[0] will create a line that
fits some part of the absorption spectrum). 'group#' is a similar list
@@ -87,7 +87,7 @@
# convert to NumPy array if we have a YTArray
if isinstance(x, YTArray):
x = x.d
-
+
#Empty dictionary for fitted lines
allSpeciesLines = {}
@@ -129,8 +129,8 @@
#Check if any flux at partner sites
if not _line_exists(speciesDict['wavelength'],
- fluxData,z,x0,xRes,fitLim):
- continue
+ fluxData,z,x0,xRes,fitLim):
+ continue
#Fit Using complex tools
newLinesP,flag=_complex_fit(xBounded,yDatBounded,yFitBounded,
@@ -139,7 +139,7 @@
#If flagged as a bad fit, species is lyman alpha,
# and it may be a saturated line, use special tools
if flag and species=='lya' and min(yDatBounded)<.1:
- newLinesP=_large_flag_fit(xBounded,yDatBounded,
+ newLinesP=_large_flag_fit(xBounded,yDatBounded,
yFitBounded,z,speciesDict,
fitLim,minError)
@@ -151,7 +151,7 @@
#Check existence of partner lines if applicable
if len(speciesDict['wavelength']) != 1:
- newLinesP = _remove_unaccepted_partners(newLinesP, x, fluxData,
+ newLinesP = _remove_unaccepted_partners(newLinesP, x, fluxData,
b, minError, x0, xRes, speciesDict)
@@ -177,15 +177,15 @@
return (allSpeciesLines,yFit)
-def _complex_fit(x, yDat, yFit, initz, minSize, errBound, speciesDict,
+def _complex_fit(x, yDat, yFit, initz, minSize, errBound, speciesDict,
initP=None):
""" Fit an absorption complex by iteratively adding and optimizing
voigt profiles.
-
+
A complex is defined as a region where some number of lines may be present,
or a region of non zero of absorption. Lines are iteratively added
and optimized until the difference between the flux generated using
- the optimized parameters has a least squares difference between the
+ the optimized parameters has a least squares difference between the
desired flux profile less than the error bound.
Parameters
@@ -199,7 +199,7 @@
array of flux profile fitted for the wavelength
space given by x already. Same size as x.
initz : float
- redshift to try putting first line at
+ redshift to try putting first line at
(maximum absorption for region)
minsize : float
minimum absorption allowed for a line to still count as a line
@@ -219,7 +219,7 @@
Returns
-------
linesP : (3,) ndarray
- Array of best parameters if a good enough fit is found in
+ Array of best parameters if a good enough fit is found in
the form [[N1,b1,z1], [N2,b2,z2],...]
flag : bool
boolean value indicating the success of the fit (True if unsuccessful)
@@ -227,8 +227,8 @@
#Setup initial line guesses
if initP is None: #Regular fit
- initP = [0,0,0]
- if min(yDat)<.01: #Large lines get larger initial guess
+ initP = [0,0,0]
+ if min(yDat)<.01: #Large lines get larger initial guess
initP[0] = speciesDict['init_N']*10**2
elif min(yDat)<.5:
initP[0] = speciesDict['init_N']*10**1
@@ -246,9 +246,9 @@
wl0 = speciesDict['wavelength'][0]
#Check if first line exists still
- if min(yDat-yFit+1)>minSize:
+ if min(yDat-yFit+1)>minSize:
return [],False
-
+
#Values to proceed through first run
errSq,prevErrSq,prevLinesP=1,10*len(x),[]
@@ -283,9 +283,9 @@
if any(linesP[:,1]==speciesDict['init_b']):
flag = True
break
-
+
#If good enough, break
- if errSq < errBound:
+ if errSq < errBound:
break
#If last fit was worse, reject the last line and revert to last fit
@@ -297,7 +297,7 @@
linesP = prevLinesP
break
- #If too many lines
+ #If too many lines
if np.shape(linesP)[0]>8 or np.size(linesP)+3>=len(x):
#If its fitable by flag tools and still bad, use flag tools
if errSq >1E2*errBound and speciesDict['name']=='HI lya':
@@ -311,11 +311,11 @@
prevLinesP = linesP
#Set up initial condition for new line
- newP = [0,0,0]
+ newP = [0,0,0]
yAdjusted = 1+yFit*yNewFit-yDat
-
- if min(yAdjusted)<.01: #Large lines get larger initial guess
+
+ if min(yAdjusted)<.01: #Large lines get larger initial guess
newP[0] = speciesDict['init_N']*10**2
elif min(yAdjusted)<.5:
newP[0] = speciesDict['init_N']*10**1
@@ -333,7 +333,7 @@
remove=[]
for i,p in enumerate(linesP):
check=_check_params(np.array([p]),speciesDict,x)
- if check:
+ if check:
remove.append(i)
linesP = np.delete(linesP,remove,axis=0)
@@ -358,7 +358,7 @@
array of flux profile fitted for the wavelength
space given by x already. Same size as x.
initz : float
- redshift to try putting first line at
+ redshift to try putting first line at
(maximum absorption for region)
speciesDict : dictionary
dictionary containing all relevant parameters needed
@@ -374,15 +374,15 @@
Returns
-------
bestP : (3,) ndarray
- array of best parameters if a good enough fit is found in
- the form [[N1,b1,z1], [N2,b2,z2],...]
+ array of best parameters if a good enough fit is found in
+ the form [[N1,b1,z1], [N2,b2,z2],...]
"""
#Set up some initial line guesses
lineTests = _get_test_lines(initz)
#Keep track of the lowest achieved error
- bestError = 1000
+ bestError = 1000
#Iterate through test line guesses
for initLines in lineTests:
@@ -403,7 +403,7 @@
bestError = errSq
bestP = linesP
- if bestError>10*errBound*len(x):
+ if bestError>10*errBound*len(x):
return []
else:
return bestP
@@ -411,11 +411,11 @@
def _get_test_lines(initz):
"""
Returns a 3d numpy array of lines to test as initial guesses for difficult
- to fit lyman alpha absorbers that are saturated.
-
+ to fit lyman alpha absorbers that are saturated.
+
The array is 3d because
the first dimension gives separate initial guesses, the second dimension
- has multiple lines for the same guess (trying a broad line plus a
+ has multiple lines for the same guess (trying a broad line plus a
saturated line) and the 3d dimension contains the 3 fit parameters (N,b,z)
Parameters
@@ -426,7 +426,7 @@
Returns
-------
testP : (,3,) ndarray
- numpy array of the form
+ numpy array of the form
[[[N1a,b1a,z1a], [N1b,b1b,z1b]], [[N2a,b2,z2a],...] ...]
"""
@@ -454,8 +454,8 @@
return testP
def _get_bounds(z, b, wl, x0, xRes):
- """
- Gets the indices of range of wavelength that the wavelength wl is in
+ """
+ Gets the indices of range of wavelength that the wavelength wl is in
with the size of some initial wavelength range.
Used for checking if species with multiple lines (as in the OVI doublet)
@@ -466,21 +466,21 @@
z : float
redshift
b : (3) ndarray/list
- initial bounds in form [i0,i1,i2] where i0 is the index of the
- minimum flux for the complex, i1 is index of the lower wavelength
+ initial bounds in form [i0,i1,i2] where i0 is the index of the
+ minimum flux for the complex, i1 is index of the lower wavelength
edge of the complex, and i2 is the index of the higher wavelength
edge of the complex.
wl : float
- unredshifted wavelength of the peak of the new region
+ unredshifted wavelength of the peak of the new region
x0 : float
wavelength of the index 0
xRes : float
difference in wavelength for two consecutive indices
-
+
Returns
-------
indices : (2) tuple
- Tuple (i1,i2) where i1 is the index of the lower wavelength bound of
+ Tuple (i1,i2) where i1 is the index of the lower wavelength bound of
the new region and i2 is the index of the higher wavelength bound of
the new region
"""
@@ -492,7 +492,7 @@
return indices
-def _remove_unaccepted_partners(linesP, x, y, b, errBound,
+def _remove_unaccepted_partners(linesP, x, y, b, errBound,
x0, xRes, speciesDict):
"""
Given a set of parameters [N,b,z] that form multiple lines for a given
@@ -507,18 +507,18 @@
Parameters
----------
linesP : (3,) ndarray
- array giving sets of line parameters in
+ array giving sets of line parameters in
form [[N1, b1, z1], ...]
x : (N) ndarray
wavelength array [nm]
y : (N) ndarray
normalized flux array of original data
b : (3) tuple/list/ndarray
- indices that give the bounds of the original region so that another
+ indices that give the bounds of the original region so that another
region of similar size can be used to determine the goodness
of fit of the other wavelengths
errBound : float
- size of the error that is appropriate for a given region,
+ size of the error that is appropriate for a given region,
adjusted to account for the size of the region.
Returns
@@ -569,7 +569,7 @@
#Remove all bad line fits
linesP = np.delete(linesP,removeLines,axis=0)
- return linesP
+ return linesP
@@ -588,12 +588,12 @@
xRes : float
difference in wavelength between consecutive cells in flux array
fluxMin : float
- maximum flux to count as a line existing.
+ maximum flux to count as a line existing.
Returns
-------
- flag : boolean
+ flag : boolean
value indicating whether all lines exist. True if all lines exist
"""
@@ -618,7 +618,7 @@
def _find_complexes(x, yDat, complexLim=.999, fitLim=.99,
minLength =3, maxLength=1000, splitLim=.99):
"""Breaks up the wavelength space into groups
- where there is some absorption.
+ where there is some absorption.
Parameters
----------
@@ -628,33 +628,33 @@
array of flux corresponding to the wavelengths given
in x. (needs to be the same size as x)
complexLim : float, optional
- Maximum flux to start the edge of an absorption complex. Different
- from fitLim because it decides extent of a complex rather than
- whether or not a complex is accepted.
+ Maximum flux to start the edge of an absorption complex. Different
+ from fitLim because it decides extent of a complex rather than
+ whether or not a complex is accepted.
fitLim : float,optional
- Maximum flux where the level of absorption will trigger
+ Maximum flux where the level of absorption will trigger
identification of the region as an absorption complex. Default = .98.
(ex: for a minSize=.98, a region where all the flux is between 1.0 and
.99 will not be separated out to be fit as an absorbing complex, but
a region that contains a point where the flux is .97 will be fit
as an absorbing complex.)
minLength : int, optional
- number of cells required for a complex to be included.
+ number of cells required for a complex to be included.
default is 3 cells.
maxLength : int, optional
number of cells required for a complex to be split up. Default
is 1000 cells.
splitLim : float, optional
if attempting to split a region for being larger than maxlength
- the point of the split must have a flux greater than splitLim
+ the point of the split must have a flux greater than splitLim
(ie: absorption greater than splitLim). Default= .99.
Returns
-------
- cBounds : (3,)
- list of bounds in the form [[i0,i1,i2],...] where i0 is the
+ cBounds : (3,)
+ list of bounds in the form [[i0,i1,i2],...] where i0 is the
index of the maximum flux for a complex, i1 is the index of the
- beginning of the complex, and i2 is the index of the end of the
+ beginning of the complex, and i2 is the index of the end of the
complex. Indexes refer to the indices of x and yDat.
"""
@@ -743,7 +743,7 @@
x : (N) ndarray
Array of wavelength
linesP: (3,) ndarray
- Array giving sets of line parameters in
+ Array giving sets of line parameters in
form [[N1, b1, z1], ...]
speciesDict : dictionary
Dictionary containing all relevant parameters needed
@@ -762,7 +762,7 @@
g=speciesDict['Gamma'][i]
wl=speciesDict['wavelength'][i]
y = y+ _gen_tau(x,p,f,g,wl)
- if firstLine:
+ if firstLine:
break
flux = np.exp(-y)
@@ -772,15 +772,15 @@
"""This calculates a flux distribution for given parameters using the yt
voigt profile generator"""
N,b,z= p
-
+
#Calculating quantities
tau_o = 1.4973614E-15*N*f*lambda_unshifted/b
a=7.95774715459E-15*Gamma*lambda_unshifted/b
x=299792.458/b*(lambda_unshifted*(1+z)/t-1)
-
+
H = np.zeros(len(x))
H = voigt(a,x)
-
+
tau = tau_o*H
return tau
@@ -797,8 +797,8 @@
Parameters
----------
- pTotal : (3,) ndarray
- Array with form [[N1, b1, z1], ...]
+ pTotal : (3,) ndarray
+ Array with form [[N1, b1, z1], ...]
x : (N) ndarray
array of wavelengths [nm]
yDat : (N) ndarray
@@ -816,7 +816,7 @@
error : (N) ndarray
the difference between the fit generated by the parameters
given in pTotal multiplied by the previous fit and the desired
- flux profile, w/ first index modified appropriately for bad
+ flux profile, w/ first index modified appropriately for bad
parameter choices and additional penalty for fitting with a lower
flux than observed.
"""
@@ -834,13 +834,13 @@
def _check_params(p, speciesDict,xb):
"""
- Check to see if any of the parameters in p fall outside the range
+ Check to see if any of the parameters in p fall outside the range
given in speciesDict or on the boundaries
Parameters
----------
p : (3,) ndarray
- array with form [[N1, b1, z1], ...]
+ array with form [[N1, b1, z1], ...]
speciesDict : dictionary
dictionary with properties giving the max and min
values appropriate for each parameter N,b, and z.
@@ -864,21 +864,21 @@
any(p[:,1] <= speciesDict['minb']) or\
any(p[:,2] >= maxz) or\
any(p[:,2] <= minz):
- check = 999
-
+ check = 999
+
return check
def _check_optimization_init(p,speciesDict,initz,xb,yDat,yFit,minSize,errorBound):
"""
Check to see if any of the parameters in p are the
- same as initial paramters and if so, attempt to
+ same as initial paramters and if so, attempt to
split the region and refit it.
Parameters
----------
p : (3,) ndarray
- array with form [[N1, b1, z1], ...]
+ array with form [[N1, b1, z1], ...]
speciesDict : dictionary
dictionary with properties giving the max and min
values appropriate for each parameter N,b, and z.
@@ -894,38 +894,38 @@
any(p[:,1] == speciesDict['init_b']) or\
any(p[:,1] == speciesDict['maxb']):
- # These are the initial bounds
- init_bounds = [yDat.argmin(),0,len(xb)-1]
+ # These are the initial bounds
+ init_bounds = [yDat.argmin(),0,len(xb)-1]
- # Gratitutous limit for splitting region
- newSplitLim = 1 - (1-min(yDat))*.5
+ # Gratitutous limit for splitting region
+ newSplitLim = 1 - (1-min(yDat))*.5
- # Attempt to split region
- split = _split_region(yDat,init_bounds,newSplitLim)
-
- # If we can't split it, just reject it. Its unphysical
- # to just keep the default parameters and we're out of
- # options at this point
- if not split:
- return []
+ # Attempt to split region
+ split = _split_region(yDat,init_bounds,newSplitLim)
- # Else set up the bounds for each region and fit separately
- b1,b2 = split[0][2], split[1][1]
+ # If we can't split it, just reject it. Its unphysical
+ # to just keep the default parameters and we're out of
+ # options at this point
+ if not split:
+ return []
- p1,flag = _complex_fit(xb[:b1], yDat[:b1], yFit[:b1],
- initz, minSize, errorBound, speciesDict)
+ # Else set up the bounds for each region and fit separately
+ b1,b2 = split[0][2], split[1][1]
- p2,flag = _complex_fit(xb[b2:], yDat[b2:], yFit[b2:],
- initz, minSize, errorBound, speciesDict)
+ p1,flag = _complex_fit(xb[:b1], yDat[:b1], yFit[:b1],
+ initz, minSize, errorBound, speciesDict)
- # Make the final line parameters. Its annoying because
- # one or both regions may have fit to nothing
- if np.size(p1)> 0 and np.size(p2)>0:
- p = np.r_[p1,p2]
- elif np.size(p1) > 0:
- p = p1
- else:
- p = p2
+ p2,flag = _complex_fit(xb[b2:], yDat[b2:], yFit[b2:],
+ initz, minSize, errorBound, speciesDict)
+
+ # Make the final line parameters. Its annoying because
+ # one or both regions may have fit to nothing
+ if np.size(p1)> 0 and np.size(p2)>0:
+ p = np.r_[p1,p2]
+ elif np.size(p1) > 0:
+ p = p1
+ else:
+ p = p2
return p
@@ -939,7 +939,7 @@
Parameters
----------
p : (3,) ndarray
- array with form [[N1, b1, z1], ...]
+ array with form [[N1, b1, z1], ...]
speciesDict : dictionary
dictionary with properties giving the max and min
values appropriate for each parameter N,b, and z.
@@ -979,7 +979,7 @@
# we just set it to no fit because we've tried
# everything else at this point. this region just sucks :(
remove_lines.append(i)
-
+
if remove_lines:
p = np.delete(p, remove_lines, axis=0)
@@ -988,7 +988,7 @@
def _output_fit(lineDic, file_name = 'spectrum_fit.h5'):
"""
This function is designed to output the parameters of the series
- of lines used to fit an absorption spectrum.
+ of lines used to fit an absorption spectrum.
The dataset contains entries in the form species/N, species/b
species/z, and species/complex. The ith entry in each of the datasets
@@ -999,11 +999,11 @@
Parameters
----------
lineDic : dictionary
- Dictionary of dictionaries representing the fit lines.
+ Dictionary of dictionaries representing the fit lines.
Top level keys are the species given in orderFits and the corresponding
- entries are dictionaries with the keys 'N','b','z', and 'group#'.
+ entries are dictionaries with the keys 'N','b','z', and 'group#'.
Each of these corresponds to a list of the parameters for every
- accepted fitted line.
+ accepted fitted line.
fileName : string, optional
Name of the file to output fit to. Default = 'spectrum_fit.h5'
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.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.spacepope.org/pipermail/yt-svn-spacepope.org/attachments/20160420/36ff654d/attachment.html>
More information about the yt-svn
mailing list