SUGAR: a SED model of SNe Ia based on spectral features

The SUGAR model is available here: sugar_template_v1.zip.

Below is an example of how to use the SUGAR template to obtain the SUGAR spectral time series for a given set of SUGAR parameters (unzip the file in order to be able to execute this code).

import numpy as np
import copy
import os
import pylab as plt

NBIN = 196
NPHASE = 21
PHASEMIN = -12
PHASEMAX = 48

# function to convert in flux
def go_to_flux(X, Y, ABmag0=48.59):
    """Convert AB mag to flux."""
    Flux_nu = 10**(-0.4 * (Y + ABmag0))
    f = X**2 / 299792458. * 1.e-10
    Flux_lambda = Flux_nu / f
    return Flux_lambda

class sugar_model(object):

        def __init__(self, rep_sugar='sugar_template_v1/'):

            arange = np.arange(NBIN*NPHASE).reshape(NPHASE, NBIN).T.reshape(-1)

            m0 = np.loadtxt(os.path.join(rep_sugar, 'sugar_template_0.dat'))
            self.alpha_1 = np.loadtxt(os.path.join(rep_sugar, 'sugar_template_1.dat'))[:,2][arange]
            self.alpha_2 = np.loadtxt(os.path.join(rep_sugar, 'sugar_template_2.dat'))[:,2][arange]
            self.alpha_3 = np.loadtxt(os.path.join(rep_sugar, 'sugar_template_3.dat'))[:,2][arange]
            self.f = np.loadtxt(os.path.join(rep_sugar, 'sugar_template_4.dat'))[:,2][arange]
            self.epoch = m0[:,0][arange]
            self.wavelength = m0[:,1][arange]
            self.m0 =m0[:,2][arange]

            self.Y = None
            self.q1 = None
            self.q2 = None
            self.q3 = None
            self.av = None
            self.grey = None

        def sugar_in_mag(self, q1, q2, q3, av, grey):

            self.q1 = q1
            self.q2 = q2
            self.q3 = q3
            self.av = av
            self.grey = grey

            # init with average template
            self.Y = copy.deepcopy(self.m0)
            # add intrinsic features (stretch, velocity, and second-strech/calcim related)
            self.Y += self.q1 * self.alpha_1
            self.Y += self.q2 * self.alpha_2
            self.Y += self.q3 * self.alpha_3
            # add color
            self.Y += self.av * self.f
            # add grey offset
            self.Y += self.grey

            arange = np.arange(NBIN*NPHASE).reshape(NBIN, NPHASE).T.reshape(-1)
            self.Y = self.Y[arange]
            self.X = self.wavelength[arange]

        def plot_sed_in_mag(self):

            cmap = plt.cm.seismic
            col = cmap(np.linspace(0,1,NPHASE))
            fig = plt.figure(figsize=(8,6))
            sugar_params = '$q_1$=%.1f, $q_2$=%.1f, $q_3$=%.1f, $A_V$=%.1f, $\Delta M_{grey}$=%.1f'%((self.q1,
                                                                                                      self.q2,
                                                                                                      self.q3,
                                                                                                      self.av,
                                                                                                      self.grey))
            plt.title('SUGAR template \n '+ sugar_params, fontsize=14)
            for i in range(NPHASE):
                plt.plot(self.X[i*NBIN:(i+1)*NBIN], self.Y[i*NBIN:(i+1)*NBIN], c=col[i])
            plt.ylabel('AB mag + const.', fontsize=16)
            plt.xlabel('wavelength $[\AA]$', fontsize=16)
            plt.gca().invert_yaxis()

            cbar_ax = fig.add_axes([0.93, 0.1, 0.02, 0.8])
            min_max = plt.mpl.colors.Normalize(vmin=PHASEMIN, vmax=PHASEMAX)
            cbar = plt.mpl.colorbar.ColorbarBase(cbar_ax, cmap=cmap, norm=min_max)
            cbar.set_label('Epoch relative to SUGAR t0 (days)', fontsize=14)

        def plot_sed_in_flux(self):

            cmap = plt.cm.seismic
            col = cmap(np.linspace(0,1,NPHASE))
            fig = plt.figure(figsize=(8,6))
            sugar_params = '$q_1$=%.1f, $q_2$=%.1f, $q_3$=%.1f, $A_V$=%.1f, $\Delta M_{grey}$=%.1f'%((self.q1,
                                                                                                      self.q2,
                                                                                                      self.q3,
                                                                                                      self.av,
                                                                                                      self.grey))
            plt.title('SUGAR template \n '+ sugar_params, fontsize=14)
            for i in range(NPHASE):
                y_flux = go_to_flux(self.X[i*NBIN:(i+1)*NBIN], self.Y[i*NBIN:(i+1)*NBIN])
                plt.plot(self.X[i*NBIN:(i+1)*NBIN], y_flux, c=col[i])
            plt.ylabel('Flux $\\times$ const.', fontsize=16)
            plt.xlabel('wavelength $[\AA]$', fontsize=16)

            cbar_ax = fig.add_axes([0.93, 0.1, 0.02, 0.8])
            min_max = plt.mpl.colors.Normalize(vmin=PHASEMIN, vmax=PHASEMAX)
            cbar = plt.mpl.colorbar.ColorbarBase(cbar_ax, cmap=cmap, norm=min_max)
            cbar.set_label('Epoch relative to SUGAR t0 (days)', fontsize=14)

sugar = sugar_model()
sugar.sugar_in_mag(-2., 1., 0.5, 0.1, 0)
sugar.plot_sed_in_mag()
sugar.plot_sed_in_flux()
plt.show()