Skip to content

Commit

Permalink
Merge pull request #269 from nkhadka21/LensClass_refactor
Browse files Browse the repository at this point in the history
Lens class refactor
  • Loading branch information
nkhadka21 authored Nov 7, 2024
2 parents ca033c9 + a819578 commit 5cc87c5
Show file tree
Hide file tree
Showing 23 changed files with 1,375 additions and 703 deletions.
3 changes: 2 additions & 1 deletion codecov.yml
Original file line number Diff line number Diff line change
Expand Up @@ -7,6 +7,7 @@ comment: # this is a top-level key
require_base: false # [true :: must have a base report to post]
require_head: true # [true :: must have a head report to post]
ignore:
- "slsim/lsst_science_pipeline.py"
- "slsim/LsstSciencePipeline/lsst_science_pipeline.py"
- "slsim/LsstSciencePipeline/opsim_pipeline.py"
- "setup.py"
- "tests/*"
Original file line number Diff line number Diff line change
Expand Up @@ -26,7 +26,8 @@
"from slsim.image_simulation import lens_image_series, sharp_rgb_image\n",
"from slsim.Plots.plot_functions import create_image_montage_from_image_list\n",
"from slsim.image_simulation import point_source_coordinate_properties\n",
"from slsim.LsstSciencePipeline import opsim_pipeline"
"from slsim.LsstSciencePipeline import opsim_pipeline\n",
"from slsim.LsstSciencePipeline.util_lsst import opsim_variable_lens_injection"
]
},
{
Expand Down Expand Up @@ -477,7 +478,7 @@
"num_pix = 200\n",
"transform_pix2angle = np.array([[0.2, 0], [0, 0.2]])\n",
"\n",
"images = opsim_pipeline.opsim_variable_lens_injection(\n",
"images = opsim_variable_lens_injection(\n",
" lens_class, bands, num_pix, transform_pix2angle, exposure_data[index]\n",
")\n",
"print(\"images.keys() : \", images.keys())\n",
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -18,10 +18,8 @@
"source": [
"import numpy as np\n",
"from slsim.Util.param_util import random_radec_string\n",
"from slsim.lsst_science_pipeline import (\n",
" multiple_variable_lens_injection,\n",
" dp0_time_series_images_data,\n",
")\n",
"from slsim.lsst_science_pipeline import dp0_time_series_images_data\n",
"from slsim.LsstSciencePipeline.util_lsst import multiple_variable_lens_injection\n",
"import lsst.daf.butler as dafButler\n",
"from astropy.cosmology import FlatLambdaCDM\n",
"from astropy.units import Quantity\n",
Expand Down
2 changes: 1 addition & 1 deletion requirements.txt
Original file line number Diff line number Diff line change
@@ -1,5 +1,5 @@
skypy
lenstronomy>=1.11.10
lenstronomy>=1.12.1
astropy>=5.2
numpy
scipy<=1.13.1
Expand Down
2 changes: 1 addition & 1 deletion slsim/Deflectors/elliptical_lens_galaxies.py
Original file line number Diff line number Diff line change
Expand Up @@ -102,7 +102,7 @@ def draw_deflector(self):
deflector["vel_disp"] = vel_disp
if deflector["e1_light"] == -1 or deflector["e2_light"] == -1:
e1_light, e2_light, e1_mass, e2_mass = elliptical_projected_eccentricity(
**deflector, **self._kwargs_mass2light
**deflector
)
deflector["e1_light"] = e1_light
deflector["e2_light"] = e2_light
Expand Down
127 changes: 4 additions & 123 deletions slsim/LsstSciencePipeline/lsst_science_pipeline.py
Original file line number Diff line number Diff line change
@@ -1,10 +1,8 @@
import numpy as np
from astropy.table import Table, vstack
from astropy.table import Column
from slsim.image_simulation import (
sharp_image,
lens_image,
lens_image_series,
)
from slsim.Util.param_util import transformmatrix_to_pixelscale
from scipy import interpolate
Expand Down Expand Up @@ -497,10 +495,11 @@ def cutout_image_psf_kernel(
calibFluxRadius=12,
):
"""This function extracts psf kernels from the dp0 cutout image at point source
image positions and deflector position. In the dp0.2 data, psf kernel vary with
coordinate and can be computed using given psf model.
image positions and deflector position. dp0 images are objects that has various
attributes. In the dp0.2 data, psf kernel vary with coordinate and can be computed
using given psf model.
:param dp0_image: cutout image from the dp0 data or any other image
:param dp0_image: cutout image from the dp0 data.
:param lens_class: class object containing all information of the lensing system
(e.g., Lens())
:param band: imaging band
Expand Down Expand Up @@ -798,124 +797,6 @@ def multiple_dp0_time_series_images_data(
return expo_data_list
return None


def variable_lens_injection(
lens_class, band, num_pix, transform_pix2angle, exposure_data
):
"""Injects variable lens to the dp0 time series data.
:param lens_class: Lens() object
:param band: imaging band
:param num_pix: number of pixels per axis
:param transform_pix2angle: transformation matrix (2x2) of pixels into coordinate
displacements
:param exposure_data: An astropy table of exposure data. It must contain calexp
images or generated noisy background image (column name should be
"time_series_images", these images are single exposure images of the same part
of the sky at different time), magnitude zero point (column name should be
"zero_point", these are zero point magnitudes for each single exposure images in
time series image) , psf kernel for each exposure (column name should be
"psf_kernel", these are pixel psf kernel for each single exposure images in time
series image), exposure time (column name should be "expo_time", these are
exposure time for each single exposure images in time series images),
observation time (column name should be "obs_time", these are observation time
in days for each single exposure images in time series images)
:return: Astropy table of injected lenses and exposure information of dp0 data
"""
# the range of observation time of single exposure images might be outside of the
# lightcurve time. So, we use random observation time from the lens class lightcurve
# time to ensure simulation of reasonable images.
observation_time = np.random.uniform(
min(lens_class.source.lightcurve_time),
max(lens_class.source.lightcurve_time),
size=len(exposure_data["obs_time"]),
)
observation_time.sort()
new_obs_time = Column(name="obs_time", data=observation_time)
exposure_data.replace_column("obs_time", new_obs_time)
lens_images = lens_image_series(
lens_class,
band=band,
mag_zero_point=exposure_data["zero_point"],
num_pix=num_pix,
psf_kernel=exposure_data["psf_kernel"],
transform_pix2angle=transform_pix2angle,
exposure_time=exposure_data["expo_time"],
t_obs=exposure_data["obs_time"],
)

final_image = []
for i in range(len(exposure_data["obs_time"])):
final_image.append(exposure_data["time_series_images"][i] + lens_images[i])
lens_col = Column(name="lens", data=lens_images)
final_image_col = Column(name="injected_lens", data=final_image)
if "lens" in exposure_data.colnames:
exposure_data.replace_column("lens", lens_col)
else:
exposure_data.add_column(lens_col)
if "injected_lens" in exposure_data.colnames:
exposure_data.replace_column("injected_lens", final_image_col)
else:
exposure_data.add_column(final_image_col)
return exposure_data


def multiple_variable_lens_injection(
lens_class_list,
band,
num_pix,
transform_matrices_list,
exposure_data_list,
output_file=None,
):
"""Injects multiple variable lenses to multiple dp0 time series data.
:param lens_class_list: list of Lens() object
:param band: imaging band
:param num_pix: number of pixels per axis
:param transform_matrices_list: list of transformation matrix (2x2) of pixels into
coordinate displacements for each exposure
:param exposure_data: list of astropy table of exposure data for each set of time
series images. It must contain calexp images or generated noisy background image
(column name should be "time_series_images", these images are single exposure
images of the same part of the sky at different time), magnitude zero point
(column name should be "zero_point", these are zero point magnitudes for each
single exposure images in time series image) , psf kernel for each exposure
(column name should be "psf_kernel", these are pixel psf kernel for each single
exposure images in time series image), exposure time (column name should be
"expo_time", these are exposure time for each single exposure images in time
series images), observation time (column name should be "obs_time", these are
observation time in days for each single exposure images in time series images)
:param output_file: path to the output FITS file where data will be saved
:return: list of astropy table of injected lenses and exposure information of dp0
data for each time series lenses. If output_file path is provided, it saves list
of these astropy table in fits file with the given name.
"""
final_images_catalog = []
for lens_class, transform_matrices, expo_data in zip(
lens_class_list, transform_matrices_list, exposure_data_list
):
variable_injected_image = variable_lens_injection(
lens_class,
band=band,
num_pix=num_pix,
transform_pix2angle=transform_matrices,
exposure_data=expo_data,
)
if output_file is None:
final_images_catalog.append(variable_injected_image)
else:
first_table = not os.path.exists(output_file)
if first_table:
variable_injected_image.write(output_file, overwrite=True)
first_table = False
else:
fits_append_table(output_file, variable_injected_image)
if len(final_images_catalog) > 1:
return final_images_catalog
return None


def measure_noise_level_in_RSP_coadd(RSP_coadd, N_pixels, plot=False):
np.random.seed(1)
"""Function to measure the noise level within a central square aperture of an RSP
Expand Down
70 changes: 0 additions & 70 deletions slsim/LsstSciencePipeline/opsim_pipeline.py
Original file line number Diff line number Diff line change
@@ -1,7 +1,5 @@
import numpy as np
from astropy.table import Table
from astropy.table import Column
from slsim.image_simulation import lens_image
import lenstronomy.Util.util as util
import lenstronomy.Util.kernel_util as kernel_util
import lenstronomy.Util.data_util as data_util
Expand Down Expand Up @@ -149,71 +147,3 @@ def opsim_time_series_images_data(

table_data_list.append(table_data)
return table_data_list


def opsim_variable_lens_injection(
lens_class, bands, num_pix, transform_pix2angle, exposure_data
):
"""Injects variable lens to the OpSim time series data (1 object).
:param lens_class: Lens() object
:param bands: list of imaging bands of interest
:param num_pix: number of pixels per axis
:param transform_pix2angle: transformation matrix (2x2) of pixels into coordinate
displacements
:param exposure_data: An astropy table of exposure data. One entry of
table_list_data generated from the opsim_time_series_images_data function. It
must contain the rms of background noise fluctuations (column name should be
"bkg_noise"), psf kernel for each exposure (column name should be "psf_kernel",
these are pixel psf kernel for each single exposure images in time series
image), observation time (column name should be "obs_time", these are
observation time in days for each single exposure images in time series images),
exposure time (column name should be "expo_time", these are exposure time for
each single exposure images in time series images), magnitude zero point (column
name should be "zero_point", these are zero point magnitudes for each single
exposure images in time series image), coordinates of the object (column name
should be "calexp_center"), these are the coordinates in (ra, dec), and the band
in which the observation is taken (column name should be "band").
:return: Astropy table of injected lenses and exposure information of dp0 data
"""

final_image = []

for obs in range(len(exposure_data["obs_time"])):

exposure_data_obs = exposure_data[obs]

if exposure_data_obs["band"] not in bands:
continue

if "bkg_noise" in exposure_data_obs.keys():
std_gaussian_noise = exposure_data_obs["bkg_noise"]
else:
std_gaussian_noise = None

lens_images = lens_image(
lens_class,
band=exposure_data_obs["band"],
mag_zero_point=exposure_data_obs["zero_point"],
num_pix=num_pix,
psf_kernel=exposure_data_obs["psf_kernel"],
transform_pix2angle=transform_pix2angle,
exposure_time=exposure_data_obs["expo_time"],
t_obs=exposure_data_obs["obs_time"],
std_gaussian_noise=std_gaussian_noise,
)

final_image.append(lens_images)

lens_col = Column(name="lens", data=final_image)
final_image_col = Column(name="injected_lens", data=final_image)

# Create a new Table with only the bands of interest
expo_bands = np.array([b for b in exposure_data["band"]])
mask = np.isin(expo_bands, bands)
exposure_data_new = exposure_data[mask]

# if len(exposure_data_new) > 0:
exposure_data_new.add_columns([lens_col, final_image_col])

return exposure_data_new
Loading

0 comments on commit 5cc87c5

Please sign in to comment.