start adding AGN-fitting mode

desihub · Dec 6, 2024 · 5e93273 · 5e93273
1 parent f1289c1
commit 5e93273
Show file tree

Hide file tree

Showing 3 changed files with 30 additions and 23 deletions.
diff --git a/py/fastspecfit/fastspecfit.py b/py/fastspecfit/fastspecfit.py
@@ -150,6 +150,7 @@ def fastspec(fastphot=False, fitstack=False, args=None, comm=None, verbose=False
         'fphotofile':        args.fphotofile,
         'fastphot':          fastphot,
         'fitstack':          fitstack,
+        'agn_continuum':     args.agn_continuum,
         'ignore_photometry': args.ignore_photometry,
         'template_file':     args.templates,
         'template_version':  args.templateversion,
@@ -335,6 +336,7 @@ def parse(options=None):
     parser.add_argument('--seed', type=int, default=1, help='Random seed for Monte Carlo reproducibility; ignored if --input-seeds is passed.')
     parser.add_argument('--nmonte', type=int, default=50, help='Number of Monte Carlo realizations.')
     parser.add_argument('--zmin', type=float, default=None, help='Override the default minimum redshift required for modeling.')
+    parser.add_argument('--agn-continuum', default=False, action='store_true', help='Fit AGN, not stellar continuum templates.')
     parser.add_argument('--no-broadlinefit', default=True, action='store_false', dest='broadlinefit',
                         help='Do not model broad Balmer and helium line-emission.')
     parser.add_argument('--ignore-photometry', default=False, action='store_true', help='Ignore the broadband photometry during model fitting.')

diff --git a/py/fastspecfit/singlecopy.py b/py/fastspecfit/singlecopy.py
@@ -22,6 +22,7 @@ def initialize(self,
                    fphotofile=None,
                    fastphot=False,
                    fitstack=False,
+                   agn_continuum=False,
                    ignore_photometry=False,
                    template_file=None,
                    template_version=None,
@@ -39,8 +40,9 @@ def initialize(self,
                                    imf=template_imf,
                                    mintemplatewave=None,
                                    maxtemplatewave=40e4,
+                                   agn_continuum=agn_continuum,
                                    fastphot=fastphot)
-        log.debug(f'Cached stellar templates {self.templates.file}')
+        log.debug(f'Cached templates {self.templates.file}')
 
         # emission line table
         self.emlines = LineTable(emlines_file)

diff --git a/py/fastspecfit/templates.py b/py/fastspecfit/templates.py
@@ -32,7 +32,7 @@ class Templates(object):
 
     def __init__(self, template_file=None, template_version=None, imf=None,
                  mintemplatewave=None, maxtemplatewave=40e4, vdisp_nominal=250.,
-                 fastphot=False, read_linefluxes=False):
+                 agn_continuum=False, fastphot=False, read_linefluxes=False):
         """"
         Read the templates into a dictionary.
 
@@ -123,27 +123,30 @@ def __init__(self, template_file=None, template_version=None, imf=None,
             #self.gamma    = dusthdr['GAMMA']
 
             # construct the AGN wavelength vector
-            iragnflux = T['AGNFLUX'].read()
-            iragnwave = T['AGNWAVE'].read()
-            #iragnflux  /= trapz(iragnflux, x=iragnwave) # should already be 1.0
-
-            trim = np.searchsorted(iragnwave, 1e4, 'left') # hack...
-            iragnflux = iragnflux[trim:]
-            iragnwave = iragnwave[trim:]
-
-            feflux = T['FEFLUX'].read()
-            fewave = T['FEWAVE'].read()
-
-            febounds = np.searchsorted(templatewave, Templates.AGN_PIXKMS_BOUNDS, 'left')
-            irbounds = np.searchsorted(templatewave, iragnwave[0], 'left')
-
-            agnwave = np.hstack((templatewave[:febounds[0]], fewave,
-                                 templatewave[febounds[1]:irbounds],
-                                 iragnwave))
-            #self.agnwave = agnwave
-
-            #agnhdr = T['AGNFLUX'].read_header()
-            #self.agntau   = agnhdr['AGNTAU']
+            if agn_continuum:
+                iragnflux = T['AGNFLUX'].read()
+                iragnwave = T['AGNWAVE'].read()
+                #iragnflux  /= trapz(iragnflux, x=iragnwave) # should already be 1.0
+
+                trim = np.searchsorted(iragnwave, 1e4, 'left') # hack...
+                iragnflux = iragnflux[trim:]
+                iragnwave = iragnwave[trim:]
+
+                feflux = T['FEFLUX'].read()
+                fewave = T['FEWAVE'].read()
+                self.fewave = feflux
+                self.feflux = fewave
+
+                febounds = np.searchsorted(templatewave, Templates.AGN_PIXKMS_BOUNDS, 'left')
+                irbounds = np.searchsorted(templatewave, iragnwave[0], 'left')
+
+                agnwave = np.hstack((templatewave[:febounds[0]], fewave,
+                                     templatewave[febounds[1]:irbounds],
+                                     iragnwave))
+                self.agnwave = agnwave
+
+                #agnhdr = T['AGNFLUX'].read_header()
+                #self.agntau   = agnhdr['AGNTAU']
         else:
             errmsg = f'Templates file {template_file} missing mandatory extensions DUSTFLUX and AGNFLUX.'
             log.critical(errmsg)