Skip to content

Commit

Permalink
address comments
Browse files Browse the repository at this point in the history
  • Loading branch information
jvshields committed Feb 19, 2024
1 parent 4443ea4 commit 589686e
Show file tree
Hide file tree
Showing 2 changed files with 31 additions and 25 deletions.
22 changes: 11 additions & 11 deletions stardis/radiation_field/opacities/opacities_solvers/base.py
Original file line number Diff line number Diff line change
Expand Up @@ -33,7 +33,7 @@


# H minus opacity
def calc_alpha_file(stellar_plasma, stellar_model, tracing_nus, spec, fpath):
def calc_alpha_file(stellar_plasma, stellar_model, tracing_nus, opacity_source, fpath):
"""
Calculates opacities when a cross-section file is provided.
Expand All @@ -43,7 +43,7 @@ def calc_alpha_file(stellar_plasma, stellar_model, tracing_nus, spec, fpath):
stellar_model : stardis.model.base.StellarModel
tracing_nus : astropy.unit.quantity.Quantity
Numpy array of frequencies used for ray tracing with units of Hz.
species : spec : str
opacity_source : str
String representing the opacity source. Used to obtain the appropriate number density of the corresponding species.
fpath : str
Path to the cross-section file.
Expand All @@ -58,11 +58,11 @@ def calc_alpha_file(stellar_plasma, stellar_model, tracing_nus, spec, fpath):
tracing_lambdas = tracing_nus.to(u.AA, u.spectral()).value

sigmas = sigma_file(
tracing_lambdas, stellar_model.temperatures.value, Path(fpath), spec
tracing_lambdas, stellar_model.temperatures.value, Path(fpath), opacity_source
)
number_density, atomic_number, ion_number = get_number_density(
stellar_plasma, spec
) # Should revisit this function to make it more general and less hacky.
stellar_plasma, opacity_source
) ###TODO: Should revisit this function to make it more general and less hacky.
return sigmas * number_density.to_numpy()[:, np.newaxis]


Expand Down Expand Up @@ -531,18 +531,18 @@ def calc_alphas(
"""

for (
spec,
opacity_source,
fpath,
) in opacity_config.file.items(): # Iterate through requested file opacities
alpha_file = calc_alpha_file(
stellar_plasma,
stellar_model,
stellar_radiation_field.frequencies,
spec,
opacity_source,
fpath,
)
stellar_radiation_field.opacities.opacities_dict[
f"alpha_file_{spec}"
f"alpha_file_{opacity_source}"
] = alpha_file

