Star.py 2.36 KB
Newer Older
Fang Yuedong's avatar
Fang Yuedong committed
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
import galsim
import os
import sys
import numpy as np
import astropy.constants as cons
from astropy.table import Table
from scipy import interpolate

from observation_sim.mock_objects._util import integrate_sed_bandpass, getNormFactorForSpecWithABMAG, getObservedSED, getABMAG, tag_sed
from observation_sim.mock_objects.MockObject import MockObject


class Star(MockObject):
    def __init__(self, param, logger=None):
        super().__init__(param, logger=logger)
        if not hasattr(self, "mu"):
            self.mu = 1.

    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