Skip to content

Commit

Permalink
[TEP005][QOL] Add MontecarloRunner.spectrum_integrated
Browse files Browse the repository at this point in the history
This is a shortcut for obtaining an integrated spectrum for the
wavelength grid defined in the configuration file. If incompatible
configuration settings are detected, an IntegrationError is raised.
  • Loading branch information
yeganer committed Jun 21, 2017
1 parent 55e0696 commit aeb3ac4
Show file tree
Hide file tree
Showing 3 changed files with 58 additions and 21 deletions.
27 changes: 7 additions & 20 deletions scripts/tardis
Original file line number Diff line number Diff line change
Expand Up @@ -67,28 +67,13 @@ tardis_config = config_reader.Configuration.from_yaml(args.config_fname)
simulation = Simulation.from_config(tardis_config)


def get_integrated_spectrum():
if (
tardis_config.plasma.disable_electron_scattering and
tardis_config.plasma.line_interaction_type == 'downbranch'
):
return simulation.runner.integrator.calculate_spectrum(
simulation.runner.spectrum_frequency[:-1],
1000)
else:
raise config_reader.ConfigurationError(
'The formal integral module currently works only for '
'plasma settings "downbranch" and disabled electron scattering'
'.\nWe got: {}'.format(str(tardis_config.plasma)))


def get_virtual_spectrum():
# Catch warning when acessing invalid spectrum_virtual
with warnings.catch_warnings(record=True) as w:
spectrum = simulation.runner.spectrum_virtual

if len(w) > 0 and w[-1]._category_name == 'UserWarning':
print(
warnings.warn(
'Virtual spectrum is not available, using the '
'real packet spectrum instead.')
spectrum = simulation.runner.spectrum
Expand All @@ -99,20 +84,22 @@ print('Saving the {} spectrum.'.format(tardis_config.spectrum.method))
# Believe it or not, that is a fancy switch-case statement
# We need functions so that only one spectrum is accessed so we can catch
# the warning properly
spectrum = {
get_spectrum = {
'real': lambda: simulation.runner.spectrum,
'virtual': get_virtual_spectrum,
'integrated': get_integrated_spectrum,
'integrated': lambda: simulation.runner.spectrum_integrated,
}[tardis_config.spectrum.method]


if args.gdb:
import os; print(os.getpid()); raw_input() # Workaround to attach gdb
import os
print(os.getpid())
raw_input() # Workaround to attach gdb
if args.profile:
import cProfile
cProfile.runctx('simulation.run()', locals(),
globals(), filename=args.profiler_log_file)
else:
simulation.run()

spectrum().to_ascii(args.spectrum_fname)
get_spectrum().to_ascii(args.spectrum_fname)
7 changes: 7 additions & 0 deletions tardis/montecarlo/base.py
Original file line number Diff line number Diff line change
Expand Up @@ -124,6 +124,13 @@ def spectrum_virtual(self):
self.spectrum_frequency,
self.montecarlo_virtual_luminosity)

@property
def spectrum_integrated(self):
return self.integrator.calculate_spectrum(
self.spectrum_frequency[:-1],
1000,
raises=True)

@property
def integrator(self):
if self._integrator is None:
Expand Down
45 changes: 44 additions & 1 deletion tardis/montecarlo/formal_integral.py
Original file line number Diff line number Diff line change
@@ -1,3 +1,4 @@
import warnings
import numpy as np
import pandas as pd
from astropy import units as u
Expand All @@ -6,17 +7,59 @@
from tardis.montecarlo.spectrum import TARDISSpectrum


class IntegrationError(Exception):
pass

class FormalIntegrator(object):

def __init__(self, model, plasma, runner):
self.model = model
self.plasma = plasma
self.runner = runner

def calculate_spectrum(self, frequency, N=1000):
def check(self, raises=False):
'''
A method that determines if the formal integral can be performed with
the current configuration settings
The function returns False if the configuration conflicts with the
required settings. If raises evaluates to True, then a
IntegrationError is raised instead
'''
def raise_or_return(message):
if raises:
raise IntegrationError(message)
else:
warnings.warn(message)
return False

for obj in (self.model, self.plasma, self.runner):
if obj is None:
return raise_or_return(
'The integrator is missing either model, plasma or '
'runner. Please make sure these are provided to the '
'FormalIntegrator.'
)

if not self.runner.line_interaction_type == 'downbranch':
return raise_or_return(
'The FormalIntegrator currently only works for '
'line_interaction_type == "downbranch"'
)

if self.runner.sigma_thomson.value > 1e-100:
return raise_or_return(
'The FormalIntegrator currently only works with '
'disable_electron_scattering turned on.'
)

return True

def calculate_spectrum(self, frequency, N=1000, raises=False):
# Very crude implementation
# The c extension needs bin centers (or something similar)
# while TARDISSpectrum needs bin edges
self.check(raises)
frequency = frequency.to('Hz', u.spectral())
luminosity = u.Quantity(
formal_integral(
Expand Down

0 comments on commit aeb3ac4

Please sign in to comment.