KBMOD Visualization and Analysis

This notebook demonstrates the basic functionality for visualizing both the input data and the results.

Setup kbmod visualization demo

Before importing, make sure you have installed kbmod using pip install . in the root directory. Also be sure you are running with python3 and using the correct notebook kernel.

import math
import matplotlib.pyplot as plt
import numpy as np
import os

from kbmod.analysis.plotting import *
from kbmod.util_functions import load_deccam_layered_image
from kbmod.search import (
    ImageStack,
    PSF,
    Trajectory,
    get_stamps,
    get_median_stamp,
    get_mean_stamp,
    get_summed_stamp,
)
from kbmod.results import Results

# Data paths
im_path = "../data/small/"
res_path = "../data/fake_results/results.ecsv"

Visualizing the Input

We can visual input files with the kbmod.plotting.plot_img() function which will display the impage and optionally normalize the image (with the norm parameter) and add a title (with the title parameter).

We can load a LayeredImage given the path and filename to the FITS file as well as the PSF for the image. We use a default psf here.

p = PSF(1.0)
im = load_deccam_layered_image(im_path + "000000.fits", p)
print(f"Loaded a {im.get_width()} by {im.get_height()} image at time {im.get_obstime()}")

We can visualize the image the plotting library’s plot_img() function which is a wrapper around matplotlib’s imshow() function. The plot_img() function can take images in a variety of formats including a numpy array of the pixel values, a RawImage object, or a LayeredImage object.

Note: The data/demo images contain a single bright object, so the majority of the image should be empty with a single bright spot.

plot_image(im.get_science().image, norm=True, title="Image 000000.fits")

We can also load a collection of images and plot them with plot_multiple_images(). This function handles the arrangement of subplots given a single columns parameter that defines the total number of columns to use.

The function can take an ImageStack or list of images. In the case of an ImageStack it automatically labels each image with its time stamp.

files = [im_path + f for f in os.listdir(im_path) if ".fits" in f]
files.sort()
imgs = [load_deccam_layered_image(f, p) for f in files]

# Load the images.
stack = ImageStack(imgs)

num_images = stack.img_count()
print(f"Loaded {num_images} images.")

plot_multiple_images(stack, columns=4)

Create and Visualize Stamps

Stamps are a critical tool for analyzing and debugging proposed detections. They can be created automatically using the stamp creation utilities. It requires a few pieces of data:

• search_stack - provides the machinery for making predictions on the image (needed to handle the various corrections).

• trajectory - Contains the information about where to place the stamps (the underlying trajectory).

• stamp_radius - The radius in pixels.

# Create the trajectory with a given parameters and then the trajectory result.
trj = Trajectory()
trj.x = 11
trj.y = 27
trj.vx = 16.0
trj.vy = 3.3

# Create the stamps around this trajectory.
stamps = get_stamps(stack, trj, 20)

Now we can display the stamps around each predicted object position using the plotting library’s plot_multiple_images() function. Note that this time the function is taking a list of RawImage objects.

plot_multiple_images(stamps)

We can also compute the sum, mean, or median stamp. The stamp functions take a list of bools corresponding to whether each time is valid. These can be extracted from the result data. You can use an empty array (as we do below) to build the stamp out of all valid times.

# The coadds requires a vector of which indices to use.
inds = [True] * stack.img_count()

plot_multiple_images(
    [
        get_summed_stamp(stack, trj, 10, inds),
        get_mean_stamp(stack, trj, 10, inds),
        get_median_stamp(stack, trj, 10, inds),
    ],
    labels=["Summed", "Mean", "Median"],
    norm=True,
    columns=2,
)

Plotting Results

We can plot the results stored in the Results object from a KBMOD run. In this example, you need results in the form of either a Results object or saved files. For the purposes of this notebook you can use that fake result data (corresponding to the data/demo images) included in data/fake_results/results.ecsv. These results were generated by the KBMOD_Demo notebook.

We load all of the result data as a Results object using the class’s read_table function. We extract the first row to use for later examples.

Note that in previous versions (<=1.1) KBMOD saved output was also saved in a series of individual files such as results_DEMO.txt for the a file of the Trajectory objects. You can still load these using the from_trajectory_file function:

results = Results.from_trajectory_file("../data/fake_results/results_DEMO.txt")

results = Results.read_table("../data/fake_results/results.ecsv")
print(f"Loaded {len(results)} results.")

row0 = results[0]

We plot the results using the same helper functions as above and accessing the “stamp” column from the results:

plot_image(row0["stamp"])

or the “all_stamps” column:

plot_multiple_images(row0["all_stamps"], columns=5)

We can visualize the time series such as: psi_curve, phi_curve, or light_curve using the plot_time_series function. If we provide timestamps (such as with the top and bottom figure), they are used to both label the x-axis and set spacing.

all_times = [
    57130.19921875,
    57130.2109375,
    57130.21875,
    57131.19921875,
    57131.2109375,
    57131.21875,
    57132.19921875,
    57132.2109375,
    57132.21875,
    57133.19921875,
]

fig1 = plt.figure(figsize=(10, 15))
(ax1, ax2, ax3) = fig1.subplots(3, 1)

plot_time_series(row0["psi_curve"], times=all_times, ax=ax1, figure=fig1, title="Psi curve with time spacing")

plot_time_series(row0["psi_curve"], ax=ax2, figure=fig1, title="Psi curve with equal spacing")

plot_time_series(
    row0["psi_curve"],
    times=all_times,
    indices=[True, True, True, False, True, True, False, True, True, False],
    ax=ax3,
    figure=fig1,
    title="Psi curve with time spacing and invalid indices",
)

The plotting library also contains a helper function to display all this information from a single row of the Results data structure. The plot_result_row() function will read data out of the “stamp”, “psi_curve”, “phi_curve”, “obs_valid”, and “all_stamps” columns.

plot_result_row(row0)

If a column does not exist, the function will simply indicate this.

results.table.remove_column("all_stamps")
plot_result_row(results[0])