alpha_bf = calc_alpha_bf(
Expand Down Expand Up @@ -583,9 +583,9 @@ def calc_alphas(
stellar_radiation_field.frequencies,
opacity_config.line,
)
stellar_radiation_field.opacities.opacities_dict[
"alpha_line_at_nu"
] = alpha_line_at_nu
stellar_radiation_field.opacities.opacities_dict["alpha_line_at_nu"] = (
alpha_line_at_nu
)
stellar_radiation_field.opacities.opacities_dict["alpha_line_at_nu_gammas"] = gammas
stellar_radiation_field.opacities.opacities_dict[
"alpha_line_at_nu_doppler_widths"
Expand Down
34 changes: 20 additions & 14 deletions stardis/radiation_field/opacities/opacities_solvers/util.py
Original file line number Diff line number Diff line change
Expand Up @@ -7,7 +7,7 @@
from scipy.interpolate import interp1d, LinearNDInterpolator


def sigma_file(tracing_lambdas, temperatures, fpath, spec=None):
def sigma_file(tracing_lambdas, temperatures, fpath, opacity_source=None):
"""
Reads and interpolates a cross-section file.
Expand All @@ -19,7 +19,7 @@ def sigma_file(tracing_lambdas, temperatures, fpath, spec=None):
Temperatures to compute the cross-section for.
fpath : str
Filepath to cross-section file.
spec: str
opacity_source: str
Opacity source. Used to identify the table to be read appropriately.
Returns
Expand All @@ -29,7 +29,7 @@ def sigma_file(tracing_lambdas, temperatures, fpath, spec=None):
for each wavelength and temperature combination.
"""
if (
spec == "H2plus_bf"
opacity_source == "H2plus_bf"
): # This section specifically ingests the Stancil 1994 h2_plus_bf_S1994.dat table found in data.
h2_plus_bf_table = pd.read_csv(fpath, delimiter="\s+", index_col=0, comment="#")
h2_plus_bf_table.replace({"-": "e-"}, regex=True, inplace=True)
Expand All @@ -51,7 +51,7 @@ def sigma_file(tracing_lambdas, temperatures, fpath, spec=None):
) # Scaling from Stancil 1994 table

elif (
spec == "Hminus_ff"
opacity_source == "Hminus_ff"
): # This section specifically ingests the Bell and Berrington 1987 h_minus_ff_B1987.dat table found in data.
h_minus_ff_table = pd.read_csv(fpath, delimiter="\s+", comment="#")
h_minus_ff_table.columns = h_minus_ff_table.columns.str.strip(",")
Expand All @@ -74,7 +74,10 @@ def sigma_file(tracing_lambdas, temperatures, fpath, spec=None):
* const.k_B.cgs.value
* temperatures[:, np.newaxis]
)
else: # This is currently used to read h_minus_bf_W1979.dat, Wishart 1979

elif (
opacity_source == "Hminus_bf"
): # This is currently used to read h_minus_bf_W1979.dat, Wishart 1979
h_minus_bf_table = pd.read_csv(
fpath, header=None, comment="#", names=["wavelength", "cross_section"]
)
Expand All @@ -89,10 +92,13 @@ def sigma_file(tracing_lambdas, temperatures, fpath, spec=None):
)
sigmas = linear_interp_1d_from_file(tracing_lambdas)

else:
raise ValueError(f"Unknown opacity_source: {opacity_source}")

return sigmas


def get_number_density(stellar_plasma, spec):
def get_number_density(stellar_plasma, opacity_source):
"""
Computes number density, atomic number, and ion number for an opacity
source provided as a string.
Expand All @@ -107,41 +113,41 @@ def get_number_density(stellar_plasma, spec):
ion_number : int
"""

if spec == "Hminus_bf":
if opacity_source == "Hminus_bf":
return stellar_plasma.h_minus_density, None, None
elif spec == "Hminus_ff":
elif opacity_source == "Hminus_ff":
return (
stellar_plasma.ion_number_density.loc[1, 0]
* stellar_plasma.electron_densities,
None,
None,
)
elif spec == "Heminus_ff":
elif opacity_source == "Heminus_ff":
return (
stellar_plasma.ion_number_density.loc[2, 0]
* stellar_plasma.electron_densities,
None,
None,
)
elif spec == "H2minus_ff":
elif opacity_source == "H2minus_ff":
return stellar_plasma.h2_density * stellar_plasma.electron_densities, None, None
elif spec == "H2plus_ff":
elif opacity_source == "H2plus_ff":
return (
stellar_plasma.ion_number_density.loc[1, 0]
* stellar_plasma.ion_number_density.loc[1, 1],
None,
None,
)
elif spec == "H2plus_bf":
elif opacity_source == "H2plus_bf":
return stellar_plasma.h2_plus_density, None, None

ion = spec[: len(spec) - 3]
ion = opacity_source[: len(opacity_source) - 3]

atomic_number, ion_number = species_string_to_tuple(ion.replace("_", " "))

number_density = 1

if spec[len(spec) - 2 :] == "ff":
if opacity_source[len(opacity_source) - 2 :] == "ff":
ion_number += 1
number_density *= stellar_plasma.electron_densities

Expand Down

0 comments on commit 589686e

Please sign in to comment.