Module pyaurorax.tools.keogram

Generate keograms.

# Copyright 2024 University of Calgary
#
# Licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at
#
#     http://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.
"""
Generate keograms.
"""

from ._create import create
from ._create_custom import create_custom

__all__ = [
    "create",
    "create_custom",
]

Functions

def create(images: numpy.ndarray, timestamp: List[datetime.datetime], axis: int = 0) ‑> Keogram

Create a keogram from a set of images.

Args

images : numpy.ndarray

A set of images. Normally this would come directly from a data read call, but can also be any arbitrary set of images. It is anticipated that the order of axes is [rows, cols, num_images] or [row, cols, channels, num_images]. If it is not, then be sure to specify the axis parameter accordingly.

timestamp : List[datetime.datetime]

A list of timestamps corresponding to each image.

axis : int

The axis to extract the keogram slice from. Default is 0, meaning the rows (or Y) axis.

Returns

A Keogram object.

Raises

ValueError

issue with supplied parameters.

Expand source code

def create(images: np.ndarray, timestamp: List[datetime.datetime], axis: int = 0) -> Keogram:
    """
    Create a keogram from a set of images.

    Args:
        images (numpy.ndarray): 
            A set of images. Normally this would come directly from a data `read` call, but can also 
            be any arbitrary set of images. It is anticipated that the order of axes is [rows, cols, num_images]
            or [row, cols, channels, num_images]. If it is not, then be sure to specify the `axis` parameter
            accordingly.

        timestamp (List[datetime.datetime]): 
            A list of timestamps corresponding to each image.

        axis (int): 
            The axis to extract the keogram slice from. Default is `0`, meaning the rows (or Y) axis.

    Returns:
        A `pyaurorax.tools.Keogram` object.

    Raises:
        ValueError: issue with supplied parameters.
    """
    # set y axis
    ccd_y = np.arange(0, images.shape[axis])

    # determine if we are single or 3 channel
    n_channels = 1
    if (len(images.shape) == 3):
        # single channel
        n_channels = 1
    elif (len(images.shape) == 4):
        # three channel
        n_channels = 3
    else:
        ValueError("Unable to determine number of channels based on the supplied images. Make sure you are supplying a " +
                   "[rows,cols,images] or [rows,cols,channels,images] sized array.")

    # initialize keogram data
    n_rows = images.shape[0]
    n_imgs = images.shape[-1]
    if (n_channels == 1):
        keo_arr = np.full([n_rows, n_imgs], 0, dtype=images.dtype)
    else:
        keo_arr = np.full([n_rows, n_imgs, n_channels], 0, dtype=images.dtype)

    # extract the keogram slices
    middle_column_idx = int(np.floor((images.shape[1]) / 2 - 1))
    for img_idx in range(0, n_imgs):
        if (n_channels == 1):
            # single channel
            frame = images[:, :, img_idx]
            frame_middle_slice = frame[:, middle_column_idx]
            keo_arr[:, img_idx] = frame_middle_slice
        else:
            # 3-channel
            frame = images[:, :, :, img_idx]
            frame_middle_slice = frame[:, middle_column_idx, :]
            keo_arr[:, img_idx, :] = frame_middle_slice

    # create the keogram object
    keo_obj = Keogram(data=keo_arr, slice_idx=middle_column_idx, timestamp=timestamp, ccd_y=ccd_y)

    # return
    return keo_obj

def create_custom(images: numpy.ndarray, timestamp: List[datetime.datetime], coordinate_system: Literal['ccd', 'mag', 'geo'], width: int, x_locs: Union[List[Union[float, int]], numpy.ndarray], y_locs: Union[List[Union[float, int]], numpy.ndarray], preview: bool = False, skymap: Optional[pyucalgarysrs.data.classes.Skymap] = None, altitude_km: Union[int, float, ForwardRef(None)] = None, metric: Literal['mean', 'median', 'sum'] = 'median') ‑> Keogram

