Star.py

import galsim
import os, sys
import numpy as np
import astropy.constants as cons
from astropy.table import Table
from scipy import interpolate

from ObservationSim.MockObject._util import integrate_sed_bandpass, getNormFactorForSpecWithABMAG, getObservedSED, getABMAG, tag_sed
from ObservationSim.MockObject.MockObject import MockObject

class Star(MockObject):
    def __init__(self, param, logger=None):
        super().__init__(param, logger=logger)

    def unload_SED(self):
        """(Test) free up SED memory
        """
        del self.sed

    def getGSObj(self, psf, g1=0, g2=0, flux=None, filt=None, tel=None, exptime=150.):
        if flux == None:
            flux = self.getElectronFluxFilt(filt, tel, exptime)
        # star = galsim.Gaussian(sigma=1.e-8, flux=1.)
        star = galsim.DeltaFunction()
        star = star.withFlux(flux)
        final = galsim.Convolve(psf, star)
        return final
        
    def getGSObj_multiband(self, tel, psf_list, bandpass_list, filt, nphotons_tot=None, g1=0, g2=0, exptime=150.):
        if len(psf_list) != len(bandpass_list):
            raise ValueError("!!!The number of PSF profiles and the number of bandpasses must be equal.")
        objs = []
        if nphotons_tot == None:
            nphotons_tot = self.getElectronFluxFilt(filt, tel, exptime)

        try:
            full = integrate_sed_bandpass(sed=self.sed, bandpass=filt.bandpass_full)
        except Exception as e:
            print(e)
            self.logger.error(e)
            return -1

        for i in range(len(bandpass_list)):
            bandpass = bandpass_list[i]
            psf = psf_list[i]
            try:
                sub = integrate_sed_bandpass(sed=self.sed, bandpass=bandpass)
            except Exception as e:
                print(e)
                self.logger.error(e)
                return -1
        
            ratio = sub/full

            if not (ratio == -1 or (ratio != ratio)):
                nphotons = ratio * nphotons_tot
            else:
                return -1
            star = galsim.DeltaFunction()
            star = star.withFlux(nphotons)
            star = galsim.Convolve(psf, star)
            objs.append(star)
        final = galsim.Sum(objs)
        return final