Module pyaurorax.tools.bounding_box.extract_metric
Extract various metrics from a given bounding box.
Functions
def azimuth(images: numpy.ndarray,
skymap: pyucalgarysrs.data.classes.Skymap,
azimuth_bounds: Sequence[float | int],
metric: Literal['mean', 'median', 'sum'] = 'median',
n_channels: int | None = None,
show_preview: bool = False) ‑> numpy.ndarray-
Expand source code
def azimuth(images: np.ndarray, skymap: Skymap, azimuth_bounds: Sequence[Union[int, float]], metric: Literal["mean", "median", "sum"] = "median", n_channels: Optional[int] = None, show_preview: bool = False) -> np.ndarray: """ Compute a metric of image data within an azimuthal boundary. 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]. skymap (pyaurorax.data.ucalgary.Skymap): The skymap corresponding to the image data. azimuth_bounds (Sequence[int, float]): A 2-element sequence specifying the azimuthal bounds from which to extract the metric. Anticipated order is [az_min, az_max]. metric (str): The name of the metric that is to be computed for the bounded area. Valid metrics are `mean`, `median`, `sum`. Default is `median`. n_channels (int): By default, function will assume the type of data passed as input - this argument can be used to manually specify the number of channels contained in image data. show_preview (bool): Plot a preview of the bounded area. Returns: A numpy.ndarray containing the metrics computed within azimuth range, for all image frames. Raises: ValueError: issue encountered with value supplied in parameter """ # Select individual azimuths from list az_0 = azimuth_bounds[0] az_1 = azimuth_bounds[1] # Determine if we are using single or 3 channel images = np.squeeze(images) if n_channels is not None: n_channels = n_channels else: n_channels = 1 if (len(images.shape) == 3): # single channel n_channels = 1 elif (len(images.shape) == 4): # three channel n_channels = 3 # Ensure that coordinates are valid if az_0 > 360 or az_0 < 0: raise ValueError("Invalid Azimuth: " + str(az_0)) elif az_1 > 360 or az_1 < 0: raise ValueError("Invalid Azimuth: " + str(az_1)) # Ensure that azimuths are properly ordered if az_0 > az_1: az_0, az_1 = az_1, az_0 # Ensure that this is a valid polygon if (az_0 == az_1): raise ValueError("Azimuth bounds defined with zero area.") # Obtain azimuth array from skymap if (skymap.full_azimuth is None): raise ValueError("Skymap 'full_azimuth' value is None. Unable to perform function") az = np.squeeze(skymap.full_azimuth) # Obtain indices into skymap within azimuth range bound_idx = np.where(np.logical_and(az > float(az_0), az < float(az_1))) # If boundaries contain no data, raise error if len(bound_idx[0]) == 0 or len(bound_idx[1]) == 0: raise ValueError("No data within desired bounds. Try a larger area.") # Slice out the bounded data, and plot preview if requested if n_channels == 1: bound_data = images[bound_idx[0], bound_idx[1], :] if show_preview: preview_img = scale_intensity(images[:, :, 0], top=230) preview_img[bound_idx[0], bound_idx[1]] = 255 plt.figure() plt.imshow(preview_img, cmap="grey", origin="lower") plt.title("Bounded Area Preview") plt.axis("off") plt.show() elif n_channels == 3: bound_data = images[bound_idx[0], bound_idx[1], :, :] if show_preview: preview_img = scale_intensity(images[:, :, :, 0], top=230) preview_img[bound_idx[0], bound_idx[1], 0] = 255 preview_img[bound_idx[0], bound_idx[1], 1:] = 0 plt.figure() plt.imshow(preview_img, origin="lower") plt.title("Bounded Area Preview") plt.axis("off") plt.show() else: raise ValueError("Unrecognized image format with shape: " + str(images.shape)) # Compute metric of interest if metric == 'median': result = np.median(bound_data, axis=0) elif metric == 'mean': result = np.mean(bound_data, axis=0) elif metric == 'sum': result = np.sum(bound_data, axis=0) else: raise ValueError("Metric " + str(metric) + " is not recognized.") return result
Compute a metric of image data within an azimuthal boundary.
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]. skymap
:Skymap
- The skymap corresponding to the image data.
azimuth_bounds
:Sequence[int, float]
- A 2-element sequence specifying the azimuthal bounds from which to extract the metric. Anticipated order is [az_min, az_max].
metric
:str
- The name of the metric that is to be computed for the bounded area. Valid metrics are
mean
,median
,sum
. Default ismedian
. n_channels
:int
- By default, function will assume the type of data passed as input - this argument can be used to manually specify the number of channels contained in image data.
show_preview
:bool
- Plot a preview of the bounded area.
Returns
A numpy.ndarray containing the metrics computed within azimuth range, for all image frames.
Raises
ValueError
- issue encountered with value supplied in parameter
def ccd(images: numpy.ndarray,
ccd_bounds: Sequence[int],
metric: Literal['mean', 'median', 'sum'] = 'median',
n_channels: int | None = None,
show_preview: bool = False) ‑> numpy.ndarray-
Expand source code
def ccd(images: np.ndarray, ccd_bounds: Sequence[int], metric: Literal["mean", "median", "sum"] = "median", n_channels: Optional[int] = None, show_preview: bool = False) -> np.ndarray: """ Compute a metric of image data within a CCD boundary. 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]. ccd_bounds (List[int]): A 4-element sequence specifying the (inclusive) CCD bounds from which to extract the metric. Anticipated order is [x_0, x_1, y_0, y_1]. metric (str): The name of the metric that is to be computed for the bounded area. Valid metrics are `mean`, `median`, `sum`. Defaults to `median`. n_channels (int): By default, function will assume the type of data passed as input - this argument can be used to manually specify the number of channels contained in image data. show_preview (bool): Plot a preview of the bounded area. Returns: A numpy.ndarray containing the metrics computed within CCD bounds, for all image frames. Raises: ValueError: issue encountered with value supplied in parameter """ # Select individual ccd x/y from list x_0 = ccd_bounds[0] x_1 = ccd_bounds[1] y_0 = ccd_bounds[2] y_1 = ccd_bounds[3] # Determine if we are using single or 3 channel images = np.squeeze(images) if n_channels is not None: n_channels = n_channels else: n_channels = 1 if (len(images.shape) == 3): # single channel n_channels = 1 elif (len(images.shape) == 4): # three channel n_channels = 3 # Ensure that coordinates are valid max_x = images.shape[0] max_y = images.shape[1] if y_0 > max_y or y_0 < 0: raise ValueError("CCD Y Coordinate " + str(y_0) + " out of range for image of shape " + str((max_y, max_x)) + ".") elif y_1 > max_y or y_1 < 0: raise ValueError("CCD Y Coordinate " + str(y_1) + " out of range for image of shape " + str((max_y, max_x)) + ".") elif x_0 > max_x or x_0 < 0: raise ValueError("CCD Y Coordinate " + str(x_0) + " out of range for image of shape " + str((max_y, max_x)) + ".") elif x_1 > max_x or x_1 < 0: raise ValueError("CCD Y Coordinate " + str(x_1) + " out of range for image of shape " + str((max_y, max_x)) + ".") # Ensure that coordinates are properly ordered if y_0 > y_1: y_0, y_1 = y_1, y_0 if x_0 > x_1: x_0, x_1 = x_1, x_0 # Ensure that this is a valid polygon if (y_0 == y_1) or (x_0 == x_1): raise ValueError("Polygon defined with zero area.") # Slice out the bounded data if n_channels == 1: bound_data = images[y_0:y_1, x_0:x_1, :] if show_preview: preview_img = scale_intensity(images[:, :, 0], top=230) preview_img[y_0:y_1, x_0:x_1] = 255 plt.figure() plt.imshow(preview_img, cmap="grey", origin="lower") plt.title("Bounded Area Preview") plt.axis("off") plt.show() elif n_channels == 3: bound_data = images[y_0:y_1, x_0:x_1, :, :] if show_preview: preview_img = scale_intensity(images[:, :, :, 0], top=230) preview_img[y_0:y_1, x_0:x_1, 0] = 255 preview_img[y_0:y_1, x_0:x_1, 1:] = 0 plt.figure() plt.imshow(preview_img, origin="lower") plt.title("Bounded Area Preview") plt.axis("off") plt.show() else: raise ValueError("Unrecognized image format with shape: " + str(images.shape)) # Compute metric of interest if metric == 'median': result = np.median(bound_data, axis=(0, 1)) elif metric == 'mean': result = np.mean(bound_data, axis=(0, 1)) elif metric == 'sum': result = np.sum(bound_data, axis=(0, 1)) else: raise ValueError("Metric " + str(metric) + " is not recognized.") return result
Compute a metric of image data within a CCD boundary.
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]. ccd_bounds
:List[int]
- A 4-element sequence specifying the (inclusive) CCD bounds from which to extract the metric. Anticipated order is [x_0, x_1, y_0, y_1].
metric
:str
- The name of the metric that is to be computed for the bounded area. Valid metrics are
mean
,median
,sum
. Defaults tomedian
. n_channels
:int
- By default, function will assume the type of data passed as input - this argument can be used to manually specify the number of channels contained in image data.
show_preview
:bool
- Plot a preview of the bounded area.
Returns
A numpy.ndarray containing the metrics computed within CCD bounds, for all image frames.
Raises
ValueError
- issue encountered with value supplied in parameter
def elevation(images: numpy.ndarray,
skymap: pyucalgarysrs.data.classes.Skymap,
elevation_bounds: Sequence[float | int],
metric: Literal['mean', 'median', 'sum'] = 'median',
n_channels: int | None = None,
show_preview: bool = False) ‑> numpy.ndarray-
Expand source code
def elevation(images: np.ndarray, skymap: Skymap, elevation_bounds: Sequence[Union[int, float]], metric: Literal["mean", "median", "sum"] = "median", n_channels: Optional[int] = None, show_preview: bool = False) -> np.ndarray: """ Compute a metric of image data within an elevation boundary. 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]. skymap (pyaurorax.data.ucalgary.Skymap): The skymap corresponding to the image data. elevation_bounds (Sequence): A 2-element sequence specifying the elevation bounds from which to extract the metric. Anticipated order is [el_min, el_max]. metric (str): The name of the metric that is to be computed for the bounded area. Valid metrics are `mean`, `median`, `sum`. Default is `median`. n_channels (int): By default, function will assume the type of data passed as input - this argument can be used to manually specify the number of channels contained in image data. show_preview (bool): Plot a preview of the bounded area. Returns: A numpy.ndarray containing the metrics computed within elevation range, for all image frames. Raises: ValueError: issue encountered with value supplied in parameter """ # Select individual elevations from list elev_0 = elevation_bounds[0] elev_1 = elevation_bounds[1] # Determine if we are using single or 3 channel images = np.squeeze(images) if n_channels is not None: n_channels = n_channels else: n_channels = 1 if (len(images.shape) == 3): # single channel n_channels = 1 elif (len(images.shape) == 4): # three channel n_channels = 3 # Ensure that coordinates are valid if elev_0 > 90 or elev_0 < 0: raise ValueError("Invalid Elevation: " + str(elev_0)) elif elev_1 > 90 or elev_1 < 0: raise ValueError("Invalid Elevation: " + str(elev_1)) # Ensure that elevations are properly ordered if elev_0 > elev_1: elev_0, elev_1 = elev_1, elev_0 # Ensure that this is a valid bounded area if (elev_0 == elev_1): raise ValueError("Elevation bounds defined with zero area.") # Obtain elevation array from skymap elev = np.squeeze(skymap.full_elevation) # Obtain indices into skymap within elevation range bound_idx = np.where(np.logical_and(elev >= float(elev_0), elev <= float(elev_1))) # If boundaries contain no data, raise error if len(bound_idx[0]) == 0 or len(bound_idx[1]) == 0: raise ValueError("No data within desired bounds. Try a larger area.") # Slice out the bounded data if n_channels == 1: bound_data = images[bound_idx[0], bound_idx[1], :] if show_preview: preview_img = scale_intensity(images[:, :, 0], top=230) preview_img[bound_idx[0], bound_idx[1]] = 255 plt.figure() plt.imshow(preview_img, cmap="grey", origin="lower") plt.title("Bounded Area Preview") plt.axis("off") plt.show() elif n_channels == 3: bound_data = images[bound_idx[0], bound_idx[1], :, :] if show_preview: preview_img = scale_intensity(images[:, :, :, 0], top=230) preview_img[bound_idx[0], bound_idx[1], 0] = 255 preview_img[bound_idx[0], bound_idx[1], 1:] = 0 plt.figure() plt.imshow(preview_img, origin="lower") plt.title("Bounded Area Preview") plt.axis("off") plt.show() else: raise ValueError("Unrecognized image format with shape: " + str(images.shape)) # Compute metric of interest if metric == 'median': result = np.median(bound_data, axis=0) elif metric == 'mean': result = np.mean(bound_data, axis=0) elif metric == 'sum': result = np.sum(bound_data, axis=0) else: raise ValueError("Metric " + str(metric) + " is not recognized.") return result
Compute a metric of image data within an elevation boundary.
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]. skymap
:Skymap
- The skymap corresponding to the image data.
elevation_bounds
:Sequence
- A 2-element sequence specifying the elevation bounds from which to extract the metric. Anticipated order is [el_min, el_max].
metric
:str
- The name of the metric that is to be computed for the bounded area. Valid metrics are
mean
,median
,sum
. Default ismedian
. n_channels
:int
- By default, function will assume the type of data passed as input - this argument can be used to manually specify the number of channels contained in image data.
show_preview
:bool
- Plot a preview of the bounded area.
Returns
A numpy.ndarray containing the metrics computed within elevation range, for all image frames.
Raises
ValueError
- issue encountered with value supplied in parameter
def geo(images: numpy.ndarray,
skymap: pyucalgarysrs.data.classes.Skymap,
altitude_km: float | int,
lonlat_bounds: Sequence[float | int],
metric: Literal['mean', 'median', 'sum'] = 'median',
n_channels: int | None = None,
show_preview: bool = False) ‑> numpy.ndarray-
Expand source code
def geo(images: np.ndarray, skymap: Skymap, altitude_km: Union[int, float], lonlat_bounds: Sequence[Union[int, float]], metric: Literal["mean", "median", "sum"] = "median", n_channels: Optional[int] = None, show_preview: bool = False) -> np.ndarray: """ Compute a metric of image data within a geographic lat/lon boundary. 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]. skymap (pyaurorax.data.ucalgary.Skymap): The skymap corresponding to the image data. altitude_km (int or float): The altitude of the image data in kilometers. lonlat_bounds (Sequence): A 4-element sequence specifying the lat/lon bounds from which to extract the metric. Anticipated order is [lon_0, lon_1, lat_0, lat_1]. metric (str): The name of the metric that is to be computed for the bounded area. Valid metrics are `mean`, `median`, `sum`. Default is `median`. n_channels (int): By default, function will assume the type of data passed as input - this argument can be used to manually specify the number of channels contained in image data. show_preview (bool): Plot a preview of the bounded area. Returns: A numpy.ndarray containing the metrics computed within elevation range, for all image frames. Raises: ValueError: issue encountered with value supplied in parameter """ # Select individual lats/lons from list lon_0 = lonlat_bounds[0] lon_1 = lonlat_bounds[1] lat_0 = lonlat_bounds[2] lat_1 = lonlat_bounds[3] # Determine if we are using single or 3 channel images = np.squeeze(images) if n_channels is not None: n_channels = n_channels else: n_channels = 1 if (len(images.shape) == 3): # single channel n_channels = 1 elif (len(images.shape) == 4): # three channel n_channels = 3 # Ensure that coordinates are valid if lat_0 > 90 or lat_0 < -90: raise ValueError("Invalid Latitude: " + str(lat_0)) elif lat_1 > 90 or lat_1 < -90: raise ValueError("Invalid Latitude: " + str(lat_1)) elif lon_0 > 360 or lon_0 < -180: raise ValueError("Invalid Longitude: " + str(lon_0)) elif lon_1 > 360 or lon_1 < -180: raise ValueError("Invalid Longitude: " + str(lon_1)) # Convert(0,360) longitudes to (-180,180) if entered as such if lon_0 > 180: lon_0 -= 360.0 if lon_1 > 180: lon_1 -= 360.0 # Ensure that coordinates are properly ordered if lat_0 > lat_1: lat_0, lat_1 = lat_1, lat_0 if lon_0 > lon_1: lon_0, lon_1 = lon_1, lon_0 # Ensure that this is a valid polygon if (lat_0 == lat_1) or (lon_0 == lon_1): raise ValueError("Polygon defined with zero area.") # Obtain lat/lon arrays from skymap interp_alts = skymap.full_map_altitude / 1000.0 if (altitude_km in interp_alts): altitude_idx = np.where(altitude_km == interp_alts) lats = np.squeeze(skymap.full_map_latitude[altitude_idx, :, :]) lons = np.squeeze(skymap.full_map_longitude[altitude_idx, :, :]) lons[np.where(lons > 180)] -= 360.0 # Fix skymap to be in (-180,180) format else: # Make sure altitude is in range that can be interpolated if (altitude_km < interp_alts[0]) or (altitude_km > interp_alts[2]): raise ValueError("Altitude " + str(altitude_km) + " outside valid range of " + str((interp_alts[0], interp_alts[2]))) # Initialze empty lat/lon arrays lats = np.full(np.squeeze(skymap.full_map_latitude[0, :, :]).shape, np.nan, dtype=skymap.full_map_latitude[0, :, :].dtype) lons = np.full(np.squeeze(skymap.full_map_latitude[0, :, :]).shape, np.nan, dtype=skymap.full_map_latitude[0, :, :].dtype) # Interpolate lats and lons at desired altitude for i in range(skymap.full_map_latitude.shape[1]): for j in range(skymap.full_map_latitude.shape[2]): pixel_lats = skymap.full_map_latitude[:, i, j] pixel_lons = skymap.full_map_longitude[:, i, j] if np.isnan(pixel_lats).any() or np.isnan(pixel_lons).any(): continue lats[i, j] = np.interp(altitude_km, interp_alts, pixel_lats) lons[i, j] = np.interp(altitude_km, interp_alts, pixel_lons) lons[np.where(lons > 180)] -= 360.0 # Fix skymap to be in (-180,180) format # Check that lat/lon range is reasonable min_skymap_lat = np.nanmin(lats) max_skymap_lat = np.nanmax(lats) min_skymap_lon = np.nanmin(lons) max_skymap_lon = np.nanmax(lons) if (lat_0 <= min_skymap_lat) or (lat_1 >= max_skymap_lat): raise ValueError(f"Latitude range supplied is outside the valid range for this skymap {(min_skymap_lat,max_skymap_lat)}.") if (lon_0 <= min_skymap_lon) or (lon_1 >= max_skymap_lon): raise ValueError(f"Longitude range supplied is outside the valid range for this skymap {(min_skymap_lon,max_skymap_lon)}.") # Obtain indices into skymap within lat/lon range bound_idx = np.where(np.logical_and.reduce((lats >= float(lat_0), lats <= float(lat_1), lons >= float(lon_0), lons <= float(lon_1)))) # If boundaries contain no data, raise error if len(bound_idx[0]) == 0 or len(bound_idx[1]) == 0: raise ValueError("No data within desired bounds. Try a larger area.") # Convert from skymap coords to image coords bound_idx = tuple(i - 1 for i in bound_idx) bound_idx = tuple(np.maximum(idx, 0) for idx in bound_idx) # Slice out the bounded data if n_channels == 1: bound_data = images[bound_idx[0], bound_idx[1], :] if show_preview: preview_img = scale_intensity(images[:, :, 0], top=230) preview_img[bound_idx[0], bound_idx[1]] = 255 plt.figure() plt.imshow(preview_img, cmap="grey", origin="lower") plt.title("Bounded Area Preview") plt.axis("off") plt.show() elif n_channels == 3: bound_data = images[bound_idx[0], bound_idx[1], :, :] if show_preview: preview_img = scale_intensity(images[:, :, :, 0], top=230) preview_img[bound_idx[0], bound_idx[1], 0] = 255 preview_img[bound_idx[0], bound_idx[1], 1:] = 0 plt.figure() plt.imshow(preview_img, origin="lower") plt.title("Bounded Area Preview") plt.axis("off") plt.show() else: raise ValueError("Unrecognized image format with shape: " + str(images.shape)) # Compute metric of interest if metric == 'median': result = np.median(bound_data, axis=0) elif metric == 'mean': result = np.mean(bound_data, axis=0) elif metric == 'sum': result = np.sum(bound_data, axis=0) else: raise ValueError("Metric " + str(metric) + " is not recognized.") return result
Compute a metric of image data within a geographic lat/lon boundary.
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]. skymap
:Skymap
- The skymap corresponding to the image data.
altitude_km
:int
orfloat
- The altitude of the image data in kilometers.
lonlat_bounds
:Sequence
- A 4-element sequence specifying the lat/lon bounds from which to extract the metric. Anticipated order is [lon_0, lon_1, lat_0, lat_1].
metric
:str
- The name of the metric that is to be computed for the bounded area. Valid metrics are
mean
,median
,sum
. Default ismedian
. n_channels
:int
- By default, function will assume the type of data passed as input - this argument can be used to manually specify the number of channels contained in image data.
show_preview
:bool
- Plot a preview of the bounded area.
Returns
A numpy.ndarray containing the metrics computed within elevation range, for all image frames.
Raises
ValueError
- issue encountered with value supplied in parameter
def mag(images: numpy.ndarray,
timestamp: datetime.datetime,
skymap: pyucalgarysrs.data.classes.Skymap,
altitude_km: float | int,
lonlat_bounds: Sequence[float | int],
metric: Literal['mean', 'median', 'sum'] = 'median',
n_channels: int | None = None,
show_preview: bool = False) ‑> numpy.ndarray-
Expand source code
def mag(images: np.ndarray, timestamp: datetime.datetime, skymap: Skymap, altitude_km: Union[int, float], lonlat_bounds: Sequence[Union[int, float]], metric: Literal["mean", "median", "sum"] = "median", n_channels: Optional[int] = None, show_preview: bool = False) -> np.ndarray: """ Compute a metric of image data within a magnetic lat/lon boundary. 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]. skymap (pyaurorax.data.ucalgary.Skymap): The skymap corresponding to the image data. altitude_km (int or float): The altitude of the image data in kilometers. lonlat_bounds (Sequence): A 4-element sequence specifying the magnetic lat/lon bounds from which to extract the metric. Anticipated order is [lon_0, lon_1, lat_0, lat_1]. metric (str): The name of the metric that is to be computed for the bounded area. Valid metrics are `mean`, `median`, `sum`. Default is `median`. n_channels (int): By default, function will assume the type of data passed as input - this argument can be used to manually specify the number of channels contained in image data. show_preview (bool): Plot a preview of the bounded area. Returns: A numpy.ndarray containing the metrics computed within elevation range, for all image frames. Raises: ValueError: issue encountered with value supplied in parameter """ # Select individual lats/lons from list lon_0 = lonlat_bounds[0] lon_1 = lonlat_bounds[1] lat_0 = lonlat_bounds[2] lat_1 = lonlat_bounds[3] # Determine if we are using single or 3 channel images = np.squeeze(images) if n_channels is not None: n_channels = n_channels else: n_channels = 1 if (len(images.shape) == 3): # single channel n_channels = 1 elif (len(images.shape) == 4): # three channel n_channels = 3 # Ensure that coordinates are valid if lat_0 > 90 or lat_0 < -90: raise ValueError("Invalid Latitude: " + str(lat_0)) elif lat_1 > 90 or lat_1 < -90: raise ValueError("Invalid Latitude: " + str(lat_1)) elif lon_0 > 360 or lon_0 < -180: raise ValueError("Invalid Longitude: " + str(lon_0)) elif lon_1 > 360 or lon_1 < -180: raise ValueError("Invalid Longitude: " + str(lon_1)) # Convert(0,360) longitudes to (-180,180) if entered as such if lon_0 > 180: lon_0 -= 360.0 if lon_1 > 180: lon_1 -= 360.0 # Ensure that coordinates are properly ordered if lat_0 > lat_1: lat_0, lat_1 = lat_1, lat_0 if lon_0 > lon_1: lon_0, lon_1 = lon_1, lon_0 # Ensure that this is a valid polygon if (lat_0 == lat_1) or (lon_0 == lon_1): raise ValueError("Polygon defined with zero area.") # Obtain lat/lon arrays from skymap interp_alts = skymap.full_map_altitude / 1000.0 if (altitude_km in interp_alts): altitude_idx = np.where(altitude_km == interp_alts) lats = np.squeeze(skymap.full_map_latitude[altitude_idx, :, :]) lons = np.squeeze(skymap.full_map_longitude[altitude_idx, :, :]) lons[np.where(lons > 180)] -= 360.0 # Fix skymap to be in (-180,180) format else: # Make sure altitude is in range that can be interpolated if (altitude_km < interp_alts[0]) or (altitude_km > interp_alts[2]): raise ValueError("Altitude " + str(altitude_km) + " outside valid range of " + str((interp_alts[0], interp_alts[2]))) # Initialze empty lat/lon arrays lats = np.full(np.squeeze(skymap.full_map_latitude[0, :, :]).shape, np.nan, dtype=skymap.full_map_latitude[0, :, :].dtype) lons = np.full(np.squeeze(skymap.full_map_latitude[0, :, :]).shape, np.nan, dtype=skymap.full_map_latitude[0, :, :].dtype) # Interpolate lats and lons at desired altitude for i in range(skymap.full_map_latitude.shape[1]): for j in range(skymap.full_map_latitude.shape[2]): pixel_lats = skymap.full_map_latitude[:, i, j] pixel_lons = skymap.full_map_longitude[:, i, j] if np.isnan(pixel_lats).any() or np.isnan(pixel_lons).any(): continue lats[i, j] = np.interp(altitude_km, interp_alts, pixel_lats) lons[i, j] = np.interp(altitude_km, interp_alts, pixel_lons) lons[np.where(lons > 180)] -= 360.0 # Fix skymap to be in (-180,180) format # Convert skymap to magnetic coords mag_lats, mag_lons, mag_alts = aacgmv2.convert_latlon_arr(lats.flatten(), lons.flatten(), (lons * 0.0).flatten(), timestamp, method_code='G2A') mag_lats = np.reshape(mag_lats, lats.shape) mag_lons = np.reshape(mag_lons, lons.shape) # Check that lat/lon range is reasonable min_skymap_lat = np.nanmin(mag_lats) max_skymap_lat = np.nanmax(mag_lats) min_skymap_lon = np.nanmin(mag_lons) max_skymap_lon = np.nanmax(mag_lons) if (lat_0 <= min_skymap_lat) or (lat_1 >= max_skymap_lat): raise ValueError(f"Latitude range supplied is outside the valid range for this skymap {(min_skymap_lat,max_skymap_lat)}.") if (lon_0 <= min_skymap_lon) or (lon_1 >= max_skymap_lon): raise ValueError(f"Latitude range supplied is outside the valid range for this skymap {(min_skymap_lon,max_skymap_lon)}.") # Obtain indices into skymap within lat/lon range bound_idx = np.where( np.logical_and.reduce(( mag_lats >= float(lat_0), mag_lats <= float(lat_1), mag_lons >= float(lon_0), mag_lons <= float(lon_1), ))) # If boundaries contain no data, raise error if len(bound_idx[0]) == 0 or len(bound_idx[1]) == 0: raise ValueError("No data within desired bounds. Try a larger area.") # Convert from skymap coords to image coords bound_idx = tuple(i - 1 for i in bound_idx) bound_idx = tuple(np.maximum(idx, 0) for idx in bound_idx) # Slice out the bounded data if n_channels == 1: bound_data = images[bound_idx[0], bound_idx[1], :] if show_preview: preview_img = scale_intensity(images[:, :, 0], top=230) preview_img[bound_idx[0], bound_idx[1]] = 255 plt.figure() plt.imshow(preview_img, cmap="grey", origin="lower") plt.title("Bounded Area Preview") plt.axis("off") plt.show() elif n_channels == 3: bound_data = images[bound_idx[0], bound_idx[1], :, :] if show_preview: preview_img = scale_intensity(images[:, :, :, 0], top=230) preview_img[bound_idx[0], bound_idx[1], 0] = 255 preview_img[bound_idx[0], bound_idx[1], 1:] = 0 plt.figure() plt.imshow(preview_img, origin="lower") plt.title("Bounded Area Preview") plt.axis("off") plt.show() else: raise ValueError("Unrecognized image format with shape: " + str(images.shape)) # Compute metric of interest if metric == 'median': result = np.median(bound_data, axis=0) elif metric == 'mean': result = np.mean(bound_data, axis=0) elif metric == 'sum': result = np.sum(bound_data, axis=0) else: raise ValueError("Metric " + str(metric) + " is not recognized.") return result
Compute a metric of image data within a magnetic lat/lon boundary.
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]. skymap
:Skymap
- The skymap corresponding to the image data.
altitude_km
:int
orfloat
- The altitude of the image data in kilometers.
lonlat_bounds
:Sequence
- A 4-element sequence specifying the magnetic lat/lon bounds from which to extract the metric. Anticipated order is [lon_0, lon_1, lat_0, lat_1].
metric
:str
- The name of the metric that is to be computed for the bounded area. Valid metrics are
mean
,median
,sum
. Default ismedian
. n_channels
:int
- By default, function will assume the type of data passed as input - this argument can be used to manually specify the number of channels contained in image data.
show_preview
:bool
- Plot a preview of the bounded area.
Returns
A numpy.ndarray containing the metrics computed within elevation range, for all image frames.
Raises
ValueError
- issue encountered with value supplied in parameter