Create a keogram, from a custom slice of a set of images. The slice used is defined by a set of points, in CCD, geographic, or geomagnetic coordinates, within the bounds of the image data. Keogram is created from the bottom up, meaning the first point will correspond to the bottom of the keogram data.

Args

images : numpy.ndarray

A set of images. Normally this would come directly from a data read call, but can also be any arbitrary set of images. It is anticipated that the order of axes is [rows, cols, num_images] or [row, cols, channels, num_images]. If it is not, then be sure to specify the axis parameter accordingly.

timestamp : List[datetime.datetime]

A list of timestamps corresponding to each image.

coordinate_system (str): The coordinate system in which input points are defined. Valid options are "ccd", "geo", or "mag".

width (int): Width of the desired keogram slice, in CCD pixel units.

x_locs (Sequence[float | int]): Sequence of points giving the x-coordinates that define a path through the image data, from which to build the keogram.

y_locs (Sequence[float | int]): Sequence of points giving the y-coordinates that define a path through the image data, from which to build the keogram.

preview (Optional[bool]): When True, the first frame in images will be displayed, with the keogram slice plotted.

skymap (Skymap): The skymap to use in georeferencing when working in geographic or magnetic coordinates.

altitude_km (float | int): The altitude of the image data, in km, to use in georeferencing when working in goegraphic or magnetic coordinates.

metric (str): The metric used to compute values for each keogram pixel. Valid options are "median", "mean", and "sum". Defaults to "median".

Returns

A Keogram object. Raises:

Expand source code

