Module diskchef.maps.radmc_lines
Expand source code
import warnings
from contextlib import redirect_stdout
import typing
import io
from datetime import timedelta
import subprocess
import time
import shutil
import os
from dataclasses import dataclass
import numpy as np
from astropy.io.fits import Header
from matplotlib.figure import Figure
import astropy.coordinates
from astropy import units as u
from astropy import constants as c
from astropy.wcs import WCS
from spectral_cube import SpectralCube
import radmc3dPy
import diskchef
from diskchef.engine.exceptions import CHEFNotImplementedError
from diskchef.engine.other import PathLike
from diskchef.lamda.line import Line, DummyLine
from diskchef.maps.radmcrt import RadMCBase, RadMCOutput
@dataclass
class RadMCRTImage(RadMCBase):
"""
Class with interface to run RadMC3D by Cornelis Dullemond
By initialization creates the basic table (`self.polar_table`) with the grid in polar RadMC3D coordinates.
Runs a generation of a continuum image, can be subclassed for line cube
"""
allowed_modes = ["image"]
velocity: u.km / u.s = 0 * u.km / u.s
mode: typing.Literal["image"] = 'image'
def __post_init__(self):
super().__post_init__()
self.check_mode()
def check_mode(self):
if self.mode not in self.allowed_modes:
self.logger.error(f"`mode` should be one of {self.allowed_modes}. Undefined behavior is expected")
@property
def command(self) -> str:
return (
f"{self.executable} {self.mode} "
)
@u.quantity_input
def run(
self, *,
inclination: u.deg = 0 * u.deg, position_angle: u.deg = 0 * u.deg,
distance: u.pc = 140 * u.pc, wav=1 * u.mm,
threads: int = 1, npix: int = 100,
) -> None:
"""Run RadMC3D after files were created with `create_files()`"""
self.logger.info("Running radmc3d")
start = time.time()
wavelength = wav.to(u.um, equivalencies=u.spectral())
command = (self.command +
f"incl {inclination.to_value(u.deg)} "
f"posang {position_angle.to_value(u.deg)} "
f"setthreads {threads} "
f"npix {npix} "
f"lambda {wavelength.value}"
)
self.logger.info("Running radmc3d for wavelength %.4e um: %s", wavelength.value, command)
proc = subprocess.run(
command,
cwd=self.folder,
text=True,
capture_output=True,
shell=True
)
self.logger.info("radmc3d finished after %s", timedelta(seconds=time.time() - start))
self.catch_radmc_messages(proc)
newname = os.path.join(self.folder, f"{wavelength.value:.4e}_um_image.out")
shutil.move(os.path.join(self.folder, "image.out"), newname)
self.radmc_to_fits(
name=newname,
line=wavelength,
distance=distance
)
def radmc_to_fits(
self, name: PathLike, line: typing.Union[Line, u.Quantity], distance,
) -> PathLike:
"""Saves RadMC3D `image.out` files as FITS files
Args:
name: PathLike -- image.out-like file to process
line: Union[Line, wavelength-like u.Quantity]
Returns:
PathLike object of a newly created fits files
"""
with redirect_stdout(io.StringIO()) as f:
im = radmc3dPy.image.readImage(fname=name)
self.logger.debug(f.getvalue())
if isinstance(line, Line):
restfreq = line.frequency.to_value(u.Hz)
else:
restfreq = line.to_value(u.Hz, equivalencies=u.spectral())
x_deg = ((im.x << u.cm) / distance).to(u.deg, equivalencies=u.dimensionless_angles())
y_deg = ((im.y << u.cm) / distance).to(u.deg, equivalencies=u.dimensionless_angles())
freq_hz = (im.freq << u.Hz)
x_len = len(x_deg)
midx_i = x_len // 2
midx_val = x_deg[midx_i]
mean_x_width = (x_deg[1:] - x_deg[:-1]).mean()
y_len = len(y_deg)
midy_i = y_len // 2
midy_val = y_deg[midy_i]
mean_y_width = (y_deg[1:] - y_deg[:-1]).mean()
freq_len = len(freq_hz)
midfreq_i = freq_len // 2
midfreq_hz = freq_hz[midfreq_i]
mean_channel_width = (freq_hz[1:] - freq_hz[:-1]).mean()
coordinate = self.coordinate
if coordinate is not None:
try:
coordinate = astropy.coordinates.SkyCoord(coordinate)
except ValueError:
coordinate = astropy.coordinates.SkyCoord.from_name(coordinate)
else:
coordinate = astropy.coordinates.SkyCoord("0h 0d")
wcs_dict = {
'CTYPE3': 'RA---TAN', 'CUNIT3': 'deg',
'CDELT3': -mean_x_width.value, 'CRPIX3': midx_i + 1,
'CRVAL3': midx_val.value + coordinate.ra.to_value(u.deg), 'NAXIS3': x_len,
'CTYPE2': 'DEC--TAN', 'CUNIT2': 'deg',
'CDELT2': mean_y_width.value, 'CRPIX2': midy_i + 1,
'CRVAL2': midy_val.value + coordinate.dec.to_value(u.deg), 'NAXIS2': y_len,
'CTYPE1': 'FREQ ', 'CUNIT1': 'Hz',
'CDELT1': mean_channel_width.value if np.isfinite(mean_channel_width.value) else 1,
'CRPIX1': midfreq_i + 1, 'CRVAL1': (midfreq_hz * (1 - self.velocity / c.c)).to_value(u.Hz),
'NAXIS1': freq_len,
}
header = WCS(wcs_dict).to_header()
header["HISTORY"] = "Created with DiskCheF package"
header["HISTORY"] = "(G.V. Smirnov-Pinchukov, https://gitlab.com/SmirnGreg/diskchef)"
header["HISTORY"] = "using RadMC3D (C. Dullemond, https://github.com/dullemond/radmc3d-2.0)"
header["RESTFREQ"] = restfreq
self._extra_header(header)
cube = SpectralCube(
data=np.rot90(im.image * self.factor, axes=[0, 1]) << self.unit,
wcs=WCS(wcs_dict), header=header
)
fitsname = name.replace(".out", self.fitsname_suffix + ".fits")
cube.write(fitsname, overwrite=True)
self.fitsfiles.append(fitsname)
self.logger.info("Saved as %s and %s", name, fitsname)
self.logger.info("Saved as %s and %s", name, fitsname)
return fitsname
@property
def unit(self) -> u.Unit:
"""Unit of the returned data"""
return u.Jy / u.sr
@property
def factor(self) -> float:
"""Factor of the returned data relative to self.unit"""
return 1e23
def _extra_header(self, header: Header, **kwargs):
"""Add extra header items into the output fits files"""
for key, value in kwargs.items():
try:
header[key] = value
except Exception as e:
self.logger.error(e)
@property
def fitsname_suffix(self) -> str:
"""Additional suffix of an output fits file"""
return ""
@dataclass
class RadMCRT(RadMCRTImage):
"""
Class with interface to run RadMC3D by Cornelis Dullemond
By initialization creates the basic table (`self.polar_table`) with the grid in polar RadMC3D coordinates.
Runs a generation of a line cube
"""
def __post_init__(self):
super().__post_init__()
self.interpolate("n(H+2H2)")
def create_files(self) -> None:
"""Creates all the files necessary to run RadMC3D"""
super().create_files()
self.lines()
self.gas_velocity()
self.gas_temperature()
for molecule in self.molecules_list:
self.numberdens(species=molecule)
self.molecule(species=molecule, lamda_file=self.lamda_files_dict[molecule])
self.logger.info("Line files written to %s", self.folder)
def gas_temperature(self, out_file: PathLike = None) -> None:
"""Writes the gas temperature files"""
if out_file is None:
out_file = os.path.join(self.folder, 'gas_temperature.inp')
self.interpolate('Gas temperature')
with open(out_file, 'w') as file:
print('1', file=file) # Typically 1 at present
print(self.nrcells, file=file)
print('\n'.join(f"{entry:.7e}" for entry in self.polar_table['Gas temperature'].to(u.K).value), file=file)
def numberdens(self, species: str, out_file: PathLike = None) -> None:
"""Writes the gas number density files"""
if out_file is None:
out_file = os.path.join(self.folder, f'numberdens_{species}.inp')
self.interpolate(species)
with open(out_file, 'w') as file:
print('1', file=file) # Typically 1 at present
print(self.nrcells, file=file)
print(
'\n'.join(
f"{entry:.7e}" for entry
in self.polar_table['n(H+2H2)'].to(u.cm ** (-3)).value * self.polar_table[species]),
file=file
)
def molecule(self, species: str, lamda_file: PathLike, out_file: PathLike = None) -> None:
"""
Copies the molecule transition files into working directory
species: str name of the molecule
"""
if out_file is None:
out_file = os.path.join(self.folder, f'molecule_{species}.inp')
shutil.copyfile(lamda_file, out_file)
def lines(self, out_file: PathLike = None) -> None:
self.molecules_list = sorted(list(set([line.molecule for line in self.line_list])))
self.lamda_files_dict = {
line.molecule: line.lamda_file for line in self.line_list
}
if out_file is None:
out_file = os.path.join(self.folder, 'lines.inp')
with open(out_file, 'w') as file:
print('2', file=file) # 2 for molecules
print(len(self.molecules_list), file=file) # unique list of molecules
for molecule in self.molecules_list:
coll_partners = next(line for line in self.line_list if line.molecule == molecule).collision_partner
print(f'{molecule} leiden 0 0 {len(coll_partners)}', file=file)
print('\n'.join(coll_partners), file=file)
def gas_velocity(self, out_file: PathLike = None) -> None:
if out_file is None:
out_file = os.path.join(self.folder, 'gas_velocity.inp')
with open(out_file, 'w') as file:
print('1', file=file) # Typically 1 at present
print(self.nrcells, file=file)
for vr, vtheta, vphi in zip(self.polar_table['Velocity R'].to(u.cm / u.s).value,
self.polar_table['Velocity Theta'].to(u.cm / u.s).value,
self.polar_table['Velocity Phi'].to(u.cm / u.s).value):
print(f"{vr:.7e}", f"{vtheta:.7e}", f"{vphi:.7e}", file=file)
@u.quantity_input
def run(
self,
inclination: u.deg = 0 * u.deg, position_angle: u.deg = 0 * u.deg,
distance: u.pc = 140 * u.pc, velocity_offset: u.km / u.s = 0 * u.km / u.s,
threads: int = 1, npix: int = 100,
) -> None:
"""Run RadMC3D after files were created with `create_files()`"""
self.logger.info("Running radmc3d")
for line in self.line_list:
self._run_single(
molecule=self.molecules_list.index(line.molecule) + 1,
line=line.transition,
inclination=inclination.to(u.deg).value, position_angle=position_angle.to(u.deg).value,
velocity_offset=velocity_offset.to(u.km / u.s).value,
name=f"{line.name}_image.out",
distance=distance.to(u.pc).value,
threads=threads,
lineobj=line,
npix=npix,
# coordinate=coordinate,
)
def _run_single(
self, molecule: int, line: int, inclination: float = 0, position_angle: float = 0,
name: PathLike = None, distance: float = 140, velocity_offset: float = 0,
n_channels: int = 100, threads: int = 1, lineobj: Line = None, npix: int = 100,
) -> None:
start = time.time()
command = (self.command +
f"imolspec {molecule} "
f"iline {line} "
f"widthkms 10 "
f"incl {inclination} "
f"posang {position_angle} "
f"vkms {velocity_offset} "
f"linenlam {n_channels} "
f"setthreads {threads} "
f"npix {npix} "
)
self.logger.info("Running radmc3d for molecule %d and transition %d: %s", molecule, line, command)
proc = subprocess.run(
command,
cwd=self.folder,
text=True,
capture_output=True,
shell=True
)
self.logger.info("radmc3d finished after %s", timedelta(seconds=time.time() - start))
self.catch_radmc_messages(proc)
output = RadMCOutput(lineobj, mode=self.mode)
if name is not None:
newname = os.path.join(self.folder, name)
shutil.move(os.path.join(self.folder, "image.out"), newname)
output.file_radmc = newname
output.file_fits = self.radmc_to_fits(newname, self.line_list[line], distance * u.pc)
self.outputs[lineobj] = output
def copy_for_propype(self, folder: PathLike = None) -> None:
"""Creates a copy of `self.fitsfiles` in a format which is ready to run them through PRODIGE pipeline"""
if folder is None:
folder = os.path.join(self.folder, "propype")
try:
os.mkdir(folder)
except FileExistsError as e:
self.logger.debug("%s", e)
for line, output in self.outputs.items():
molecule = line.molecule
fitsfile = output.file_fits
transition = output.line.transition
shutil.copyfile(fitsfile, os.path.join(folder, f"m-Line-00-{molecule}_{transition}+D.fits"))
@dataclass
class RadMCRTLines(RadMCRT):
"""
Subclass of RadMCRT to run RadMC only once for all the required frequencies
Usage:
>>> physics = diskchef.physics.WilliamsBest2014(star_mass=0.52 * u.solMass, radial_bins=10, vertical_bins=10)
>>> chem = diskchef.chemistry.NonzeroChemistryWB2014(physics)
>>> chem.run_chemistry()
>>> mapping = RadMCRTLines(chemistry=chem, line_list=[
... Line(name='CO J=2-1', transition=1, molecule='CO'),
... Line(name='CO J=3-2', transition=2, molecule='CO'),
... ])
>>> mapping.wavelength_multiple_lines(channels_per_line=5, window_width=2 * u.km/u.s) # doctest: +NORMALIZE_WHITESPACE
array([2600.76630869, 2600.76197107, 2600.75763346, 2600.75329588,
2600.7489583 , 1300.40799349, 1300.40582464, 1300.4036558 ,
1300.40148696, 1300.39931813])
"""
def frequency_centers(self) -> typing.List[u.Quantity]:
"""Fetches frequencies of lines from `self.line_list`"""
for line in self.line_list:
line.parse_lamda()
self.ordered_line_list = sorted(self.line_list, key=lambda line: line.frequency)
frequency_centers = [line.frequency for line in self.ordered_line_list]
return frequency_centers
def wavelength_multiple_lines(self, channels_per_line: int, window_width: u.km / u.s) -> np.ndarray:
"""
Return a sorted array of wavelength for all the lines, in um (unitless)
"""
frequencies_list = []
self.channels_per_line = []
for frequency_center in self.frequency_centers():
width = np.abs(window_width.to(u.GHz, equivalencies=u.doppler_radio(frequency_center)) - frequency_center)
frequencies_list.append(
np.linspace(frequency_center - width / 2, frequency_center + width / 2, channels_per_line)
)
self.channels_per_line.append(channels_per_line)
frequencies = np.hstack(frequencies_list)
wavelengths = frequencies.to(u.um, equivalencies=u.spectral()).value
return wavelengths
@u.quantity_input
def camera_wavelength_micron(
self,
out_file: PathLike = None,
window_width: u.km / u.s = 6 * u.km / u.s,
channels_per_line: int = 200
) -> None:
"""Creates a `camera_wavelength_micron.inp` file with all the frequencies for all the lines"""
wavelengths = self.wavelength_multiple_lines(channels_per_line, window_width)
if out_file is None:
out_file = os.path.join(self.folder, 'camera_wavelength_micron.inp')
with open(out_file, 'w') as file:
print(len(wavelengths), file=file)
print('\n'.join(f"{entry:.10e}" for entry in wavelengths), file=file)
def create_files(self,
window_width: u.km / u.s = 6 * u.km / u.s,
channels_per_line: int = 200
) -> None:
super().create_files()
self.camera_wavelength_micron(window_width=window_width, channels_per_line=channels_per_line)
@u.quantity_input
def run(
self,
inclination: u.deg = 0 * u.deg, position_angle: u.deg = 0 * u.deg,
distance: u.pc = 140 * u.pc, velocity_offset: u.km / u.s = 0 * u.km / u.s,
threads: int = 1, npix: int = 100
) -> None:
"""
Run radmc3d to create the line emission cubes
Args:
inclination:
position_angle:
distance:
velocity_offset:
threads:
npix:
"""
self.logger.info("Running radmc3d")
start = time.time()
command = (self.command +
f"incl {inclination.to(u.deg).value} "
f"posang {position_angle.to(u.deg).value} "
f"setthreads {threads} "
f"npix {npix} "
"loadlambda "
)
self.logger.info("Running radmc3d for all transition at once: %s", command)
proc = subprocess.run(
command,
cwd=self.folder,
text=True,
capture_output=True,
shell=True
)
self.logger.info("radmc3d finished after %s", timedelta(seconds=time.time() - start))
self.catch_radmc_messages(proc)
names = [os.path.join(self.folder, f"{line.name}_image.out") for line in self.ordered_line_list]
self.split(names=names)
for line, name in zip(self.ordered_line_list, names):
self.outputs[line] = RadMCOutput(line, file_radmc=name, distance=distance, mode=self.mode)
self.outputs[line].file_fits = self.radmc_to_fits(name, line, distance)
def channel_maps(
self,
filename: PathLike = None,
extension: typing.Union[None, str] = 'png',
distance: u.pc = None,
**kwargs
) -> typing.List[Figure]:
"""
Create channel maps for each line or just a given fits file
This method is in unstable state and will be modified and fixed in future versions
Args:
filename: input fits file to create channel maps from
extension: extension for the output figure, e.g. png, pdf, etc. Does not create any file if None
distance (u.pc): distance to the source. If None, the distance will be searched with Simbad
Returns:
list of matplotlib figures containing the channel maps
"""
if filename is None:
outputs = self.outputs
else:
line = DummyLine(name=filename, molecule=filename)
outputs: typing.Dict[DummyLine, RadMCOutput] = {
line: RadMCOutput(line=line, file_fits=filename, distance=distance)
}
figures = []
for output in outputs.values():
fig = output.plot_channel_map(**kwargs)
if extension is not None:
fig.savefig(os.path.join(self.folder, f"{output.line.name}.{extension}"))
figures.append(fig)
return figures
def split(
self, filename=None,
names: typing.List[PathLike] = None
):
"""Split the `image.out` file with all lines into multiple single-frequency files"""
lengths = self.channels_per_line
if filename is None:
filename = os.path.join(self.folder, 'image.out')
with open(filename, 'r') as bigfile:
files = [open(name, 'w') for name in names]
header = [next(bigfile) for _ in range(4)]
px_number = [int(entry) for entry in header[1].split()]
self.logger.info("Splitting %s in %s", filename, names)
for file, length in zip(files, lengths):
file.write(header[0])
file.write(header[1])
file.write(str(length) + "\n")
file.write(header[3])
for line in range(length):
file.write(next(bigfile))
for file, length in zip(files, lengths):
for line in range(length * (px_number[0] * px_number[1] + 1)):
file.write(next(bigfile))
for file in files:
file.close()
@dataclass
class RadMCRTSingleCall(RadMCRTLines):
"""Deprecated in favor of RadMCRTLines"""
def __post_init__(self):
super().__post_init__()
self.logger.warning("Deprecation warning. Renamed to RadMCRTLines!")
@dataclass
class RadMCRTLinesTraceTau(RadMCRTLines):
"""Class to calculate optical depth map in lines"""
allowed_modes = ["image tracetau"]
mode: typing.Literal["image tracetau"] = "image tracetau"
@property
def unit(self) -> u.Unit:
"""Unit of the returned data"""
return u.dimensionless_unscaled
@property
def factor(self) -> float:
"""Factor of the returned data relative to self.unit"""
return 1.
def _extra_header(self, header: Header, **kwargs):
"""Add extra header items into the output fits files"""
super()._extra_header(header, **{"HISTORY": "Optical depth map"}, **kwargs)
@property
def fitsname_suffix(self) -> str:
"""Additional suffix of an output fits file"""
return "_tracetau"
@dataclass
class RadMCRTLinesTauSurf(RadMCRTLines):
"""Class to calculate optical depth map in lines"""
allowed_modes = ["tausurf"]
mode: typing.Literal["image tracetau"] = "tausurf"
tau: float = 1.
@property
def unit(self) -> u.Unit:
"""Unit of the returned data"""
return u.au
@property
def command(self) -> str:
return super().command + f"{self.tau} "
@property
def factor(self) -> float:
"""Factor of the returned data relative to self.unit"""
return (u.cm / u.au).si.scale
def _extra_header(self, header: Header, **kwargs):
"""Add extra header items into the output fits files"""
super()._extra_header(header, **{"HISTORY": f"Optical depth tau={self.tau: .2e} map"}, **kwargs)
@property
def fitsname_suffix(self) -> str:
"""Additional suffix of an output fits file"""
return f"_tau_{self.tau: .2e}"
Classes
class RadMCRT (chemistry: Union[ChemistryBase, PhysicsBase] = None, line_list: List[Line] = None, executable: Union[str, os.PathLike] = 'radmc3d', folder: Union[str, os.PathLike] = 'radmc', radii_bins: Optional[None] = None, theta_bins: Optional[None] = None, outer_radius: Optional[None] = None, verbosity: int = 0, wavelengths: astropy.units.quantity.Quantity = <Quantity [9.12000000e-02, 1.00185062e-01, 1.10055336e-01, 1.20898034e-01, 1.32808958e-01, 1.45893352e-01, 1.60266826e-01, 1.76056380e-01, 1.93401528e-01, 2.12455528e-01, 2.33386735e-01, 2.56380093e-01, 2.81638767e-01, 3.09385935e-01, 3.39866766e-01, 3.73350582e-01, 4.10133237e-01, 4.50539734e-01, 4.94927096e-01, 5.43687520e-01, 5.97251841e-01, 6.56093341e-01, 7.20731931e-01, 7.91738741e-01, 8.69741171e-01, 9.55428434e-01, 1.04955764e+00, 1.15296050e+00, 1.26655065e+00, 1.39133174e+00, 1.52840633e+00, 1.67898555e+00, 1.84439991e+00, 2.02611096e+00, 2.22572425e+00, 2.44500354e+00, 2.68588630e+00, 2.95050094e+00, 3.24118552e+00, 3.56050847e+00, 3.91129125e+00, 4.29663329e+00, 4.71993939e+00, 5.18494979e+00, 5.69577320e+00, 6.25692315e+00, 6.87335783e+00, 7.55052390e+00, 8.29440466e+00, 9.11157287e+00, 1.00092488e+01, 1.09953642e+01, 1.20786320e+01, 1.32686239e+01, 1.45758543e+01, 1.60118735e+01, 1.75893699e+01, 1.93222820e+01, 2.12259213e+01, 2.33171079e+01, 2.56143191e+01, 2.81378525e+01, 3.09100054e+01, 3.39552720e+01, 3.73005596e+01, 4.09754262e+01, 4.50123423e+01, 4.94469770e+01, 5.43185138e+01, 5.96699964e+01, 6.55487093e+01, 7.20065955e+01, 7.91007152e+01, 8.68937507e+01, 9.54545592e+01, 1.04858782e+02, 1.15189513e+02, 1.26538032e+02, 1.39004611e+02, 1.52699404e+02, 1.67743413e+02, 1.84269563e+02, 2.02423878e+02, 2.22366762e+02, 2.44274428e+02, 2.68340447e+02, 2.94777460e+02, 3.23819058e+02, 3.55721846e+02, 3.90767711e+02, 4.29266308e+02, 4.71557804e+02, 5.18015875e+02, 5.69051015e+02, 6.25114158e+02, 6.86700665e+02, 7.54354701e+02, 8.28674040e+02, 9.10315352e+02, 1.00000000e+03] um>, modified_random_walk: bool = True, scattering_mode_max: int = None, nphot_therm: int = None, nphot_mono: int = None, nphot_spec: int = None, coordinate: Union[str, astropy.coordinates.sky_coordinate.SkyCoord] = None, external_source_type: Literal['Draine1978', 'WeingartnerDraine2001', 'None'] = None, external_source_multiplier: float = 1, velocity: Unit("km / s") = <Quantity 0. km / s>, mode: Literal['image'] = 'image')
-
Class with interface to run RadMC3D by Cornelis Dullemond
By initialization creates the basic table (
self.polar_table
) with the grid in polar RadMC3D coordinates.Runs a generation of a line cube
Expand source code
@dataclass class RadMCRT(RadMCRTImage): """ Class with interface to run RadMC3D by Cornelis Dullemond By initialization creates the basic table (`self.polar_table`) with the grid in polar RadMC3D coordinates. Runs a generation of a line cube """ def __post_init__(self): super().__post_init__() self.interpolate("n(H+2H2)") def create_files(self) -> None: """Creates all the files necessary to run RadMC3D""" super().create_files() self.lines() self.gas_velocity() self.gas_temperature() for molecule in self.molecules_list: self.numberdens(species=molecule) self.molecule(species=molecule, lamda_file=self.lamda_files_dict[molecule]) self.logger.info("Line files written to %s", self.folder) def gas_temperature(self, out_file: PathLike = None) -> None: """Writes the gas temperature files""" if out_file is None: out_file = os.path.join(self.folder, 'gas_temperature.inp') self.interpolate('Gas temperature') with open(out_file, 'w') as file: print('1', file=file) # Typically 1 at present print(self.nrcells, file=file) print('\n'.join(f"{entry:.7e}" for entry in self.polar_table['Gas temperature'].to(u.K).value), file=file) def numberdens(self, species: str, out_file: PathLike = None) -> None: """Writes the gas number density files""" if out_file is None: out_file = os.path.join(self.folder, f'numberdens_{species}.inp') self.interpolate(species) with open(out_file, 'w') as file: print('1', file=file) # Typically 1 at present print(self.nrcells, file=file) print( '\n'.join( f"{entry:.7e}" for entry in self.polar_table['n(H+2H2)'].to(u.cm ** (-3)).value * self.polar_table[species]), file=file ) def molecule(self, species: str, lamda_file: PathLike, out_file: PathLike = None) -> None: """ Copies the molecule transition files into working directory species: str name of the molecule """ if out_file is None: out_file = os.path.join(self.folder, f'molecule_{species}.inp') shutil.copyfile(lamda_file, out_file) def lines(self, out_file: PathLike = None) -> None: self.molecules_list = sorted(list(set([line.molecule for line in self.line_list]))) self.lamda_files_dict = { line.molecule: line.lamda_file for line in self.line_list } if out_file is None: out_file = os.path.join(self.folder, 'lines.inp') with open(out_file, 'w') as file: print('2', file=file) # 2 for molecules print(len(self.molecules_list), file=file) # unique list of molecules for molecule in self.molecules_list: coll_partners = next(line for line in self.line_list if line.molecule == molecule).collision_partner print(f'{molecule} leiden 0 0 {len(coll_partners)}', file=file) print('\n'.join(coll_partners), file=file) def gas_velocity(self, out_file: PathLike = None) -> None: if out_file is None: out_file = os.path.join(self.folder, 'gas_velocity.inp') with open(out_file, 'w') as file: print('1', file=file) # Typically 1 at present print(self.nrcells, file=file) for vr, vtheta, vphi in zip(self.polar_table['Velocity R'].to(u.cm / u.s).value, self.polar_table['Velocity Theta'].to(u.cm / u.s).value, self.polar_table['Velocity Phi'].to(u.cm / u.s).value): print(f"{vr:.7e}", f"{vtheta:.7e}", f"{vphi:.7e}", file=file) @u.quantity_input def run( self, inclination: u.deg = 0 * u.deg, position_angle: u.deg = 0 * u.deg, distance: u.pc = 140 * u.pc, velocity_offset: u.km / u.s = 0 * u.km / u.s, threads: int = 1, npix: int = 100, ) -> None: """Run RadMC3D after files were created with `create_files()`""" self.logger.info("Running radmc3d") for line in self.line_list: self._run_single( molecule=self.molecules_list.index(line.molecule) + 1, line=line.transition, inclination=inclination.to(u.deg).value, position_angle=position_angle.to(u.deg).value, velocity_offset=velocity_offset.to(u.km / u.s).value, name=f"{line.name}_image.out", distance=distance.to(u.pc).value, threads=threads, lineobj=line, npix=npix, # coordinate=coordinate, ) def _run_single( self, molecule: int, line: int, inclination: float = 0, position_angle: float = 0, name: PathLike = None, distance: float = 140, velocity_offset: float = 0, n_channels: int = 100, threads: int = 1, lineobj: Line = None, npix: int = 100, ) -> None: start = time.time() command = (self.command + f"imolspec {molecule} " f"iline {line} " f"widthkms 10 " f"incl {inclination} " f"posang {position_angle} " f"vkms {velocity_offset} " f"linenlam {n_channels} " f"setthreads {threads} " f"npix {npix} " ) self.logger.info("Running radmc3d for molecule %d and transition %d: %s", molecule, line, command) proc = subprocess.run( command, cwd=self.folder, text=True, capture_output=True, shell=True ) self.logger.info("radmc3d finished after %s", timedelta(seconds=time.time() - start)) self.catch_radmc_messages(proc) output = RadMCOutput(lineobj, mode=self.mode) if name is not None: newname = os.path.join(self.folder, name) shutil.move(os.path.join(self.folder, "image.out"), newname) output.file_radmc = newname output.file_fits = self.radmc_to_fits(newname, self.line_list[line], distance * u.pc) self.outputs[lineobj] = output def copy_for_propype(self, folder: PathLike = None) -> None: """Creates a copy of `self.fitsfiles` in a format which is ready to run them through PRODIGE pipeline""" if folder is None: folder = os.path.join(self.folder, "propype") try: os.mkdir(folder) except FileExistsError as e: self.logger.debug("%s", e) for line, output in self.outputs.items(): molecule = line.molecule fitsfile = output.file_fits transition = output.line.transition shutil.copyfile(fitsfile, os.path.join(folder, f"m-Line-00-{molecule}_{transition}+D.fits"))
Ancestors
Subclasses
Class variables
var velocity : Unit("km / s")
Methods
def copy_for_propype(self, folder: Union[str, os.PathLike] = None) ‑> None
-
Creates a copy of
self.fitsfiles
in a format which is ready to run them through PRODIGE pipelineExpand source code
def copy_for_propype(self, folder: PathLike = None) -> None: """Creates a copy of `self.fitsfiles` in a format which is ready to run them through PRODIGE pipeline""" if folder is None: folder = os.path.join(self.folder, "propype") try: os.mkdir(folder) except FileExistsError as e: self.logger.debug("%s", e) for line, output in self.outputs.items(): molecule = line.molecule fitsfile = output.file_fits transition = output.line.transition shutil.copyfile(fitsfile, os.path.join(folder, f"m-Line-00-{molecule}_{transition}+D.fits"))
def gas_temperature(self, out_file: Union[str, os.PathLike] = None) ‑> None
-
Writes the gas temperature files
Expand source code
def gas_temperature(self, out_file: PathLike = None) -> None: """Writes the gas temperature files""" if out_file is None: out_file = os.path.join(self.folder, 'gas_temperature.inp') self.interpolate('Gas temperature') with open(out_file, 'w') as file: print('1', file=file) # Typically 1 at present print(self.nrcells, file=file) print('\n'.join(f"{entry:.7e}" for entry in self.polar_table['Gas temperature'].to(u.K).value), file=file)
def gas_velocity(self, out_file: Union[str, os.PathLike] = None) ‑> None
-
Expand source code
def gas_velocity(self, out_file: PathLike = None) -> None: if out_file is None: out_file = os.path.join(self.folder, 'gas_velocity.inp') with open(out_file, 'w') as file: print('1', file=file) # Typically 1 at present print(self.nrcells, file=file) for vr, vtheta, vphi in zip(self.polar_table['Velocity R'].to(u.cm / u.s).value, self.polar_table['Velocity Theta'].to(u.cm / u.s).value, self.polar_table['Velocity Phi'].to(u.cm / u.s).value): print(f"{vr:.7e}", f"{vtheta:.7e}", f"{vphi:.7e}", file=file)
def lines(self, out_file: Union[str, os.PathLike] = None) ‑> None
-
Expand source code
def lines(self, out_file: PathLike = None) -> None: self.molecules_list = sorted(list(set([line.molecule for line in self.line_list]))) self.lamda_files_dict = { line.molecule: line.lamda_file for line in self.line_list } if out_file is None: out_file = os.path.join(self.folder, 'lines.inp') with open(out_file, 'w') as file: print('2', file=file) # 2 for molecules print(len(self.molecules_list), file=file) # unique list of molecules for molecule in self.molecules_list: coll_partners = next(line for line in self.line_list if line.molecule == molecule).collision_partner print(f'{molecule} leiden 0 0 {len(coll_partners)}', file=file) print('\n'.join(coll_partners), file=file)
def molecule(self, species: str, lamda_file: Union[str, os.PathLike], out_file: Union[str, os.PathLike] = None) ‑> None
-
Copies the molecule transition files into working directory
species: str name of the molecule
Expand source code
def molecule(self, species: str, lamda_file: PathLike, out_file: PathLike = None) -> None: """ Copies the molecule transition files into working directory species: str name of the molecule """ if out_file is None: out_file = os.path.join(self.folder, f'molecule_{species}.inp') shutil.copyfile(lamda_file, out_file)
def numberdens(self, species: str, out_file: Union[str, os.PathLike] = None) ‑> None
-
Writes the gas number density files
Expand source code
def numberdens(self, species: str, out_file: PathLike = None) -> None: """Writes the gas number density files""" if out_file is None: out_file = os.path.join(self.folder, f'numberdens_{species}.inp') self.interpolate(species) with open(out_file, 'w') as file: print('1', file=file) # Typically 1 at present print(self.nrcells, file=file) print( '\n'.join( f"{entry:.7e}" for entry in self.polar_table['n(H+2H2)'].to(u.cm ** (-3)).value * self.polar_table[species]), file=file )
Inherited members
class RadMCRTImage (chemistry: Union[ChemistryBase, PhysicsBase] = None, line_list: List[Line] = None, executable: Union[str, os.PathLike] = 'radmc3d', folder: Union[str, os.PathLike] = 'radmc', radii_bins: Optional[None] = None, theta_bins: Optional[None] = None, outer_radius: Optional[None] = None, verbosity: int = 0, wavelengths: astropy.units.quantity.Quantity = <Quantity [9.12000000e-02, 1.00185062e-01, 1.10055336e-01, 1.20898034e-01, 1.32808958e-01, 1.45893352e-01, 1.60266826e-01, 1.76056380e-01, 1.93401528e-01, 2.12455528e-01, 2.33386735e-01, 2.56380093e-01, 2.81638767e-01, 3.09385935e-01, 3.39866766e-01, 3.73350582e-01, 4.10133237e-01, 4.50539734e-01, 4.94927096e-01, 5.43687520e-01, 5.97251841e-01, 6.56093341e-01, 7.20731931e-01, 7.91738741e-01, 8.69741171e-01, 9.55428434e-01, 1.04955764e+00, 1.15296050e+00, 1.26655065e+00, 1.39133174e+00, 1.52840633e+00, 1.67898555e+00, 1.84439991e+00, 2.02611096e+00, 2.22572425e+00, 2.44500354e+00, 2.68588630e+00, 2.95050094e+00, 3.24118552e+00, 3.56050847e+00, 3.91129125e+00, 4.29663329e+00, 4.71993939e+00, 5.18494979e+00, 5.69577320e+00, 6.25692315e+00, 6.87335783e+00, 7.55052390e+00, 8.29440466e+00, 9.11157287e+00, 1.00092488e+01, 1.09953642e+01, 1.20786320e+01, 1.32686239e+01, 1.45758543e+01, 1.60118735e+01, 1.75893699e+01, 1.93222820e+01, 2.12259213e+01, 2.33171079e+01, 2.56143191e+01, 2.81378525e+01, 3.09100054e+01, 3.39552720e+01, 3.73005596e+01, 4.09754262e+01, 4.50123423e+01, 4.94469770e+01, 5.43185138e+01, 5.96699964e+01, 6.55487093e+01, 7.20065955e+01, 7.91007152e+01, 8.68937507e+01, 9.54545592e+01, 1.04858782e+02, 1.15189513e+02, 1.26538032e+02, 1.39004611e+02, 1.52699404e+02, 1.67743413e+02, 1.84269563e+02, 2.02423878e+02, 2.22366762e+02, 2.44274428e+02, 2.68340447e+02, 2.94777460e+02, 3.23819058e+02, 3.55721846e+02, 3.90767711e+02, 4.29266308e+02, 4.71557804e+02, 5.18015875e+02, 5.69051015e+02, 6.25114158e+02, 6.86700665e+02, 7.54354701e+02, 8.28674040e+02, 9.10315352e+02, 1.00000000e+03] um>, modified_random_walk: bool = True, scattering_mode_max: int = None, nphot_therm: int = None, nphot_mono: int = None, nphot_spec: int = None, coordinate: Union[str, astropy.coordinates.sky_coordinate.SkyCoord] = None, external_source_type: Literal['Draine1978', 'WeingartnerDraine2001', 'None'] = None, external_source_multiplier: float = 1, velocity: Unit("km / s") = <Quantity 0. km / s>, mode: Literal['image'] = 'image')
-
Class with interface to run RadMC3D by Cornelis Dullemond
By initialization creates the basic table (
self.polar_table
) with the grid in polar RadMC3D coordinates.Runs a generation of a continuum image, can be subclassed for line cube
Expand source code
@dataclass class RadMCRTImage(RadMCBase): """ Class with interface to run RadMC3D by Cornelis Dullemond By initialization creates the basic table (`self.polar_table`) with the grid in polar RadMC3D coordinates. Runs a generation of a continuum image, can be subclassed for line cube """ allowed_modes = ["image"] velocity: u.km / u.s = 0 * u.km / u.s mode: typing.Literal["image"] = 'image' def __post_init__(self): super().__post_init__() self.check_mode() def check_mode(self): if self.mode not in self.allowed_modes: self.logger.error(f"`mode` should be one of {self.allowed_modes}. Undefined behavior is expected") @property def command(self) -> str: return ( f"{self.executable} {self.mode} " ) @u.quantity_input def run( self, *, inclination: u.deg = 0 * u.deg, position_angle: u.deg = 0 * u.deg, distance: u.pc = 140 * u.pc, wav=1 * u.mm, threads: int = 1, npix: int = 100, ) -> None: """Run RadMC3D after files were created with `create_files()`""" self.logger.info("Running radmc3d") start = time.time() wavelength = wav.to(u.um, equivalencies=u.spectral()) command = (self.command + f"incl {inclination.to_value(u.deg)} " f"posang {position_angle.to_value(u.deg)} " f"setthreads {threads} " f"npix {npix} " f"lambda {wavelength.value}" ) self.logger.info("Running radmc3d for wavelength %.4e um: %s", wavelength.value, command) proc = subprocess.run( command, cwd=self.folder, text=True, capture_output=True, shell=True ) self.logger.info("radmc3d finished after %s", timedelta(seconds=time.time() - start)) self.catch_radmc_messages(proc) newname = os.path.join(self.folder, f"{wavelength.value:.4e}_um_image.out") shutil.move(os.path.join(self.folder, "image.out"), newname) self.radmc_to_fits( name=newname, line=wavelength, distance=distance ) def radmc_to_fits( self, name: PathLike, line: typing.Union[Line, u.Quantity], distance, ) -> PathLike: """Saves RadMC3D `image.out` files as FITS files Args: name: PathLike -- image.out-like file to process line: Union[Line, wavelength-like u.Quantity] Returns: PathLike object of a newly created fits files """ with redirect_stdout(io.StringIO()) as f: im = radmc3dPy.image.readImage(fname=name) self.logger.debug(f.getvalue()) if isinstance(line, Line): restfreq = line.frequency.to_value(u.Hz) else: restfreq = line.to_value(u.Hz, equivalencies=u.spectral()) x_deg = ((im.x << u.cm) / distance).to(u.deg, equivalencies=u.dimensionless_angles()) y_deg = ((im.y << u.cm) / distance).to(u.deg, equivalencies=u.dimensionless_angles()) freq_hz = (im.freq << u.Hz) x_len = len(x_deg) midx_i = x_len // 2 midx_val = x_deg[midx_i] mean_x_width = (x_deg[1:] - x_deg[:-1]).mean() y_len = len(y_deg) midy_i = y_len // 2 midy_val = y_deg[midy_i] mean_y_width = (y_deg[1:] - y_deg[:-1]).mean() freq_len = len(freq_hz) midfreq_i = freq_len // 2 midfreq_hz = freq_hz[midfreq_i] mean_channel_width = (freq_hz[1:] - freq_hz[:-1]).mean() coordinate = self.coordinate if coordinate is not None: try: coordinate = astropy.coordinates.SkyCoord(coordinate) except ValueError: coordinate = astropy.coordinates.SkyCoord.from_name(coordinate) else: coordinate = astropy.coordinates.SkyCoord("0h 0d") wcs_dict = { 'CTYPE3': 'RA---TAN', 'CUNIT3': 'deg', 'CDELT3': -mean_x_width.value, 'CRPIX3': midx_i + 1, 'CRVAL3': midx_val.value + coordinate.ra.to_value(u.deg), 'NAXIS3': x_len, 'CTYPE2': 'DEC--TAN', 'CUNIT2': 'deg', 'CDELT2': mean_y_width.value, 'CRPIX2': midy_i + 1, 'CRVAL2': midy_val.value + coordinate.dec.to_value(u.deg), 'NAXIS2': y_len, 'CTYPE1': 'FREQ ', 'CUNIT1': 'Hz', 'CDELT1': mean_channel_width.value if np.isfinite(mean_channel_width.value) else 1, 'CRPIX1': midfreq_i + 1, 'CRVAL1': (midfreq_hz * (1 - self.velocity / c.c)).to_value(u.Hz), 'NAXIS1': freq_len, } header = WCS(wcs_dict).to_header() header["HISTORY"] = "Created with DiskCheF package" header["HISTORY"] = "(G.V. Smirnov-Pinchukov, https://gitlab.com/SmirnGreg/diskchef)" header["HISTORY"] = "using RadMC3D (C. Dullemond, https://github.com/dullemond/radmc3d-2.0)" header["RESTFREQ"] = restfreq self._extra_header(header) cube = SpectralCube( data=np.rot90(im.image * self.factor, axes=[0, 1]) << self.unit, wcs=WCS(wcs_dict), header=header ) fitsname = name.replace(".out", self.fitsname_suffix + ".fits") cube.write(fitsname, overwrite=True) self.fitsfiles.append(fitsname) self.logger.info("Saved as %s and %s", name, fitsname) self.logger.info("Saved as %s and %s", name, fitsname) return fitsname @property def unit(self) -> u.Unit: """Unit of the returned data""" return u.Jy / u.sr @property def factor(self) -> float: """Factor of the returned data relative to self.unit""" return 1e23 def _extra_header(self, header: Header, **kwargs): """Add extra header items into the output fits files""" for key, value in kwargs.items(): try: header[key] = value except Exception as e: self.logger.error(e) @property def fitsname_suffix(self) -> str: """Additional suffix of an output fits file""" return ""
Ancestors
Subclasses
Class variables
var allowed_modes
var velocity : Unit("km / s")
Instance variables
var command : str
-
Expand source code
@property def command(self) -> str: return ( f"{self.executable} {self.mode} " )
var factor : float
-
Factor of the returned data relative to self.unit
Expand source code
@property def factor(self) -> float: """Factor of the returned data relative to self.unit""" return 1e23
var fitsname_suffix : str
-
Additional suffix of an output fits file
Expand source code
@property def fitsname_suffix(self) -> str: """Additional suffix of an output fits file""" return ""
var unit : astropy.units.core.Unit
-
Unit of the returned data
Expand source code
@property def unit(self) -> u.Unit: """Unit of the returned data""" return u.Jy / u.sr
Methods
def check_mode(self)
-
Expand source code
def check_mode(self): if self.mode not in self.allowed_modes: self.logger.error(f"`mode` should be one of {self.allowed_modes}. Undefined behavior is expected")
def radmc_to_fits(self, name: Union[str, os.PathLike], line: Union[Line, astropy.units.quantity.Quantity], distance) ‑> Union[str, os.PathLike]
-
Saves RadMC3D
image.out
files as FITS filesArgs
name
- PathLike – image.out-like file to process
line
- Union[Line, wavelength-like u.Quantity]
Returns
PathLike object of a newly created fits files
Expand source code
def radmc_to_fits( self, name: PathLike, line: typing.Union[Line, u.Quantity], distance, ) -> PathLike: """Saves RadMC3D `image.out` files as FITS files Args: name: PathLike -- image.out-like file to process line: Union[Line, wavelength-like u.Quantity] Returns: PathLike object of a newly created fits files """ with redirect_stdout(io.StringIO()) as f: im = radmc3dPy.image.readImage(fname=name) self.logger.debug(f.getvalue()) if isinstance(line, Line): restfreq = line.frequency.to_value(u.Hz) else: restfreq = line.to_value(u.Hz, equivalencies=u.spectral()) x_deg = ((im.x << u.cm) / distance).to(u.deg, equivalencies=u.dimensionless_angles()) y_deg = ((im.y << u.cm) / distance).to(u.deg, equivalencies=u.dimensionless_angles()) freq_hz = (im.freq << u.Hz) x_len = len(x_deg) midx_i = x_len // 2 midx_val = x_deg[midx_i] mean_x_width = (x_deg[1:] - x_deg[:-1]).mean() y_len = len(y_deg) midy_i = y_len // 2 midy_val = y_deg[midy_i] mean_y_width = (y_deg[1:] - y_deg[:-1]).mean() freq_len = len(freq_hz) midfreq_i = freq_len // 2 midfreq_hz = freq_hz[midfreq_i] mean_channel_width = (freq_hz[1:] - freq_hz[:-1]).mean() coordinate = self.coordinate if coordinate is not None: try: coordinate = astropy.coordinates.SkyCoord(coordinate) except ValueError: coordinate = astropy.coordinates.SkyCoord.from_name(coordinate) else: coordinate = astropy.coordinates.SkyCoord("0h 0d") wcs_dict = { 'CTYPE3': 'RA---TAN', 'CUNIT3': 'deg', 'CDELT3': -mean_x_width.value, 'CRPIX3': midx_i + 1, 'CRVAL3': midx_val.value + coordinate.ra.to_value(u.deg), 'NAXIS3': x_len, 'CTYPE2': 'DEC--TAN', 'CUNIT2': 'deg', 'CDELT2': mean_y_width.value, 'CRPIX2': midy_i + 1, 'CRVAL2': midy_val.value + coordinate.dec.to_value(u.deg), 'NAXIS2': y_len, 'CTYPE1': 'FREQ ', 'CUNIT1': 'Hz', 'CDELT1': mean_channel_width.value if np.isfinite(mean_channel_width.value) else 1, 'CRPIX1': midfreq_i + 1, 'CRVAL1': (midfreq_hz * (1 - self.velocity / c.c)).to_value(u.Hz), 'NAXIS1': freq_len, } header = WCS(wcs_dict).to_header() header["HISTORY"] = "Created with DiskCheF package" header["HISTORY"] = "(G.V. Smirnov-Pinchukov, https://gitlab.com/SmirnGreg/diskchef)" header["HISTORY"] = "using RadMC3D (C. Dullemond, https://github.com/dullemond/radmc3d-2.0)" header["RESTFREQ"] = restfreq self._extra_header(header) cube = SpectralCube( data=np.rot90(im.image * self.factor, axes=[0, 1]) << self.unit, wcs=WCS(wcs_dict), header=header ) fitsname = name.replace(".out", self.fitsname_suffix + ".fits") cube.write(fitsname, overwrite=True) self.fitsfiles.append(fitsname) self.logger.info("Saved as %s and %s", name, fitsname) self.logger.info("Saved as %s and %s", name, fitsname) return fitsname
Inherited members
class RadMCRTLines (chemistry: Union[ChemistryBase, PhysicsBase] = None, line_list: List[Line] = None, executable: Union[str, os.PathLike] = 'radmc3d', folder: Union[str, os.PathLike] = 'radmc', radii_bins: Optional[None] = None, theta_bins: Optional[None] = None, outer_radius: Optional[None] = None, verbosity: int = 0, wavelengths: astropy.units.quantity.Quantity = <Quantity [9.12000000e-02, 1.00185062e-01, 1.10055336e-01, 1.20898034e-01, 1.32808958e-01, 1.45893352e-01, 1.60266826e-01, 1.76056380e-01, 1.93401528e-01, 2.12455528e-01, 2.33386735e-01, 2.56380093e-01, 2.81638767e-01, 3.09385935e-01, 3.39866766e-01, 3.73350582e-01, 4.10133237e-01, 4.50539734e-01, 4.94927096e-01, 5.43687520e-01, 5.97251841e-01, 6.56093341e-01, 7.20731931e-01, 7.91738741e-01, 8.69741171e-01, 9.55428434e-01, 1.04955764e+00, 1.15296050e+00, 1.26655065e+00, 1.39133174e+00, 1.52840633e+00, 1.67898555e+00, 1.84439991e+00, 2.02611096e+00, 2.22572425e+00, 2.44500354e+00, 2.68588630e+00, 2.95050094e+00, 3.24118552e+00, 3.56050847e+00, 3.91129125e+00, 4.29663329e+00, 4.71993939e+00, 5.18494979e+00, 5.69577320e+00, 6.25692315e+00, 6.87335783e+00, 7.55052390e+00, 8.29440466e+00, 9.11157287e+00, 1.00092488e+01, 1.09953642e+01, 1.20786320e+01, 1.32686239e+01, 1.45758543e+01, 1.60118735e+01, 1.75893699e+01, 1.93222820e+01, 2.12259213e+01, 2.33171079e+01, 2.56143191e+01, 2.81378525e+01, 3.09100054e+01, 3.39552720e+01, 3.73005596e+01, 4.09754262e+01, 4.50123423e+01, 4.94469770e+01, 5.43185138e+01, 5.96699964e+01, 6.55487093e+01, 7.20065955e+01, 7.91007152e+01, 8.68937507e+01, 9.54545592e+01, 1.04858782e+02, 1.15189513e+02, 1.26538032e+02, 1.39004611e+02, 1.52699404e+02, 1.67743413e+02, 1.84269563e+02, 2.02423878e+02, 2.22366762e+02, 2.44274428e+02, 2.68340447e+02, 2.94777460e+02, 3.23819058e+02, 3.55721846e+02, 3.90767711e+02, 4.29266308e+02, 4.71557804e+02, 5.18015875e+02, 5.69051015e+02, 6.25114158e+02, 6.86700665e+02, 7.54354701e+02, 8.28674040e+02, 9.10315352e+02, 1.00000000e+03] um>, modified_random_walk: bool = True, scattering_mode_max: int = None, nphot_therm: int = None, nphot_mono: int = None, nphot_spec: int = None, coordinate: Union[str, astropy.coordinates.sky_coordinate.SkyCoord] = None, external_source_type: Literal['Draine1978', 'WeingartnerDraine2001', 'None'] = None, external_source_multiplier: float = 1, velocity: Unit("km / s") = <Quantity 0. km / s>, mode: Literal['image'] = 'image')
-
Subclass of RadMCRT to run RadMC only once for all the required frequencies
Usage:
>>> physics = diskchef.physics.WilliamsBest2014(star_mass=0.52 * u.solMass, radial_bins=10, vertical_bins=10) >>> chem = diskchef.chemistry.NonzeroChemistryWB2014(physics) >>> chem.run_chemistry() >>> mapping = RadMCRTLines(chemistry=chem, line_list=[ ... Line(name='CO J=2-1', transition=1, molecule='CO'), ... Line(name='CO J=3-2', transition=2, molecule='CO'), ... ]) >>> mapping.wavelength_multiple_lines(channels_per_line=5, window_width=2 * u.km/u.s) # doctest: +NORMALIZE_WHITESPACE array([2600.76630869, 2600.76197107, 2600.75763346, 2600.75329588, 2600.7489583 , 1300.40799349, 1300.40582464, 1300.4036558 , 1300.40148696, 1300.39931813])
Expand source code
@dataclass class RadMCRTLines(RadMCRT): """ Subclass of RadMCRT to run RadMC only once for all the required frequencies Usage: >>> physics = diskchef.physics.WilliamsBest2014(star_mass=0.52 * u.solMass, radial_bins=10, vertical_bins=10) >>> chem = diskchef.chemistry.NonzeroChemistryWB2014(physics) >>> chem.run_chemistry() >>> mapping = RadMCRTLines(chemistry=chem, line_list=[ ... Line(name='CO J=2-1', transition=1, molecule='CO'), ... Line(name='CO J=3-2', transition=2, molecule='CO'), ... ]) >>> mapping.wavelength_multiple_lines(channels_per_line=5, window_width=2 * u.km/u.s) # doctest: +NORMALIZE_WHITESPACE array([2600.76630869, 2600.76197107, 2600.75763346, 2600.75329588, 2600.7489583 , 1300.40799349, 1300.40582464, 1300.4036558 , 1300.40148696, 1300.39931813]) """ def frequency_centers(self) -> typing.List[u.Quantity]: """Fetches frequencies of lines from `self.line_list`""" for line in self.line_list: line.parse_lamda() self.ordered_line_list = sorted(self.line_list, key=lambda line: line.frequency) frequency_centers = [line.frequency for line in self.ordered_line_list] return frequency_centers def wavelength_multiple_lines(self, channels_per_line: int, window_width: u.km / u.s) -> np.ndarray: """ Return a sorted array of wavelength for all the lines, in um (unitless) """ frequencies_list = [] self.channels_per_line = [] for frequency_center in self.frequency_centers(): width = np.abs(window_width.to(u.GHz, equivalencies=u.doppler_radio(frequency_center)) - frequency_center) frequencies_list.append( np.linspace(frequency_center - width / 2, frequency_center + width / 2, channels_per_line) ) self.channels_per_line.append(channels_per_line) frequencies = np.hstack(frequencies_list) wavelengths = frequencies.to(u.um, equivalencies=u.spectral()).value return wavelengths @u.quantity_input def camera_wavelength_micron( self, out_file: PathLike = None, window_width: u.km / u.s = 6 * u.km / u.s, channels_per_line: int = 200 ) -> None: """Creates a `camera_wavelength_micron.inp` file with all the frequencies for all the lines""" wavelengths = self.wavelength_multiple_lines(channels_per_line, window_width) if out_file is None: out_file = os.path.join(self.folder, 'camera_wavelength_micron.inp') with open(out_file, 'w') as file: print(len(wavelengths), file=file) print('\n'.join(f"{entry:.10e}" for entry in wavelengths), file=file) def create_files(self, window_width: u.km / u.s = 6 * u.km / u.s, channels_per_line: int = 200 ) -> None: super().create_files() self.camera_wavelength_micron(window_width=window_width, channels_per_line=channels_per_line) @u.quantity_input def run( self, inclination: u.deg = 0 * u.deg, position_angle: u.deg = 0 * u.deg, distance: u.pc = 140 * u.pc, velocity_offset: u.km / u.s = 0 * u.km / u.s, threads: int = 1, npix: int = 100 ) -> None: """ Run radmc3d to create the line emission cubes Args: inclination: position_angle: distance: velocity_offset: threads: npix: """ self.logger.info("Running radmc3d") start = time.time() command = (self.command + f"incl {inclination.to(u.deg).value} " f"posang {position_angle.to(u.deg).value} " f"setthreads {threads} " f"npix {npix} " "loadlambda " ) self.logger.info("Running radmc3d for all transition at once: %s", command) proc = subprocess.run( command, cwd=self.folder, text=True, capture_output=True, shell=True ) self.logger.info("radmc3d finished after %s", timedelta(seconds=time.time() - start)) self.catch_radmc_messages(proc) names = [os.path.join(self.folder, f"{line.name}_image.out") for line in self.ordered_line_list] self.split(names=names) for line, name in zip(self.ordered_line_list, names): self.outputs[line] = RadMCOutput(line, file_radmc=name, distance=distance, mode=self.mode) self.outputs[line].file_fits = self.radmc_to_fits(name, line, distance) def channel_maps( self, filename: PathLike = None, extension: typing.Union[None, str] = 'png', distance: u.pc = None, **kwargs ) -> typing.List[Figure]: """ Create channel maps for each line or just a given fits file This method is in unstable state and will be modified and fixed in future versions Args: filename: input fits file to create channel maps from extension: extension for the output figure, e.g. png, pdf, etc. Does not create any file if None distance (u.pc): distance to the source. If None, the distance will be searched with Simbad Returns: list of matplotlib figures containing the channel maps """ if filename is None: outputs = self.outputs else: line = DummyLine(name=filename, molecule=filename) outputs: typing.Dict[DummyLine, RadMCOutput] = { line: RadMCOutput(line=line, file_fits=filename, distance=distance) } figures = [] for output in outputs.values(): fig = output.plot_channel_map(**kwargs) if extension is not None: fig.savefig(os.path.join(self.folder, f"{output.line.name}.{extension}")) figures.append(fig) return figures def split( self, filename=None, names: typing.List[PathLike] = None ): """Split the `image.out` file with all lines into multiple single-frequency files""" lengths = self.channels_per_line if filename is None: filename = os.path.join(self.folder, 'image.out') with open(filename, 'r') as bigfile: files = [open(name, 'w') for name in names] header = [next(bigfile) for _ in range(4)] px_number = [int(entry) for entry in header[1].split()] self.logger.info("Splitting %s in %s", filename, names) for file, length in zip(files, lengths): file.write(header[0]) file.write(header[1]) file.write(str(length) + "\n") file.write(header[3]) for line in range(length): file.write(next(bigfile)) for file, length in zip(files, lengths): for line in range(length * (px_number[0] * px_number[1] + 1)): file.write(next(bigfile)) for file in files: file.close()
Ancestors
Subclasses
Class variables
var velocity : Unit("km / s")
Methods
def camera_wavelength_micron(self, out_file: Union[str, os.PathLike] = None, window_width: Unit("km / s") = <Quantity 6. km / s>, channels_per_line: int = 200) ‑> None
-
Creates a
camera_wavelength_micron.inp
file with all the frequencies for all the linesExpand source code
@u.quantity_input def camera_wavelength_micron( self, out_file: PathLike = None, window_width: u.km / u.s = 6 * u.km / u.s, channels_per_line: int = 200 ) -> None: """Creates a `camera_wavelength_micron.inp` file with all the frequencies for all the lines""" wavelengths = self.wavelength_multiple_lines(channels_per_line, window_width) if out_file is None: out_file = os.path.join(self.folder, 'camera_wavelength_micron.inp') with open(out_file, 'w') as file: print(len(wavelengths), file=file) print('\n'.join(f"{entry:.10e}" for entry in wavelengths), file=file)
def channel_maps(self, filename: Union[str, os.PathLike] = None, extension: Optional[None] = 'png', distance: Unit("pc") = None, **kwargs) ‑> List[matplotlib.figure.Figure]
-
Create channel maps for each line or just a given fits file
This method is in unstable state and will be modified and fixed in future versions
Args
filename
- input fits file to create channel maps from
extension
- extension for the output figure, e.g. png, pdf, etc. Does not create any file if None
distance
:u.pc
- distance to the source. If None, the distance will be searched with Simbad
Returns
list of matplotlib figures containing the channel maps
Expand source code
def channel_maps( self, filename: PathLike = None, extension: typing.Union[None, str] = 'png', distance: u.pc = None, **kwargs ) -> typing.List[Figure]: """ Create channel maps for each line or just a given fits file This method is in unstable state and will be modified and fixed in future versions Args: filename: input fits file to create channel maps from extension: extension for the output figure, e.g. png, pdf, etc. Does not create any file if None distance (u.pc): distance to the source. If None, the distance will be searched with Simbad Returns: list of matplotlib figures containing the channel maps """ if filename is None: outputs = self.outputs else: line = DummyLine(name=filename, molecule=filename) outputs: typing.Dict[DummyLine, RadMCOutput] = { line: RadMCOutput(line=line, file_fits=filename, distance=distance) } figures = [] for output in outputs.values(): fig = output.plot_channel_map(**kwargs) if extension is not None: fig.savefig(os.path.join(self.folder, f"{output.line.name}.{extension}")) figures.append(fig) return figures
def frequency_centers(self) ‑> List[astropy.units.quantity.Quantity]
-
Fetches frequencies of lines from
self.line_list
Expand source code
def frequency_centers(self) -> typing.List[u.Quantity]: """Fetches frequencies of lines from `self.line_list`""" for line in self.line_list: line.parse_lamda() self.ordered_line_list = sorted(self.line_list, key=lambda line: line.frequency) frequency_centers = [line.frequency for line in self.ordered_line_list] return frequency_centers
def run(self, inclination: Unit("deg") = <Quantity 0. deg>, position_angle: Unit("deg") = <Quantity 0. deg>, distance: Unit("pc") = <Quantity 140. pc>, velocity_offset: Unit("km / s") = <Quantity 0. km / s>, threads: int = 1, npix: int = 100) ‑> None
-
Run radmc3d to create the line emission cubes
Args
inclination: position_angle: distance: velocity_offset: threads: npix:
Expand source code
@u.quantity_input def run( self, inclination: u.deg = 0 * u.deg, position_angle: u.deg = 0 * u.deg, distance: u.pc = 140 * u.pc, velocity_offset: u.km / u.s = 0 * u.km / u.s, threads: int = 1, npix: int = 100 ) -> None: """ Run radmc3d to create the line emission cubes Args: inclination: position_angle: distance: velocity_offset: threads: npix: """ self.logger.info("Running radmc3d") start = time.time() command = (self.command + f"incl {inclination.to(u.deg).value} " f"posang {position_angle.to(u.deg).value} " f"setthreads {threads} " f"npix {npix} " "loadlambda " ) self.logger.info("Running radmc3d for all transition at once: %s", command) proc = subprocess.run( command, cwd=self.folder, text=True, capture_output=True, shell=True ) self.logger.info("radmc3d finished after %s", timedelta(seconds=time.time() - start)) self.catch_radmc_messages(proc) names = [os.path.join(self.folder, f"{line.name}_image.out") for line in self.ordered_line_list] self.split(names=names) for line, name in zip(self.ordered_line_list, names): self.outputs[line] = RadMCOutput(line, file_radmc=name, distance=distance, mode=self.mode) self.outputs[line].file_fits = self.radmc_to_fits(name, line, distance)
def split(self, filename=None, names: List[Union[str, os.PathLike]] = None)
-
Split the
image.out
file with all lines into multiple single-frequency filesExpand source code
def split( self, filename=None, names: typing.List[PathLike] = None ): """Split the `image.out` file with all lines into multiple single-frequency files""" lengths = self.channels_per_line if filename is None: filename = os.path.join(self.folder, 'image.out') with open(filename, 'r') as bigfile: files = [open(name, 'w') for name in names] header = [next(bigfile) for _ in range(4)] px_number = [int(entry) for entry in header[1].split()] self.logger.info("Splitting %s in %s", filename, names) for file, length in zip(files, lengths): file.write(header[0]) file.write(header[1]) file.write(str(length) + "\n") file.write(header[3]) for line in range(length): file.write(next(bigfile)) for file, length in zip(files, lengths): for line in range(length * (px_number[0] * px_number[1] + 1)): file.write(next(bigfile)) for file in files: file.close()
def wavelength_multiple_lines(self, channels_per_line: int, window_width: Unit("km / s")) ‑> numpy.ndarray
-
Return a sorted array of wavelength for all the lines, in um (unitless)
Expand source code
def wavelength_multiple_lines(self, channels_per_line: int, window_width: u.km / u.s) -> np.ndarray: """ Return a sorted array of wavelength for all the lines, in um (unitless) """ frequencies_list = [] self.channels_per_line = [] for frequency_center in self.frequency_centers(): width = np.abs(window_width.to(u.GHz, equivalencies=u.doppler_radio(frequency_center)) - frequency_center) frequencies_list.append( np.linspace(frequency_center - width / 2, frequency_center + width / 2, channels_per_line) ) self.channels_per_line.append(channels_per_line) frequencies = np.hstack(frequencies_list) wavelengths = frequencies.to(u.um, equivalencies=u.spectral()).value return wavelengths
Inherited members
class RadMCRTLinesTauSurf (chemistry: Union[ChemistryBase, PhysicsBase] = None, line_list: List[Line] = None, executable: Union[str, os.PathLike] = 'radmc3d', folder: Union[str, os.PathLike] = 'radmc', radii_bins: Optional[None] = None, theta_bins: Optional[None] = None, outer_radius: Optional[None] = None, verbosity: int = 0, wavelengths: astropy.units.quantity.Quantity = <Quantity [9.12000000e-02, 1.00185062e-01, 1.10055336e-01, 1.20898034e-01, 1.32808958e-01, 1.45893352e-01, 1.60266826e-01, 1.76056380e-01, 1.93401528e-01, 2.12455528e-01, 2.33386735e-01, 2.56380093e-01, 2.81638767e-01, 3.09385935e-01, 3.39866766e-01, 3.73350582e-01, 4.10133237e-01, 4.50539734e-01, 4.94927096e-01, 5.43687520e-01, 5.97251841e-01, 6.56093341e-01, 7.20731931e-01, 7.91738741e-01, 8.69741171e-01, 9.55428434e-01, 1.04955764e+00, 1.15296050e+00, 1.26655065e+00, 1.39133174e+00, 1.52840633e+00, 1.67898555e+00, 1.84439991e+00, 2.02611096e+00, 2.22572425e+00, 2.44500354e+00, 2.68588630e+00, 2.95050094e+00, 3.24118552e+00, 3.56050847e+00, 3.91129125e+00, 4.29663329e+00, 4.71993939e+00, 5.18494979e+00, 5.69577320e+00, 6.25692315e+00, 6.87335783e+00, 7.55052390e+00, 8.29440466e+00, 9.11157287e+00, 1.00092488e+01, 1.09953642e+01, 1.20786320e+01, 1.32686239e+01, 1.45758543e+01, 1.60118735e+01, 1.75893699e+01, 1.93222820e+01, 2.12259213e+01, 2.33171079e+01, 2.56143191e+01, 2.81378525e+01, 3.09100054e+01, 3.39552720e+01, 3.73005596e+01, 4.09754262e+01, 4.50123423e+01, 4.94469770e+01, 5.43185138e+01, 5.96699964e+01, 6.55487093e+01, 7.20065955e+01, 7.91007152e+01, 8.68937507e+01, 9.54545592e+01, 1.04858782e+02, 1.15189513e+02, 1.26538032e+02, 1.39004611e+02, 1.52699404e+02, 1.67743413e+02, 1.84269563e+02, 2.02423878e+02, 2.22366762e+02, 2.44274428e+02, 2.68340447e+02, 2.94777460e+02, 3.23819058e+02, 3.55721846e+02, 3.90767711e+02, 4.29266308e+02, 4.71557804e+02, 5.18015875e+02, 5.69051015e+02, 6.25114158e+02, 6.86700665e+02, 7.54354701e+02, 8.28674040e+02, 9.10315352e+02, 1.00000000e+03] um>, modified_random_walk: bool = True, scattering_mode_max: int = None, nphot_therm: int = None, nphot_mono: int = None, nphot_spec: int = None, coordinate: Union[str, astropy.coordinates.sky_coordinate.SkyCoord] = None, external_source_type: Literal['Draine1978', 'WeingartnerDraine2001', 'None'] = None, external_source_multiplier: float = 1, velocity: Unit("km / s") = <Quantity 0. km / s>, mode: Literal['image tracetau'] = 'tausurf', tau: float = 1.0)
-
Class to calculate optical depth map in lines
Expand source code
@dataclass class RadMCRTLinesTauSurf(RadMCRTLines): """Class to calculate optical depth map in lines""" allowed_modes = ["tausurf"] mode: typing.Literal["image tracetau"] = "tausurf" tau: float = 1. @property def unit(self) -> u.Unit: """Unit of the returned data""" return u.au @property def command(self) -> str: return super().command + f"{self.tau} " @property def factor(self) -> float: """Factor of the returned data relative to self.unit""" return (u.cm / u.au).si.scale def _extra_header(self, header: Header, **kwargs): """Add extra header items into the output fits files""" super()._extra_header(header, **{"HISTORY": f"Optical depth tau={self.tau: .2e} map"}, **kwargs) @property def fitsname_suffix(self) -> str: """Additional suffix of an output fits file""" return f"_tau_{self.tau: .2e}"
Ancestors
Class variables
var allowed_modes
var tau : float
Instance variables
var command : str
-
Expand source code
@property def command(self) -> str: return super().command + f"{self.tau} "
Inherited members
RadMCRTLines
:amr_grid
camera_wavelength_micron
catch_radmc_messages
channel_maps
copy_for_propype
create_files
external_source
factor
fitsname_suffix
frequency_centers
gas_temperature
interpolate
interpolate_back
molecule
numberdens
radmc3d
radmc_to_fits
run
split
unit
wavelength_micron
wavelength_multiple_lines
RadMCBase
:
class RadMCRTLinesTraceTau (chemistry: Union[ChemistryBase, PhysicsBase] = None, line_list: List[Line] = None, executable: Union[str, os.PathLike] = 'radmc3d', folder: Union[str, os.PathLike] = 'radmc', radii_bins: Optional[None] = None, theta_bins: Optional[None] = None, outer_radius: Optional[None] = None, verbosity: int = 0, wavelengths: astropy.units.quantity.Quantity = <Quantity [9.12000000e-02, 1.00185062e-01, 1.10055336e-01, 1.20898034e-01, 1.32808958e-01, 1.45893352e-01, 1.60266826e-01, 1.76056380e-01, 1.93401528e-01, 2.12455528e-01, 2.33386735e-01, 2.56380093e-01, 2.81638767e-01, 3.09385935e-01, 3.39866766e-01, 3.73350582e-01, 4.10133237e-01, 4.50539734e-01, 4.94927096e-01, 5.43687520e-01, 5.97251841e-01, 6.56093341e-01, 7.20731931e-01, 7.91738741e-01, 8.69741171e-01, 9.55428434e-01, 1.04955764e+00, 1.15296050e+00, 1.26655065e+00, 1.39133174e+00, 1.52840633e+00, 1.67898555e+00, 1.84439991e+00, 2.02611096e+00, 2.22572425e+00, 2.44500354e+00, 2.68588630e+00, 2.95050094e+00, 3.24118552e+00, 3.56050847e+00, 3.91129125e+00, 4.29663329e+00, 4.71993939e+00, 5.18494979e+00, 5.69577320e+00, 6.25692315e+00, 6.87335783e+00, 7.55052390e+00, 8.29440466e+00, 9.11157287e+00, 1.00092488e+01, 1.09953642e+01, 1.20786320e+01, 1.32686239e+01, 1.45758543e+01, 1.60118735e+01, 1.75893699e+01, 1.93222820e+01, 2.12259213e+01, 2.33171079e+01, 2.56143191e+01, 2.81378525e+01, 3.09100054e+01, 3.39552720e+01, 3.73005596e+01, 4.09754262e+01, 4.50123423e+01, 4.94469770e+01, 5.43185138e+01, 5.96699964e+01, 6.55487093e+01, 7.20065955e+01, 7.91007152e+01, 8.68937507e+01, 9.54545592e+01, 1.04858782e+02, 1.15189513e+02, 1.26538032e+02, 1.39004611e+02, 1.52699404e+02, 1.67743413e+02, 1.84269563e+02, 2.02423878e+02, 2.22366762e+02, 2.44274428e+02, 2.68340447e+02, 2.94777460e+02, 3.23819058e+02, 3.55721846e+02, 3.90767711e+02, 4.29266308e+02, 4.71557804e+02, 5.18015875e+02, 5.69051015e+02, 6.25114158e+02, 6.86700665e+02, 7.54354701e+02, 8.28674040e+02, 9.10315352e+02, 1.00000000e+03] um>, modified_random_walk: bool = True, scattering_mode_max: int = None, nphot_therm: int = None, nphot_mono: int = None, nphot_spec: int = None, coordinate: Union[str, astropy.coordinates.sky_coordinate.SkyCoord] = None, external_source_type: Literal['Draine1978', 'WeingartnerDraine2001', 'None'] = None, external_source_multiplier: float = 1, velocity: Unit("km / s") = <Quantity 0. km / s>, mode: Literal['image tracetau'] = 'image tracetau')
-
Class to calculate optical depth map in lines
Expand source code
@dataclass class RadMCRTLinesTraceTau(RadMCRTLines): """Class to calculate optical depth map in lines""" allowed_modes = ["image tracetau"] mode: typing.Literal["image tracetau"] = "image tracetau" @property def unit(self) -> u.Unit: """Unit of the returned data""" return u.dimensionless_unscaled @property def factor(self) -> float: """Factor of the returned data relative to self.unit""" return 1. def _extra_header(self, header: Header, **kwargs): """Add extra header items into the output fits files""" super()._extra_header(header, **{"HISTORY": "Optical depth map"}, **kwargs) @property def fitsname_suffix(self) -> str: """Additional suffix of an output fits file""" return "_tracetau"
Ancestors
Class variables
var allowed_modes
Inherited members
RadMCRTLines
:amr_grid
camera_wavelength_micron
catch_radmc_messages
channel_maps
copy_for_propype
create_files
external_source
factor
fitsname_suffix
frequency_centers
gas_temperature
interpolate
interpolate_back
molecule
numberdens
radmc3d
radmc_to_fits
run
split
unit
wavelength_micron
wavelength_multiple_lines
RadMCBase
:
class RadMCRTSingleCall (chemistry: Union[ChemistryBase, PhysicsBase] = None, line_list: List[Line] = None, executable: Union[str, os.PathLike] = 'radmc3d', folder: Union[str, os.PathLike] = 'radmc', radii_bins: Optional[None] = None, theta_bins: Optional[None] = None, outer_radius: Optional[None] = None, verbosity: int = 0, wavelengths: astropy.units.quantity.Quantity = <Quantity [9.12000000e-02, 1.00185062e-01, 1.10055336e-01, 1.20898034e-01, 1.32808958e-01, 1.45893352e-01, 1.60266826e-01, 1.76056380e-01, 1.93401528e-01, 2.12455528e-01, 2.33386735e-01, 2.56380093e-01, 2.81638767e-01, 3.09385935e-01, 3.39866766e-01, 3.73350582e-01, 4.10133237e-01, 4.50539734e-01, 4.94927096e-01, 5.43687520e-01, 5.97251841e-01, 6.56093341e-01, 7.20731931e-01, 7.91738741e-01, 8.69741171e-01, 9.55428434e-01, 1.04955764e+00, 1.15296050e+00, 1.26655065e+00, 1.39133174e+00, 1.52840633e+00, 1.67898555e+00, 1.84439991e+00, 2.02611096e+00, 2.22572425e+00, 2.44500354e+00, 2.68588630e+00, 2.95050094e+00, 3.24118552e+00, 3.56050847e+00, 3.91129125e+00, 4.29663329e+00, 4.71993939e+00, 5.18494979e+00, 5.69577320e+00, 6.25692315e+00, 6.87335783e+00, 7.55052390e+00, 8.29440466e+00, 9.11157287e+00, 1.00092488e+01, 1.09953642e+01, 1.20786320e+01, 1.32686239e+01, 1.45758543e+01, 1.60118735e+01, 1.75893699e+01, 1.93222820e+01, 2.12259213e+01, 2.33171079e+01, 2.56143191e+01, 2.81378525e+01, 3.09100054e+01, 3.39552720e+01, 3.73005596e+01, 4.09754262e+01, 4.50123423e+01, 4.94469770e+01, 5.43185138e+01, 5.96699964e+01, 6.55487093e+01, 7.20065955e+01, 7.91007152e+01, 8.68937507e+01, 9.54545592e+01, 1.04858782e+02, 1.15189513e+02, 1.26538032e+02, 1.39004611e+02, 1.52699404e+02, 1.67743413e+02, 1.84269563e+02, 2.02423878e+02, 2.22366762e+02, 2.44274428e+02, 2.68340447e+02, 2.94777460e+02, 3.23819058e+02, 3.55721846e+02, 3.90767711e+02, 4.29266308e+02, 4.71557804e+02, 5.18015875e+02, 5.69051015e+02, 6.25114158e+02, 6.86700665e+02, 7.54354701e+02, 8.28674040e+02, 9.10315352e+02, 1.00000000e+03] um>, modified_random_walk: bool = True, scattering_mode_max: int = None, nphot_therm: int = None, nphot_mono: int = None, nphot_spec: int = None, coordinate: Union[str, astropy.coordinates.sky_coordinate.SkyCoord] = None, external_source_type: Literal['Draine1978', 'WeingartnerDraine2001', 'None'] = None, external_source_multiplier: float = 1, velocity: Unit("km / s") = <Quantity 0. km / s>, mode: Literal['image'] = 'image')
-
Deprecated in favor of RadMCRTLines
Expand source code
@dataclass class RadMCRTSingleCall(RadMCRTLines): """Deprecated in favor of RadMCRTLines""" def __post_init__(self): super().__post_init__() self.logger.warning("Deprecation warning. Renamed to RadMCRTLines!")
Ancestors
Class variables
var velocity : Unit("km / s")
Inherited members
RadMCRTLines
:amr_grid
camera_wavelength_micron
catch_radmc_messages
channel_maps
copy_for_propype
create_files
external_source
factor
fitsname_suffix
frequency_centers
gas_temperature
interpolate
interpolate_back
molecule
numberdens
radmc3d
radmc_to_fits
run
split
unit
wavelength_micron
wavelength_multiple_lines
RadMCBase
: