Add Notebook

SL_AGN/v1.2/lib/inject.py

@@ -96,6 +96,26 @@ def calexp_inject_stamp(calexp, inj_catalog):
 
     return injected_exposure
 
+#--------------------------------------
+def visit_inject_stamp(visit, inj_catalog):
+
+    psf = visit.getPsf()
+    photo_calib = visit.getPhotoCalib()
+    wcs = visit.getWcs()
+
+    inject_config = VisitInjectConfig()
+    inject_task = VisitInjectTask(config=inject_config)
+
+    injected_output = inject_task.run(
+        injection_catalogs=inj_catalog,
+        input_exposure=visit.clone(),
+        psf=psf,
+        photo_calib=photo_calib,
+        wcs=wcs,
+    )
+    injected_exposure = injected_output.output_exposure
+
+    return injected_exposure
 
 #--------------------------------------
 def template_inject_stamp(template, inj_catalog):
@@ -117,3 +137,31 @@ def template_inject_stamp(template, inj_catalog):
     injected_exposure = injected_output.output_exposure
 
     return injected_exposure
+
+#-------------------
+# If you want to generate a catalog of stars, simply specify source_type = 'star' and provide magnitude
+'''
+def get_visit_catalogs(num_visits, num_sources=10, source_type, ra=37.93, dec=6.93, extension=0.03, mag_list=[10], mag_min=10, mag_max=15, filename="", n=[2],q=[0.9], beta=[31.0], hlr=[5.0]):
+    source_list = []
+    for i in range(num_sources):
+            source_list.append((np.random.uniform(ra-extension, ra+extension), np.random.uniform(dec-extension, dec+extension), np.random.uniform(10,15, size=num_visits)))
+    inj_catalogs = []
+    for i in range(num_visits):
+        inj_catalogs.append(Table())
+    if(source_type.lower() == 'star' or source_type.lower() == 'sersic'):
+        for time in range(num_times):
+            for source in range(num_sources):
+                inj_catalogs[time] = vstack([inj_catalogs[time], Table(
+                    {'injection_id': [source],
+                        'ra': [source_list[source][0]],
+                        'dec': [source_list[source][1]],
+                        'source_type': ['Star'],
+                        'mag': [source_list[source][2][time]],
+                        'n': n,
+                        'q': q,
+                        'beta': beta,
+                        'half_light_radius': hlr,
+                    }
+                )])
+
+                '''

SL_AGN/v1.2/lib/stamp.py

@@ -2,6 +2,9 @@
 import pandas as pd
 import numpy as np
 from scipy.ndimage import gaussian_filter
+import lsst.afw.math as afwMath
+import lsst.afw.geom as afwGeom
+import lsst.geom as geom
 import os
 from astropy.table import Table
 from astropy.io import fits
@@ -10,8 +13,7 @@
 matplotlib.use("Agg")
 import matplotlib.pyplot as plt
 from lib.tools import *
-
-
+import lsst.afw.image as afwImage
 
 #============================
 MAG0 = 28.17
@@ -244,3 +246,26 @@ def add_wcs(filename):
     fits.writeto("%s_wcs.fits"%tag, data, hdr_new, overwrite=True)
 
     return 0
+
+def rotate_exposure(exp, n_degrees):
+    n_degrees = n_degrees % 360
+
+    wcs = exp.getWcs()
+
+    warper = afwMath.Warper('lanczos4')
+
+    affine_rot_transform = geom.AffineTransform.makeRotation(n_degrees*geom.degrees)
+    transform_p2top2 = afwGeom.makeTransform(affine_rot_transform)
+    rotated_wcs = afwGeom.makeModifiedWcs(transform_p2top2, wcs, False)
+
+    rotated_exp = warper.warpExposure(rotated_wcs, exp)
+    return rotated_exp
+
+def make_rotated_stamp(visit, stamp_filename):
+    add_wcs(stamp_filename)
+    wcs_stamp_filename = "%s_wcs.fits"%stamp_filename[:-5]
+    stamp_img_orig = afwImage.ExposureF.readFits(wcs_stamp_filename)
+    stamp_img_rotated = rotate_exposure(stamp_img_orig, -1*visit.visitInfo.getBoresightRotAngle().asDegrees())
+    stamp_img_rotated.image.array[np.where(np.isnan(stamp_img_rotated.image.array))] = 0.0
+    rot_stamp_filename = wcs_stamp_filename[:4]+'rotated_'+wcs_stamp_filename[4:]
+    stamp_img_rotated.writeFits(rot_stamp_filename)