def create_custom(
    images: np.ndarray,
    timestamp: List[datetime.datetime],
    coordinate_system: Literal["ccd", "geo", "mag"],
    width: int,
    x_locs: Union[List[Union[float, int]], np.ndarray],
    y_locs: Union[List[Union[float, int]], np.ndarray],
    preview: bool = False,
    skymap: Optional[Skymap] = None,
    altitude_km: Optional[Union[float, int]] = None,
    metric: Literal["mean", "median", "sum"] = "median",
) -> Keogram:
    """
    Create a keogram, from a custom slice of a set of images. The slice used is defined by a set of points, 
    in CCD, geographic, or geomagnetic coordinates, within the bounds of the image data. Keogram is created
    from the bottom up, meaning the first point will correspond to the bottom of the keogram data.

    Args:
        images (numpy.ndarray): 
            A set of images. Normally this would come directly from a data `read` call, but can also 
            be any arbitrary set of images. It is anticipated that the order of axes is [rows, cols, num_images]
            or [row, cols, channels, num_images]. If it is not, then be sure to specify the `axis` parameter
            accordingly.

        timestamp (List[datetime.datetime]): 
            A list of timestamps corresponding to each image.
        
        coordinate_system (str):
            The coordinate system in which input points are defined. Valid options are "ccd", "geo", or "mag".
        
        width (int):
            Width of the desired keogram slice, in CCD pixel units.

        x_locs (Sequence[float | int]):
            Sequence of points giving the x-coordinates that define a path through the image data, from
            which to build the keogram.

        y_locs (Sequence[float | int]):
            Sequence of points giving the y-coordinates that define a path through the image data, from
            which to build the keogram.

        preview (Optional[bool]):
            When True, the first frame in images will be displayed, with the keogram slice plotted.

        skymap (Skymap):
            The skymap to use in georeferencing when working in geographic or magnetic coordinates.

        altitude_km (float | int):
            The altitude of the image data, in km, to use in georeferencing when working in goegraphic
            or magnetic coordinates.

        metric (str):
            The metric used to compute values for each keogram pixel. Valid options are "median", "mean",
            and "sum". Defaults to "median".
        
    Returns:
        A `pyaurorax.tools.Keogram` object.

    Raises:
    """

    # If using CCD coordinates we don't need a skymao or altitude
    if (coordinate_system == 'ccd') and (skymap is not None or altitude_km is not None):
        raise ValueError("Confliction in passing a Skymap in when working in CCD coordinates. Skymap is obsolete.")

    # convert any lists to np.arrays  and check shape
    x_locs = np.array(x_locs)
    y_locs = np.array(y_locs)

    # determine if we are single or 3 channel
    n_channels = 1
    if (len(images.shape) == 3):
        # single channel
        n_channels = 1
    elif (len(images.shape) == 4):
        # three channel
        n_channels = 3
    else:
        ValueError("Unable to determine number of channels based on the supplied images. Make sure you are supplying a " +
                   "[rows,cols,images] or [rows,cols,channels,images] sized array.")

    # Initialize empty keogram array
    keo_arr = np.squeeze(np.full((x_locs.shape[0] - 1, len(timestamp), n_channels), 0))

    if len(x_locs.shape) != 1:
        raise ValueError(f"X coordinates may not be multidimensional. Sequence passed with shape {x_locs.shape}")
    if len(y_locs.shape) != 1:
        raise ValueError(f"Y coordinates may not be multidimensional. Sequence passed with shape {y_locs.shape}")
    if len(x_locs.shape) != len(y_locs.shape):
        raise ValueError(f"X and Y coordinates must have same length. Sequences passed with shapes {x_locs.shape} and {y_locs.shape}")

    # Convert lat/lon coordinates to CCD
    if coordinate_system == 'mag':
        if (skymap is None or altitude_km is None):
            raise ValueError("When magnetic coordinates, a Skymap object and Altitude must be passed in through the skymap argument.")
        x_locs, y_locs = __convert_latlon_to_ccd(x_locs, y_locs, timestamp, skymap, altitude_km, magnetic=True)
    elif coordinate_system == 'geo':
        if (skymap is None or altitude_km is None):
            raise ValueError("When geographic coordinates, a Skymap object and Altitude must be passed in through the skymap argument.")
        x_locs, y_locs = __convert_latlon_to_ccd(x_locs, y_locs, timestamp, skymap, altitude_km, magnetic=False)

    # We will use the first image as a preview
    if n_channels == 1:
        preview_img = images.copy()[:, :, 0]
    else:
        preview_img = images.copy()[:, :, :, 0]

    # Now working in CCD Coordinates
    x_max = images.shape[1] - 1
    y_max = images.shape[0] - 1

    # Remove any points that are not within the image CCD
    parsed_x_locs = []
    parsed_y_locs = []
    for i in range(x_locs.shape[0]):
        x = x_locs[i]
        y = y_locs[i]

        if x < 0 or x > x_max:
            continue
        if y < 0 or y > y_max:
            continue

        parsed_x_locs.append(x)
        parsed_y_locs.append(y)
    x_locs = np.array(parsed_x_locs)
    y_locs = np.array(parsed_y_locs)

    # Make sure all supplied points are within image bounds # NOTE: This should be good to remove as it shouldn't hit, testing needed
    if len(np.where(np.logical_or(x_locs < 0, x_locs > x_max))[0]) != 0:
        raise ValueError("The following CCD coordinates passed in through x_locs are outside of the CCD image range: " +
                         str(x_locs[np.where(np.logical_or(x_locs < 0, x_locs > x_max))]))
    if len(np.where(np.logical_or(y_locs < 0, y_locs > y_max))[0]) != 0:
        raise ValueError("The following CCD coordinates passed in through y_locs are outside of the CCD image range: " +
                         str(y_locs[np.where(np.logical_or(y_locs < 0, y_locs > y_max))]))

    # Iterate points in pairs of two
    path_counter = 0
    for i in range(x_locs.shape[0] - 1):

        # Points of concern for this iteration
        x_0 = x_locs[i]
        x_1 = x_locs[i + 1]
        y_0 = y_locs[i]
        y_1 = y_locs[i + 1]

        # Compute the unit vector between the two points
        dx = x_1 - x_0
        dy = y_1 - y_0
        length = np.sqrt(dx**2 + dy**2)
        if length == 0:
            continue

        dx /= length
        dy /= length

        # Compute orthogonal unit vector
        perp_dx = -dy
        perp_dy = dx

        # Calculate (+/-) offsets for each perpendicular direction
        offset1_x = perp_dx * width / 2
        offset1_y = perp_dy * width / 2
        offset2_x = -perp_dx * width / 2
        offset2_y = -perp_dy * width / 2

        # Calculate vertices in correct order for this polygon
        vertex1 = (int(x_0 + offset1_x), int(y_0 + offset1_y))
        vertex2 = (int(x_1 + offset1_x), int(y_1 + offset1_y))
        vertex3 = (int(x_1 + offset2_x), int(y_1 + offset2_y))
        vertex4 = (int(x_0 + offset2_x), int(y_0 + offset2_y))

        # Append vertices in the correct order to form a closed polygon
        vertices = [vertex1, vertex2, vertex3, vertex4]

        # Obtain the indexes into the image of this polygon
        indices_inside = __indices_in_polygon(vertices, (images.shape[0], images.shape[1]))

        if np.any(np.array(indices_inside.shape) == 0):
            continue

        row_idx, col_idx = zip(*indices_inside)

        if n_channels == 1:
            # Update the preview image
            preview_img[row_idx, col_idx] = np.iinfo(preview_img.dtype).max

            # Extract metric from all images and add to keogram
            if metric == 'median':
                pixel_keogram = np.median(images[row_idx, col_idx, :], axis=0)
            elif metric == 'mean':
                pixel_keogram = np.mean(images[row_idx, col_idx, :], axis=0)
            elif metric == 'sum':
                pixel_keogram = np.sum(images[row_idx, col_idx, :], axis=0)
            keo_arr[i, :] = pixel_keogram
        elif n_channels == 3:
            # Update the preview image
            preview_img[row_idx, col_idx, :] = np.iinfo(preview_img.dtype).max

            # Extract metric from all images and add to keogram
            if metric == 'median':
                r_pixel_keogram = np.floor(np.median(images[row_idx, col_idx, 0, :], axis=0))
                g_pixel_keogram = np.floor(np.median(images[row_idx, col_idx, 1, :], axis=0))
                b_pixel_keogram = np.floor(np.median(images[row_idx, col_idx, 2, :], axis=0))
            elif metric == 'mean':
                r_pixel_keogram = np.floor(np.mean(images[row_idx, col_idx, 0, :], axis=0))
                g_pixel_keogram = np.floor(np.mean(images[row_idx, col_idx, 1, :], axis=0))
                b_pixel_keogram = np.floor(np.mean(images[row_idx, col_idx, 2, :], axis=0))
            elif metric == 'sum':
                r_pixel_keogram = np.floor(np.sum(images[row_idx, col_idx, 0, :], axis=0))
                g_pixel_keogram = np.floor(np.sum(images[row_idx, col_idx, 1, :], axis=0))
                b_pixel_keogram = np.floor(np.sum(images[row_idx, col_idx, 2, :], axis=0))

            keo_arr[i, :, 0] = r_pixel_keogram
            keo_arr[i, :, 1] = g_pixel_keogram
            keo_arr[i, :, 2] = b_pixel_keogram

        path_counter += 1

    if path_counter == 0:
        raise ValueError("Could not form keogram path. First ensure that coordinates are within image range. Then " +
                         "try increasing 'width' or decreasing number of points in input coordinates.")

    # Create keogram object
    keo_obj = Keogram(data=keo_arr, timestamp=timestamp)

    if preview:
        plt.figure()
        plt.imshow(preview_img, cmap='gray', origin='lower')
        plt.axis("off")
        plt.title("Keogram Domain Preview")
        plt.show()

    return keo_obj