IceFloeTracker.jl
Overview
IceFloeTracker.jl is a collection of routines and tools for processing remote sensing imagery, identifying sea ice floes, and tracking the displacement and rotation of ice floes across multiple images. It can be used either standalone to create custom processing pathways or with the Ice Floe Tracker Pipeline.
Algorithm components
The Ice Floe Tracker (IFT) package includes the core functions for the three main steps of the algorithm. These functions can be used independently and can be customized for specific use cases.
Preprocessing
IFT operates on optical satellite imagery. The main functions are designed with "true color" and "false color" imagery in mind, and have thus far primarily been tested on imagery from the Moderate Resolution Imaging Spectroradiometer (MODIS) from the NASA Aqua and Terra satellites. The preprocessing routines mask land and cloud features, and aim to adjust and sharpen the remainder of the images to amplify the contrast along the edges of sea ice floes. (TBD: Link to main preprocessing page)
Segmentation
The IFT segmentation functions include functions for semantic segmentation (pixel-by-pixel assignment into predefined categories) and object-based segmentation (groupings of pixels into distinct objects). The semantic segmentation steps use $k$-means to group pixels into water and ice regions. A combination of watershed functions, morphological operations, and further applications of $k$-means are used to identify candidate ice floes. (TBD: Link to main segmentation page)
Tracking
Ice floe tracking is carried out by comparing the shapes produced in the segmentation step. Shapes with similar area are rotated until the difference in surface area is minimized, and then the edge shapes are compared using a Ѱ-s curve. If thresholds for correlation and area differences are met, then the floe with the best correlation and smallest area differences are considered matches and the objects are assigned the same label. In the end, trajectories for individual floes are recorded in a dataframe.
Developers
IceFloeTracker.jl is a product of the Wilhelmus Lab at Brown University, led by Monica M. Wilhelmus. The original algorithm was developed by Rosalinda Lopez-Acosta during her PhD work at University of California Riverside, advised by Dr. Wilhelmus. The translation of the original Matlab code into the current modular, open source Julia package has been carried out in conjunction with the Center for Computing and Visualization at Brown University. Contributors include Daniel Watkins, Maria Isabel Restrepo, Carlos Paniagua, Tim DiVoll, John Holland, and Bradford Roarr.
Citing
If you use IceFloeTracker.jl in research, teaching, or elsewhere, please mention the IceFloeTracker package and cite our journal article outlining the algorithm:
Lopez-Acosta et al., (2019). Ice Floe Tracker: An algorithm to automatically retrieve Lagrangian trajectories via feature matching from moderate-resolution visual imagery. Remote Sensing of Environment, 234(111406), doi:10.1016/j.rse.2019.111406.
Papers using Ice Floe Tracker
- Manucharyan, Lopez-Acosta, and Wilhelmus (2022)*. Spinning ice floes reveal intensification of mesoscale eddies in the western Arctic Ocean. Scientific Reports, 12(7070), doi:10.1038/s41598-022-10712-z
- Watkins, Bliss, Hutchings, and Wilhelmus (2023)*. Evidence of Abrupt Transitions Between Sea Ice Dynamical Regimes in the East Greenland Marginal Ice Zone. Geophysical Research Letters, 50(e2023GL103558), pp. 1-10, doi:10.1029/2023GL103558
*Papers using data from the Matlab implementation of Ice Floe Tracker.
Functions
Base.isequal
— Methodisequal(matchedpairs1::MatchedPairs, matchedpairs2::MatchedPairs)
Return true
if matchedpairs1
and matchedpairs2
are equal, false
otherwise.
Base.sort!
— Methodsort!(tracked::Tracked)
Sort the floes in tracked
by area in descending order.
IceFloeTracker._adjust_histogram
— Method_adjust_histogram(masked_view, nbins, rblocks, cblocks, clip)
Perform adaptive histogram equalization to a masked image. To be invoked within imsharpen
.
Arguments
masked_view
: input image in truecolor
See imsharpen
for a description of the remaining arguments
IceFloeTracker._bin9todec
— Method_bin9todec(v)
Get decimal representation of a bit vector v
with the leading bit at its leftmost posistion.
Example
julia> _bin9todec([0 0 0 0 0 0 0 0 0])
+IceFloeTracker.jl · IceFloeTracker.jl IceFloeTracker.jl
Overview
IceFloeTracker.jl is a collection of routines and tools for processing remote sensing imagery, identifying sea ice floes, and tracking the displacement and rotation of ice floes across multiple images. It can be used either standalone to create custom processing pathways or with the Ice Floe Tracker Pipeline.
Algorithm components
The Ice Floe Tracker (IFT) package includes the core functions for the three main steps of the algorithm. These functions can be used independently and can be customized for specific use cases.
Preprocessing
IFT operates on optical satellite imagery. The main functions are designed with "true color" and "false color" imagery in mind, and have thus far primarily been tested on imagery from the Moderate Resolution Imaging Spectroradiometer (MODIS) from the NASA Aqua and Terra satellites. The preprocessing routines mask land and cloud features, and aim to adjust and sharpen the remainder of the images to amplify the contrast along the edges of sea ice floes. (TBD: Link to main preprocessing page)
Segmentation
The IFT segmentation functions include functions for semantic segmentation (pixel-by-pixel assignment into predefined categories) and object-based segmentation (groupings of pixels into distinct objects). The semantic segmentation steps use $k$-means to group pixels into water and ice regions. A combination of watershed functions, morphological operations, and further applications of $k$-means are used to identify candidate ice floes. (TBD: Link to main segmentation page)
Tracking
Ice floe tracking is carried out by comparing the shapes produced in the segmentation step. Shapes with similar area are rotated until the difference in surface area is minimized, and then the edge shapes are compared using a Ѱ-s curve. If thresholds for correlation and area differences are met, then the floe with the best correlation and smallest area differences are considered matches and the objects are assigned the same label. In the end, trajectories for individual floes are recorded in a dataframe.
Developers
IceFloeTracker.jl is a product of the Wilhelmus Lab at Brown University, led by Monica M. Wilhelmus. The original algorithm was developed by Rosalinda Lopez-Acosta during her PhD work at University of California Riverside, advised by Dr. Wilhelmus. The translation of the original Matlab code into the current modular, open source Julia package has been carried out in conjunction with the Center for Computing and Visualization at Brown University. Contributors include Daniel Watkins, Maria Isabel Restrepo, Carlos Paniagua, Tim DiVoll, John Holland, and Bradford Roarr.
Citing
If you use IceFloeTracker.jl in research, teaching, or elsewhere, please mention the IceFloeTracker package and cite our journal article outlining the algorithm:
Lopez-Acosta et al., (2019). Ice Floe Tracker: An algorithm to automatically retrieve Lagrangian trajectories via feature matching from moderate-resolution visual imagery. Remote Sensing of Environment, 234(111406), doi:10.1016/j.rse.2019.111406.
Papers using Ice Floe Tracker
- Manucharyan, Lopez-Acosta, and Wilhelmus (2022)*. Spinning ice floes reveal intensification of mesoscale eddies in the western Arctic Ocean. Scientific Reports, 12(7070), doi:10.1038/s41598-022-10712-z
- Watkins, Bliss, Hutchings, and Wilhelmus (2023)*. Evidence of Abrupt Transitions Between Sea Ice Dynamical Regimes in the East Greenland Marginal Ice Zone. Geophysical Research Letters, 50(e2023GL103558), pp. 1-10, doi:10.1029/2023GL103558
*Papers using data from the Matlab implementation of Ice Floe Tracker.
Functions
Base.isequal
— Methodisequal(matchedpairs1::MatchedPairs, matchedpairs2::MatchedPairs)
Return true
if matchedpairs1
and matchedpairs2
are equal, false
otherwise.
sourceBase.sort!
— Methodsort!(tracked::Tracked)
Sort the floes in tracked
by area in descending order.
sourceIceFloeTracker._adjust_histogram
— Method_adjust_histogram(masked_view, nbins, rblocks, cblocks, clip)
Perform adaptive histogram equalization to a masked image. To be invoked within imsharpen
.
Arguments
masked_view
: input image in truecolor
See imsharpen
for a description of the remaining arguments
sourceIceFloeTracker._bin9todec
— Method_bin9todec(v)
Get decimal representation of a bit vector v
with the leading bit at its leftmost posistion.
Example
julia> _bin9todec([0 0 0 0 0 0 0 0 0])
0
julia> _bin9todec([1 1 1 1 1 1 1 1 1])
-511
sourceIceFloeTracker._branch_filter
— Function_branch_filter(
+511
sourceIceFloeTracker._branch_filter
— Function_branch_filter(
img::AbstractArray{Bool},
func1::Function=branch_candidates_func,
-func2::Function=connected_background_count,)
Filter img
with _operator_lut
using lut1
and lut2
.
sourceIceFloeTracker._generate_se!
— Method_generate_se!(se)
Generate a structuring element by leveraging symmetry (mirroring and inverting) a given initial structuring element.
sourceIceFloeTracker._operator_lut
— Method_operator_lut(I, img, nhood, lut1, lut2)
Look up the neighborhood nhood
in lookup tables lut1
and lut2
.
Handles cases when the center of nhood
is on the edge of img
using data in I
.
sourceIceFloeTracker._swap_last_values!
— Method_swap_last_values!(df)
Swap the last two values of the area_mismatch
and corr
columns for each group in df
. For bookkeeping purposes for goodness of fit data during the tracking process.
sourceIceFloeTracker.absdiffmeanratio
— Methodabsdiffmeanratio(x, y)
Calculate the absolute difference between x
and y
divided by the mean of x
and y
.
sourceIceFloeTracker.add_padding
— Methodadd_padding(img, style)
Extrapolate the image img
according to the style
specifications type. Returns the extrapolated image.
Arguments
img
: Image to be padded.style
: A supported type (such as Pad
or Fill
) representing the extrapolation style. See the relevant documentation for details.
See also remove_padding
sourceIceFloeTracker.add_passtimes!
— Methodadd_passtimes!(props, passtimes)
Add a column passtime
to each DataFrame in props
containing the time of the image in which the floes were captured.
Arguments
props
: array of DataFrames containing floe properties.passtimes
: array of DateTime
objects containing the time of the image in which the floes were captured.
sourceIceFloeTracker.addfloemasks!
— Methodaddfloearrays(props::DataFrame, floeimg::BitMatrix)
Add a column to props
called floearray
containing the cropped floe masks from floeimg
.
sourceIceFloeTracker.addlatlon!
— Methodaddlatlon(pairedfloesdf::DataFrame, refimage::AbstractString)
Add columns latitude
, longitude
, and pixel coordinates x
, y
to pairedfloesdf
.
Arguments
pairedfloesdf
: dataframe containing floe tracking data.refimage
: path to reference image.
sourceIceFloeTracker.addmatch!
— Methodaddmatch!(matched_pairs, newmatch)
Add newmatch
to matched_pairs
.
sourceIceFloeTracker.adduuid!
— Methodadduuid!(props)
Assign a unique ID to each floe in each table of floe properties.
sourceIceFloeTracker.appendrows!
— Methodappendrows!(df::MatchingProps, props::T, ratios, idx::Int64, dist::Float64) where {T<:DataFrameRow}
Append a row to df.props
and df.ratios
with the values of props
and ratios
respectively.
sourceIceFloeTracker.apply_cloudmask
— Methodapply_cloudmask(false_color_image, cloudmask)
Zero out pixels containing clouds where clouds and ice are not discernable. Arguments should be of the same size.
Arguments
false_color_image
: reference image, e.g. corrected reflectance false color image bands [7,2,1] or grayscalecloudmask
: binary cloudmask with clouds = 0, else = 1
sourceIceFloeTracker.apply_landmask
— Methodapply_landmask(input_image, landmask_binary)
Zero out pixels in all channels of the input image using the binary landmask.
Arguments
input_image
: truecolor RGB imagelandmask_binary
: binary landmask with 1=land, 0=water/ice
sourceIceFloeTracker.atan2
— Methodatan2(y,x)
Wrapper of Base.atan
that returns the angle of vector (x,y) in the range [0, 2π).
sourceIceFloeTracker.binarize_landmask
— Methodbinarize_landmask(landmask_image)
Convert a 3-channel RGB land mask image to a 1-channel binary matrix; land = 0, ocean = 1.
Arguments
landmask_image
: RGB land mask image from fetchdata
sourceIceFloeTracker.branch
— Methodbranch(img::AbstractArray{Bool})
Find branch points in skeletonized image img
according to Definition 3 of [1].
[1] Arcelli, Carlo, and Gabriella Sanniti di Baja. "Skeletons of planar patterns." Machine Intelligence and Pattern Recognition. Vol. 19. North-Holland, 1996. 99-143.
sourceIceFloeTracker.branch_candidates_func
— Methodbranch_candidates_func(nhood)
Filter nhood
as candidate for branch point.
To be passed to the make_lut
function.
sourceIceFloeTracker.bridge
— Methodbridge(bw)
Set 0-valued pixels to 1 if they have two nonzero neighbors that are not connected. Note the following exceptions:
0 0 0 0 0 0 1 0 1 becomes 1 1 1 0 0 0 0 0 0
1 0 1 1 1 1 0 0 0 becomes 0 0 0 0 0 0 0 0 0
The same applies to all their corresponding rotations.
Examples
julia> bw = [0 0 0; 0 0 0; 1 0 1]
+func2::Function=connected_background_count,)
Filter img
with _operator_lut
using lut1
and lut2
.
sourceIceFloeTracker._generate_se!
— Method_generate_se!(se)
Generate a structuring element by leveraging symmetry (mirroring and inverting) a given initial structuring element.
sourceIceFloeTracker._operator_lut
— Method_operator_lut(I, img, nhood, lut1, lut2)
Look up the neighborhood nhood
in lookup tables lut1
and lut2
.
Handles cases when the center of nhood
is on the edge of img
using data in I
.
sourceIceFloeTracker._swap_last_values!
— Method_swap_last_values!(df)
Swap the last two values of the area_mismatch
and corr
columns for each group in df
. For bookkeeping purposes for goodness of fit data during the tracking process.
sourceIceFloeTracker.absdiffmeanratio
— Methodabsdiffmeanratio(x, y)
Calculate the absolute difference between x
and y
divided by the mean of x
and y
.
sourceIceFloeTracker.add_padding
— Methodadd_padding(img, style)
Extrapolate the image img
according to the style
specifications type. Returns the extrapolated image.
Arguments
img
: Image to be padded.style
: A supported type (such as Pad
or Fill
) representing the extrapolation style. See the relevant documentation for details.
See also remove_padding
sourceIceFloeTracker.add_passtimes!
— Methodadd_passtimes!(props, passtimes)
Add a column passtime
to each DataFrame in props
containing the time of the image in which the floes were captured.
Arguments
props
: array of DataFrames containing floe properties.passtimes
: array of DateTime
objects containing the time of the image in which the floes were captured.
sourceIceFloeTracker.addfloemasks!
— Methodaddfloearrays(props::DataFrame, floeimg::BitMatrix)
Add a column to props
called floearray
containing the cropped floe masks from floeimg
.
sourceIceFloeTracker.addlatlon!
— Methodaddlatlon(pairedfloesdf::DataFrame, refimage::AbstractString)
Add columns latitude
, longitude
, and pixel coordinates x
, y
to pairedfloesdf
.
Arguments
pairedfloesdf
: dataframe containing floe tracking data.refimage
: path to reference image.
sourceIceFloeTracker.addmatch!
— Methodaddmatch!(matched_pairs, newmatch)
Add newmatch
to matched_pairs
.
sourceIceFloeTracker.adduuid!
— Methodadduuid!(df)
Assign a unique ID to each floe in a table of floe properties.
sourceIceFloeTracker.adduuid!
— Methodadduuid!(dfs)
Assign a unique ID to each floe in an vector of tables of floe properties.
sourceIceFloeTracker.appendrows!
— Methodappendrows!(df::MatchingProps, props::T, ratios, idx::Int64, dist::Float64) where {T<:DataFrameRow}
Append a row to df.props
and df.ratios
with the values of props
and ratios
respectively.
sourceIceFloeTracker.apply_cloudmask
— Methodapply_cloudmask(false_color_image, cloudmask)
Zero out pixels containing clouds where clouds and ice are not discernable. Arguments should be of the same size.
Arguments
false_color_image
: reference image, e.g. corrected reflectance false color image bands [7,2,1] or grayscalecloudmask
: binary cloudmask with clouds = 0, else = 1
sourceIceFloeTracker.apply_landmask
— Methodapply_landmask(input_image, landmask_binary)
Zero out pixels in all channels of the input image using the binary landmask.
Arguments
input_image
: truecolor RGB imagelandmask_binary
: binary landmask with 1=land, 0=water/ice
sourceIceFloeTracker.atan2
— Methodatan2(y,x)
Wrapper of Base.atan
that returns the angle of vector (x,y) in the range [0, 2π).
sourceIceFloeTracker.binarize_landmask
— Methodbinarize_landmask(landmask_image)
Convert a 3-channel RGB land mask image to a 1-channel binary matrix; land = 0, ocean = 1.
Arguments
landmask_image
: RGB land mask image from fetchdata
sourceIceFloeTracker.branch
— Methodbranch(img::AbstractArray{Bool})
Find branch points in skeletonized image img
according to Definition 3 of [1].
[1] Arcelli, Carlo, and Gabriella Sanniti di Baja. "Skeletons of planar patterns." Machine Intelligence and Pattern Recognition. Vol. 19. North-Holland, 1996. 99-143.
sourceIceFloeTracker.branch_candidates_func
— Methodbranch_candidates_func(nhood)
Filter nhood
as candidate for branch point.
To be passed to the make_lut
function.
sourceIceFloeTracker.bridge
— Methodbridge(bw)
Set 0-valued pixels to 1 if they have two nonzero neighbors that are not connected. Note the following exceptions:
0 0 0 0 0 0 1 0 1 becomes 1 1 1 0 0 0 0 0 0
1 0 1 1 1 1 0 0 0 becomes 0 0 0 0 0 0 0 0 0
The same applies to all their corresponding rotations.
Examples
julia> bw = [0 0 0; 0 0 0; 1 0 1]
3×3 Matrix{Int64}:
0 0 0
0 0 0
@@ -40,7 +40,7 @@
3×3 BitMatrix:
1 1 1
1 1 1
- 1 1 1
sourceIceFloeTracker.bump_tile
— Methodbump_tile(tile::Tuple{UnitRange{Int64}, UnitRange{Int64}}, dims::Tuple{Int,Int})::Tuple{UnitRange{Int}, UnitRange{Int}}
Adjust the tile dimensions by adding extra rows and columns.
Arguments
tile::Tuple{Int,Int,Int,Int}
: A tuple representing the tile dimensions (a, b, c, d).dims::Tuple{Int,Int}
: A tuple representing the extra rows and columns to add (extrarows, extracols).
Returns
Tuple{UnitRange{Int}, UnitRange{Int}}
: A tuple of ranges representing the new tile dimensions.
Examples
```julia julia> bump_tile((1:3, 1:4), (1, 1)) (1:4, 1:5)
sourceIceFloeTracker.bwareamaxfilt
— Functionbwareamaxfilt(bwimg::AbstractArray{Bool}, conn)
Filter the smaller (by area) connected components in bwimg
keeping the (assumed unique) largest.
Uses 8-pixel connectivity by default (conn=8
). Use conn=4
for 4-pixel connectivity.
sourceIceFloeTracker.bwareamaxfilt!
— Functionsource IceFloeTracker.bwdist
— Methodbwdist(bwimg)
Distance transform for binary image bwdist
.
sourceIceFloeTracker.bwperim
— Methodbwperim(bwimg)
Locate the pixels at the boundary of objects in an binary image bwimg
using 8-pixel connectivity.
Arguments
bwimg
: Binary (black/white – 1/0) image
Examples
julia> A = zeros(Bool, 13, 16); A[2:6, 2:6] .= 1; A[4:8, 7:10] .= 1; A[10:12,13:15] .= 1; A[10:12,3:6] .= 1;
+ 1 1 1
sourceIceFloeTracker.bump_tile
— Methodbump_tile(tile::Tuple{UnitRange{Int64}, UnitRange{Int64}}, dims::Tuple{Int,Int})::Tuple{UnitRange{Int}, UnitRange{Int}}
Adjust the tile dimensions by adding extra rows and columns.
Arguments
tile::Tuple{Int,Int,Int,Int}
: A tuple representing the tile dimensions (a, b, c, d).dims::Tuple{Int,Int}
: A tuple representing the extra rows and columns to add (extrarows, extracols).
Returns
Tuple{UnitRange{Int}, UnitRange{Int}}
: A tuple of ranges representing the new tile dimensions.
Examples
```julia julia> bump_tile((1:3, 1:4), (1, 1)) (1:4, 1:5)
sourceIceFloeTracker.bwareamaxfilt
— Functionbwareamaxfilt(bwimg::AbstractArray{Bool}, conn)
Filter the smaller (by area) connected components in bwimg
keeping the (assumed unique) largest.
Uses 8-pixel connectivity by default (conn=8
). Use conn=4
for 4-pixel connectivity.
sourceIceFloeTracker.bwareamaxfilt!
— Functionsource IceFloeTracker.bwdist
— Methodbwdist(bwimg)
Distance transform for binary image bwdist
.
sourceIceFloeTracker.bwperim
— Methodbwperim(bwimg)
Locate the pixels at the boundary of objects in an binary image bwimg
using 8-pixel connectivity.
Arguments
bwimg
: Binary (black/white – 1/0) image
Examples
julia> A = zeros(Bool, 13, 16); A[2:6, 2:6] .= 1; A[4:8, 7:10] .= 1; A[10:12,13:15] .= 1; A[10:12,3:6] .= 1;
julia> A
13×16 Matrix{Bool}:
@@ -72,7 +72,7 @@
0 0 1 1 1 1 0 0 0 0 0 0 1 1 1 0
0 0 1 0 0 1 0 0 0 0 0 0 1 0 1 0
0 0 1 1 1 1 0 0 0 0 0 0 1 1 1 0
- 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
sourceIceFloeTracker.bwtraceboundary
— Methodbwtraceboundary(image::Union{Matrix{Int64},Matrix{Float64},T};
+ 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
sourceIceFloeTracker.bwtraceboundary
— Methodbwtraceboundary(image::Union{Matrix{Int64},Matrix{Float64},T};
P0::Union{Tuple{Int,Int},CartesianIndex{2},Nothing}=nothing,
closed::Bool=true) where T<:AbstractMatrix{Bool}
Trace the boundary of objects in image
Background pixels are represented as zero. The algorithm traces the boundary counterclockwise and an initial point P0
can be specified. If more than one boundary is detected and an initial point is provided, the boundary that contains this point is returned as a vector of CartesianIndex types. Otherwise an array of vectors is returned with all the detected boundaries in image
.
Arguments
image
: image, preferably binary with one single object, whose objects' boundaries are to be traced.P0
: initial point of a target boundary.closed
: if true
(default) makes the inital point of a boundary equal to the last point.
Example
julia> A = zeros(Int, 13, 16); A[2:6, 2:6] .= 1; A[4:8, 7:10] .= 1; A[10:12,13:15] .= 1; A[10:12,3:6] .= 1;
@@ -104,20 +104,20 @@
CartesianIndex(11, 15)
CartesianIndex(10, 15)
CartesianIndex(10, 14)
- CartesianIndex(10, 13)
sourceIceFloeTracker.callmatchcorr
— Methodcallmatchcorr(conditions)
Condition to decide whether match_corr should be called.
sourceIceFloeTracker.check_fname
— Functioncheck_fname(fname)
Checks fname
does not exist in current directory; throws an assertion if this condition is false.
Arguments
fname
: String object or Symbol to a reference to a String representing a path.
sourceIceFloeTracker.clockwise
— MethodMake a clockwise turn from the dir
direction
sourceIceFloeTracker.compute_ratios
— Methodcompute_ratios((props_day1, r), (props_day2,s))
Compute the ratios of the floe properties between the r
th floe in props_day1
and the s
th floe in props_day2
. Return a tuple of the ratios.
Arguments
props_day1
: floe properties for day 1r
: index of floe in props_day1
props_day2
: floe properties for day 2s
: index of floe in props_day2
sourceIceFloeTracker.compute_ratios_conditions
— Methodcompute_ratios_conditions((props_day1, r), (props_day2, s), delta_time, t)
Compute the conditions for a match between the r
th floe in props_day1
and the s
th floe in props_day2
. Return a tuple of the conditions.
Arguments
props_day1
: floe properties for day 1r
: index of floe in props_day1
props_day2
: floe properties for day 2s
: index of floe in props_day2
delta_time
: time elapsed from image day 1 to image day 2t
: tuple of thresholds for elapsed time and distance. See pair_floes
for details.
sourceIceFloeTracker.conditional_histeq
— Functionconditional_histeq(
+ CartesianIndex(10, 13)
sourceIceFloeTracker.callmatchcorr
— Methodcallmatchcorr(conditions)
Condition to decide whether match_corr should be called.
sourceIceFloeTracker.check_fname
— Functioncheck_fname(fname)
Checks fname
does not exist in current directory; throws an assertion if this condition is false.
Arguments
fname
: String object or Symbol to a reference to a String representing a path.
sourceIceFloeTracker.clockwise
— MethodMake a clockwise turn from the dir
direction
sourceIceFloeTracker.compute_ratios
— Methodcompute_ratios((props_day1, r), (props_day2,s))
Compute the ratios of the floe properties between the r
th floe in props_day1
and the s
th floe in props_day2
. Return a tuple of the ratios.
Arguments
props_day1
: floe properties for day 1r
: index of floe in props_day1
props_day2
: floe properties for day 2s
: index of floe in props_day2
sourceIceFloeTracker.compute_ratios_conditions
— Methodcompute_ratios_conditions((props_day1, r), (props_day2, s), delta_time, t)
Compute the conditions for a match between the r
th floe in props_day1
and the s
th floe in props_day2
. Return a tuple of the conditions.
Arguments
props_day1
: floe properties for day 1r
: index of floe in props_day1
props_day2
: floe properties for day 2s
: index of floe in props_day2
delta_time
: time elapsed from image day 1 to image day 2t
: tuple of thresholds for elapsed time and distance. See pair_floes
for details.
sourceIceFloeTracker.conditional_histeq
— Functionconditional_histeq(
true_color_image,
clouds_red,
side_length::Int,
entropy_threshold::AbstractFloat=4.0,
white_threshold::AbstractFloat=25.5,
-white_fraction_threshold::AbstractFloat=0.4,
)
Performs conditional histogram equalization on a true color image using tiles of approximately sidelength size side_length
. If a perfect tiling is not possible, the tiling on the egde of the image is adjusted to ensure that the tiles are as close to side_length
as possible. See get_tiles(array, side_length)
for more details.
sourceIceFloeTracker.conditional_histeq
— Functionconditional_histeq(
+white_fraction_threshold::AbstractFloat=0.4,
)
Performs conditional histogram equalization on a true color image using tiles of approximately sidelength size side_length
. If a perfect tiling is not possible, the tiling on the egde of the image is adjusted to ensure that the tiles are as close to side_length
as possible. See get_tiles(array, side_length)
for more details.
sourceIceFloeTracker.conditional_histeq
— Functionconditional_histeq(
true_color_image,
clouds_red,
rblocks::Int,
cblocks::Int,
entropy_threshold::AbstractFloat=4.0,
white_threshold::AbstractFloat=25.5,
-white_fraction_threshold::AbstractFloat=0.4,
)
Performs conditional histogram equalization on a true color image.
Arguments
true_color_image
: The true color image to be equalized.clouds_red
: The land/cloud masked red channel of the false color image.rblocks
: The number of row-blocks to divide the image into for histogram equalization. Default is 8.cblocks
: The number of column-blocks to divide the image into for histogram equalization. Default is 6.entropy_threshold
: The entropy threshold used to determine if a block should be equalized. Default is 4.0.white_threshold
: The white threshold used to determine if a pixel should be considered white. Default is 25.5.white_fraction_threshold
: The white fraction threshold used to determine if a block should be equalized. Default is 0.4.
Returns
The equalized true color image.
sourceIceFloeTracker.connected_background_count
— Methodconnected_background_count(nhood)
Second lut generator for neighbor transform with diamond strel (4-neighborhood).
To be passed to the make_lut
function.
sourceIceFloeTracker.consolidate_matched_pairs
— Methodconsolidate_matched_pairs(matched_pairs::MatchedPairs)
Consolidate the floe properties and similarity ratios of the matched pairs in matched_pairs
into a single dataframe. Return the consolidated dataframe. Used in iteration 0
.
sourceIceFloeTracker.convertcentroid!
— Methodconvertcentroid!(propdf, latlondata, colstodrop)
Convert the centroid coordinates from row and column to latitude and longitude dropping unwanted columns specified in colstodrop
for the output data structure. Addionally, add columns x
and y
with the pixel coordinates of the centroid.
sourceIceFloeTracker.converttounits!
— Methodconverttounits!(propdf, latlondata, colstodrop)
Convert the floe properties from pixels to kilometers and square kilometers where appropiate. Also drop the columns specified in colstodrop
.
sourceIceFloeTracker.corr
— Methodcorr(f1,f2)
Return the correlation between the psi-s curves p1
and p2
.
sourceIceFloeTracker.corr
— Methodcorr(f1,f2)
Return the normalized cross-correlation between the psi-s curves p1
and p2
.
sourceIceFloeTracker.counterclockwise
— MethodMake a counterclockwise turn from the dir
direction
sourceIceFloeTracker.create_cloudmask
— Methodcreate_cloudmask(false_color_image; prelim_threshold, band_7_threshold, band_2_threshold, ratio_lower, ratio_upper)
Convert a 3-channel false color reflectance image to a 1-channel binary matrix; clouds = 0, else = 1. Default thresholds are defined in the published Ice Floe Tracker article: Remote Sensing of the Environment 234 (2019) 111406.
Arguments
false_color_image
: corrected reflectance false color image - bands [7,2,1]prelim_threshold
: threshold value used to identify clouds in band 7, N0f8(RGB intensity/255)band_7_threshold
: threshold value used to identify cloud-ice in band 7, N0f8(RGB intensity/255)band_2_threshold
: threshold value used to identify cloud-ice in band 2, N0f8(RGB intensity/255)ratio_lower
: threshold value used to set lower ratio of cloud-ice in bands 7 and 2ratio_upper
: threshold value used to set upper ratio of cloud-ice in bands 7 and 2
sourceIceFloeTracker.create_landmask
— Methodcreate_landmask(landmask_image, struct_elem, fill_value_lower, fill_value_upper)
Convert a 3-channel RGB land mask image to a 1-channel binary matrix, and use a structuring element to extend a buffer to mask complex coastal features. In the resulting mask, land = 0 and ocean = 1. Returns a named tuple with the dilated and non-dilated landmask.
Arguments
landmask_image
: RGB land mask image from fetchdata
struct_elem
: structuring element for dilation (optional)fill_value_lower
: fill holes having at least these many pixels (optional)fill_value_upper
: fill holes having at most these many pixels (optional)
sourceIceFloeTracker.cropfloe
— Methodcropfloe(floesimg, props, i)
Crops the floe delimited by the bounding box data in props
at index i
from the floe image floesimg
.
If the dataframe has bounding box data min_row
, min_col
, max_row
, max_col
, but no label
, then returns the largest contiguous component.
If the dataframe has bounding box data min_row
, min_col
, max_row
, max_col
, and a label
, then returns the component with the label. In this case, floesimg
must be an Array{Int}.
If the dataframe has only a label
and no bounding box data, then returns the component with the label, padded by one cell of zeroes on all sides. In this case, floesimg
must be an Array{Int}.
sourceIceFloeTracker.cropfloe
— Methodcropfloe(floesimg, min_row, min_col, max_row, max_col)
Crops the floe delimited by min_row
, min_col
, max_row
, max_col
, from the floe image floesimg
.
sourceIceFloeTracker.cropfloe
— Methodcropfloe(floesimg, min_row, min_col, max_row, max_col, label)
Crops the floe from floesimg
with the label label
, returning the region bounded by min_row
, min_col
, max_row
, max_col
, and converting to a BitMatrix.
sourceIceFloeTracker.crosscorr
— Methodr, lags = crosscorr(u::Vector{T},
+white_fraction_threshold::AbstractFloat=0.4,
)
Performs conditional histogram equalization on a true color image.
Arguments
true_color_image
: The true color image to be equalized.clouds_red
: The land/cloud masked red channel of the false color image.rblocks
: The number of row-blocks to divide the image into for histogram equalization. Default is 8.cblocks
: The number of column-blocks to divide the image into for histogram equalization. Default is 6.entropy_threshold
: The entropy threshold used to determine if a block should be equalized. Default is 4.0.white_threshold
: The white threshold used to determine if a pixel should be considered white. Default is 25.5.white_fraction_threshold
: The white fraction threshold used to determine if a block should be equalized. Default is 0.4.
Returns
The equalized true color image.
sourceIceFloeTracker.connected_background_count
— Methodconnected_background_count(nhood)
Second lut generator for neighbor transform with diamond strel (4-neighborhood).
To be passed to the make_lut
function.
sourceIceFloeTracker.consolidate_matched_pairs
— Methodconsolidate_matched_pairs(matched_pairs::MatchedPairs)
Consolidate the floe properties and similarity ratios of the matched pairs in matched_pairs
into a single dataframe. Return the consolidated dataframe. Used in iteration 0
.
sourceIceFloeTracker.convertcentroid!
— Methodconvertcentroid!(propdf, latlondata, colstodrop)
Convert the centroid coordinates from row and column to latitude and longitude dropping unwanted columns specified in colstodrop
for the output data structure. Addionally, add columns x
and y
with the pixel coordinates of the centroid.
sourceIceFloeTracker.converttounits!
— Methodconverttounits!(propdf, latlondata, colstodrop)
Convert the floe properties from pixels to kilometers and square kilometers where appropiate. Also drop the columns specified in colstodrop
.
sourceIceFloeTracker.corr
— Methodcorr(f1,f2)
Return the correlation between the psi-s curves p1
and p2
.
sourceIceFloeTracker.corr
— Methodcorr(f1,f2)
Return the normalized cross-correlation between the psi-s curves p1
and p2
.
sourceIceFloeTracker.counterclockwise
— MethodMake a counterclockwise turn from the dir
direction
sourceIceFloeTracker.create_cloudmask
— Methodcreate_cloudmask(false_color_image; prelim_threshold, band_7_threshold, band_2_threshold, ratio_lower, ratio_upper)
Convert a 3-channel false color reflectance image to a 1-channel binary matrix; clouds = 0, else = 1. Default thresholds are defined in the published Ice Floe Tracker article: Remote Sensing of the Environment 234 (2019) 111406.
Arguments
false_color_image
: corrected reflectance false color image - bands [7,2,1]prelim_threshold
: threshold value used to identify clouds in band 7, N0f8(RGB intensity/255)band_7_threshold
: threshold value used to identify cloud-ice in band 7, N0f8(RGB intensity/255)band_2_threshold
: threshold value used to identify cloud-ice in band 2, N0f8(RGB intensity/255)ratio_lower
: threshold value used to set lower ratio of cloud-ice in bands 7 and 2ratio_upper
: threshold value used to set upper ratio of cloud-ice in bands 7 and 2
sourceIceFloeTracker.create_landmask
— Methodcreate_landmask(landmask_image, struct_elem, fill_value_lower, fill_value_upper)
Convert a 3-channel RGB land mask image to a 1-channel binary matrix, and use a structuring element to extend a buffer to mask complex coastal features. In the resulting mask, land = 0 and ocean = 1. Returns a named tuple with the dilated and non-dilated landmask.
Arguments
landmask_image
: RGB land mask image from fetchdata
struct_elem
: structuring element for dilation (optional)fill_value_lower
: fill holes having at least these many pixels (optional)fill_value_upper
: fill holes having at most these many pixels (optional)
sourceIceFloeTracker.cropfloe
— Methodcropfloe(floesimg, props, i)
Crops the floe delimited by the bounding box data in props
at index i
from the floe image floesimg
.
If the dataframe has bounding box data min_row
, min_col
, max_row
, max_col
, but no label
, then returns the largest contiguous component.
If the dataframe has bounding box data min_row
, min_col
, max_row
, max_col
, and a label
, then returns the component with the label. In this case, floesimg
must be an Array{Int}.
If the dataframe has only a label
and no bounding box data, then returns the component with the label, padded by one cell of zeroes on all sides. In this case, floesimg
must be an Array{Int}.
sourceIceFloeTracker.cropfloe
— Methodcropfloe(floesimg, min_row, min_col, max_row, max_col)
Crops the floe delimited by min_row
, min_col
, max_row
, max_col
, from the floe image floesimg
.
sourceIceFloeTracker.cropfloe
— Methodcropfloe(floesimg, min_row, min_col, max_row, max_col, label)
Crops the floe from floesimg
with the label label
, returning the region bounded by min_row
, min_col
, max_row
, max_col
, and converting to a BitMatrix.
sourceIceFloeTracker.crosscorr
— Methodr, lags = crosscorr(u::Vector{T},
v::Vector{T};
normalize::Bool=false,
padmode::Symbol=:longest)
Wrapper of DSP.xcorr with normalization (see https://docs.juliadsp.org/stable/convolutions/#DSP.xcorr)
Returns the pair (r, lags)
with the cross correlation scores r
and corresponding lags
according to padmode
.
Arguments
u,v
: real vectors which could have unequal length.normalize
: return normalized correlation scores (false
by default).padmode
: either :longest
(default) or :none
to control padding of shorter vector with zeros.
Examples
The example below builds two vectors, one a shifted version of the other, and computes various cross correlation scores.
julia> n = 1:5;
@@ -180,7 +180,7 @@
0.477211 4.0
0.229061 5.0
0.105402 6.0
-0.0411191 7.0
sourceIceFloeTracker.deleteallbut!
— Methoddeleteallbut!(matched_pairs, idxs, keeper)
Delete all rows in matched_pairs
except for the row with index keeper
in idxs
.
sourceIceFloeTracker.detect_move!
— MethodWorkhorse function: Get all pixel coords for detected border.
sourceIceFloeTracker.discriminate_ice_water
— Methoddiscriminate_ice_water(
+0.0411191 7.0
sourceIceFloeTracker.deleteallbut!
— Methoddeleteallbut!(matched_pairs, idxs, keeper)
Delete all rows in matched_pairs
except for the row with index keeper
in idxs
.
sourceIceFloeTracker.detect_move!
— MethodWorkhorse function: Get all pixel coords for detected border.
sourceIceFloeTracker.discriminate_ice_water
— Methoddiscriminate_ice_water(
falsecolor_image::Matrix{RGB{Float64}},
normalized_image::Matrix{Gray{Float64}},
landmask_bitmatrix::T,
@@ -195,11 +195,11 @@
st_dev_thresh_upper::Float64=Float64(98.9 / 255),
clouds_ratio_threshold::Float64=0.02,
differ_threshold::Float64=0.6,
-nbins::Real=155
)
Generates an image with ice floes apparent after filtering and combining previously processed versions of falsecolor and truecolor images from the same region of interest. Returns an image ready for segmentation to isolate floes.
Note: This function mutates the landmask object to avoid unnecessary memory allocation. If you need the original landmask, make a copy before passing it to this function. Example: discriminate_ice_water(falsecolor_image, normalized_image, copy(landmask_bitmatrix), cloudmask_bitmatrix)
Arguments
falsecolor_image
: input image in false color reflectancenormalized_image
: normalized version of true color imagelandmask_bitmatrix
: landmask for region of interestcloudmask_bitmatrix
: cloudmask for region of interestfloes_threshold
: heuristic applied to original false color imagemask_clouds_lower
: lower heuristic applied to mask out cloudsmask_clouds_upper
: upper heuristic applied to mask out cloudskurt_thresh_lower
: lower heuristic used to set pixel value threshold based on kurtosis in histogramkurt_thresh_upper
: upper heuristic used to set pixel value threshold based on kurtosis in histogramskew_thresh
: heuristic used to set pixel value threshold based on skewness in histogramst_dev_thresh_lower
: lower heuristic used to set pixel value threshold based on standard deviation in histogramst_dev_thresh_upper
: upper heuristic used to set pixel value threshold based on standard deviation in histogramclouds2_threshold
: heuristic used to set pixel value threshold based on ratio of cloudsdiffer_threshold
: heuristic used to calculate proportional intensity in histogramnbins
: number of bins during histogram build
sourceIceFloeTracker.dist
— Methoddist(p1, p2)
Return the distance between the points p1
and p2
.
sourceIceFloeTracker.filt_except_label!
— Methodfilt_except_label!(labeled_arr::Array{Int64, 2}, label::Int64)
In-place version of filt_except_label
.
See also filt_except_label
sourceIceFloeTracker.filt_except_label
— Methodfilt_except_label(labeled_arr::Array{Int64, 2}, label::Int64)
Make 0 all values in labeled_arr
that are not equal to label
.
See also filt_except_label!
sourceIceFloeTracker.find_floe_matches
— Methodfind_floe_matches(
+nbins::Real=155
)
Generates an image with ice floes apparent after filtering and combining previously processed versions of falsecolor and truecolor images from the same region of interest. Returns an image ready for segmentation to isolate floes.
Note: This function mutates the landmask object to avoid unnecessary memory allocation. If you need the original landmask, make a copy before passing it to this function. Example: discriminate_ice_water(falsecolor_image, normalized_image, copy(landmask_bitmatrix), cloudmask_bitmatrix)
Arguments
falsecolor_image
: input image in false color reflectancenormalized_image
: normalized version of true color imagelandmask_bitmatrix
: landmask for region of interestcloudmask_bitmatrix
: cloudmask for region of interestfloes_threshold
: heuristic applied to original false color imagemask_clouds_lower
: lower heuristic applied to mask out cloudsmask_clouds_upper
: upper heuristic applied to mask out cloudskurt_thresh_lower
: lower heuristic used to set pixel value threshold based on kurtosis in histogramkurt_thresh_upper
: upper heuristic used to set pixel value threshold based on kurtosis in histogramskew_thresh
: heuristic used to set pixel value threshold based on skewness in histogramst_dev_thresh_lower
: lower heuristic used to set pixel value threshold based on standard deviation in histogramst_dev_thresh_upper
: upper heuristic used to set pixel value threshold based on standard deviation in histogramclouds2_threshold
: heuristic used to set pixel value threshold based on ratio of cloudsdiffer_threshold
: heuristic used to calculate proportional intensity in histogramnbins
: number of bins during histogram build
sourceIceFloeTracker.dist
— Methoddist(p1, p2)
Return the distance between the points p1
and p2
.
sourceIceFloeTracker.filt_except_label!
— Methodfilt_except_label!(labeled_arr::Array{Int64, 2}, label::Int64)
In-place version of filt_except_label
.
See also filt_except_label
sourceIceFloeTracker.filt_except_label
— Methodfilt_except_label(labeled_arr::Array{Int64, 2}, label::Int64)
Make 0 all values in labeled_arr
that are not equal to label
.
See also filt_except_label!
sourceIceFloeTracker.find_floe_matches
— Methodfind_floe_matches(
tracked,
candidate_props,
condition_thresholds,
-mc_thresholds
)
Find matches for floes in tracked
from floes in candidate_props
.
Arguments
tracked
: dataframe containing floe trajectories.candidate_props
: dataframe containing floe candidate properties.condition_thresholds
: thresholds for deciding whether to match floe i
from tracked to floe j from candidate_props
mc_thresholds
: thresholds for area mismatch and psi-s shape correlation
sourceIceFloeTracker.find_ice_labels
— Methodfind_ice_labels(falsecolor_image, landmask; band_7_threshold, band_2_threshold, band_1_threshold, band_7_relaxed_threshold, band_1_relaxed_threshold, possible_ice_threshold)
Locate the pixels of likely ice from false color reflectance image. Returns a binary mask with ice floes contrasted from background. Default thresholds are defined in the published Ice Floe Tracker article: Remote Sensing of the Environment 234 (2019) 111406.
Arguments
falsecolor_image
: corrected reflectance false color image - bands [7,2,1]landmask
: bitmatrix landmask for region of interestband_7_threshold
: threshold value used to identify ice in band 7, N0f8(RGB intensity/255)band_2_threshold
: threshold value used to identify ice in band 2, N0f8(RGB intensity/255)band_1_threshold
: threshold value used to identify ice in band 2, N0f8(RGB intensity/255)band_7_relaxed_threshold
: threshold value used to identify ice in band 7 if not found on first pass, N0f8(RGB intensity/255)band_1_relaxed_threshold
: threshold value used to identify ice in band 1 if not found on first pass, N0f8(RGB intensity/255)
sourceIceFloeTracker.find_reflectance_peaks
— Methodfind_reflectance_peaks(reflectance_channel, possible_ice_threshold;)
Find histogram peaks in single channels of a reflectance image and return the second greatest peak. If needed, edges can be returned as the first object from build_histogram
. Similarly, peak values can be returned as the second object from findmaxima
.
Arguments
reflectance_channel
: either band 2 or band 1 of false-color reflectance imagepossible_ice_threshold
: threshold value used to identify ice if not found on first or second pass
sourceIceFloeTracker.fixzeroindexing!
— Methodfixzeroindexing!(props::DataFrame, props_to_fix::Vector{T}) where T<:Union{Symbol,String}
Fix the zero-indexing of the props_to_fix
columns in props
by adding 1 to each element.
sourceIceFloeTracker.fname_ext_splice
— Methodfname_ext_splice(fname, ext)
Join "fname"
and "ext"
with '.'
.
sourceIceFloeTracker.fname_ext_split
— Methodfname_ext_split(fname)
Split "fname.ext"
into "fname"
and "ext"
.
sourceIceFloeTracker.from_to
— MethodGet index in Moore neigborhood representing the direction from the from
pixel coords to the to
pixel coords (see definition of dir_delta below).
Clockwise Moore neighborhood.
dir_delta = [CartesianIndex(-1, 0), CartesianIndex(-1, 1), CartesianIndex(0, 1), CartesianIndex(1, 1), CartesianIndex(1, 0), CartesianIndex(1, -1), CartesianIndex(0, -1), CartesianIndex(-1,-1)]
sourceIceFloeTracker.get_area_missed
— Methodget_area_missed(side_length::Int, dims::Tuple{Int,Int})::Float64
Calculate the proportion of the area that is not covered by tiles of a given side length.
Arguments
side_length::Int
: The side length of the tile.dims::Tuple{Int,Int}
: A tuple representing the dimensions (width, height).
Returns
Float64
: The proportion of the area that is not covered by the tiles.
Examples
``` julia> getareamissed(5, (10, 20)) 0.0
julia> getareamissed(7, (10, 20)) 0.51
sourceIceFloeTracker.get_areas
— Methodget_areas(labeled_arr::Array{T, 2})::Dict{T, Int} where T
Get the "areas" (count of pixels of a given label) of the connected components in labeled_arr
.
Return a dictionary with the frequency distribution: label => countoflabel.
sourceIceFloeTracker.get_brighten_mask
— Methodget_brighten_mask(equalized_gray_reconstructed_img, gamma_green)
Arguments
equalized_gray_reconstructed_img
: The equalized gray reconstructed image (uint8 in Matlab).gamma_green
: The gamma value for the green channel (also uint8).
Returns
Difference equalizedgrayreconstructedimg - gammagreen clamped between 0 and 255.
sourceIceFloeTracker.get_dt
— Methodget_dt(props1, r, props2, s)
Return the time difference between the r
th floe in props1
and the s
th floe in props2
in minutes.
sourceIceFloeTracker.get_final
— Functionget_final(img, label, segment_mask, se_erosion, se_dilation)
Final processing following the tiling workflow.
Arguments
img
: The input image.label
: Mode of most common label in the findicelabels workflow.segment_mask
: The segment mask.se_erosion
: structuring element for erosion.se_dilation
: structuring element for dilation.apply_segment_mask=true
: Whether to filter img
the segment mask.
sourceIceFloeTracker.get_ice_masks
— Methodget_ice_masks(
+mc_thresholds
)
Find matches for floes in tracked
from floes in candidate_props
.
Arguments
tracked
: dataframe containing floe trajectories.candidate_props
: dataframe containing floe candidate properties.condition_thresholds
: thresholds for deciding whether to match floe i
from tracked to floe j from candidate_props
mc_thresholds
: thresholds for area mismatch and psi-s shape correlation
sourceIceFloeTracker.find_ice_labels
— Methodfind_ice_labels(falsecolor_image, landmask; band_7_threshold, band_2_threshold, band_1_threshold, band_7_relaxed_threshold, band_1_relaxed_threshold, possible_ice_threshold)
Locate the pixels of likely ice from false color reflectance image. Returns a binary mask with ice floes contrasted from background. Default thresholds are defined in the published Ice Floe Tracker article: Remote Sensing of the Environment 234 (2019) 111406.
Arguments
falsecolor_image
: corrected reflectance false color image - bands [7,2,1]landmask
: bitmatrix landmask for region of interestband_7_threshold
: threshold value used to identify ice in band 7, N0f8(RGB intensity/255)band_2_threshold
: threshold value used to identify ice in band 2, N0f8(RGB intensity/255)band_1_threshold
: threshold value used to identify ice in band 2, N0f8(RGB intensity/255)band_7_relaxed_threshold
: threshold value used to identify ice in band 7 if not found on first pass, N0f8(RGB intensity/255)band_1_relaxed_threshold
: threshold value used to identify ice in band 1 if not found on first pass, N0f8(RGB intensity/255)
sourceIceFloeTracker.find_reflectance_peaks
— Methodfind_reflectance_peaks(reflectance_channel, possible_ice_threshold;)
Find histogram peaks in single channels of a reflectance image and return the second greatest peak. If needed, edges can be returned as the first object from build_histogram
. Similarly, peak values can be returned as the second object from findmaxima
.
Arguments
reflectance_channel
: either band 2 or band 1 of false-color reflectance imagepossible_ice_threshold
: threshold value used to identify ice if not found on first or second pass
sourceIceFloeTracker.fixzeroindexing!
— Methodfixzeroindexing!(props::DataFrame, props_to_fix::Vector{T}) where T<:Union{Symbol,String}
Fix the zero-indexing of the props_to_fix
columns in props
by adding 1 to each element.
sourceIceFloeTracker.fname_ext_splice
— Methodfname_ext_splice(fname, ext)
Join "fname"
and "ext"
with '.'
.
sourceIceFloeTracker.fname_ext_split
— Methodfname_ext_split(fname)
Split "fname.ext"
into "fname"
and "ext"
.
sourceIceFloeTracker.from_to
— MethodGet index in Moore neigborhood representing the direction from the from
pixel coords to the to
pixel coords (see definition of dir_delta below).
Clockwise Moore neighborhood.
dir_delta = [CartesianIndex(-1, 0), CartesianIndex(-1, 1), CartesianIndex(0, 1), CartesianIndex(1, 1), CartesianIndex(1, 0), CartesianIndex(1, -1), CartesianIndex(0, -1), CartesianIndex(-1,-1)]
sourceIceFloeTracker.get_area_missed
— Methodget_area_missed(side_length::Int, dims::Tuple{Int,Int})::Float64
Calculate the proportion of the area that is not covered by tiles of a given side length.
Arguments
side_length::Int
: The side length of the tile.dims::Tuple{Int,Int}
: A tuple representing the dimensions (width, height).
Returns
Float64
: The proportion of the area that is not covered by the tiles.
Examples
``` julia> getareamissed(5, (10, 20)) 0.0
julia> getareamissed(7, (10, 20)) 0.51
sourceIceFloeTracker.get_areas
— Methodget_areas(labeled_arr::Array{T, 2})::Dict{T, Int} where T
Get the "areas" (count of pixels of a given label) of the connected components in labeled_arr
.
Return a dictionary with the frequency distribution: label => countoflabel.
sourceIceFloeTracker.get_brighten_mask
— Methodget_brighten_mask(equalized_gray_reconstructed_img, gamma_green)
Arguments
equalized_gray_reconstructed_img
: The equalized gray reconstructed image (uint8 in Matlab).gamma_green
: The gamma value for the green channel (also uint8).
Returns
Difference equalizedgrayreconstructedimg - gammagreen clamped between 0 and 255.
sourceIceFloeTracker.get_dt
— Methodget_dt(props1, r, props2, s)
Return the time difference between the r
th floe in props1
and the s
th floe in props2
in minutes.
sourceIceFloeTracker.get_final
— Functionget_final(img, label, segment_mask, se_erosion, se_dilation)
Final processing following the tiling workflow.
Arguments
img
: The input image.label
: Mode of most common label in the findicelabels workflow.segment_mask
: The segment mask.se_erosion
: structuring element for erosion.se_dilation
: structuring element for dilation.apply_segment_mask=true
: Whether to filter img
the segment mask.
sourceIceFloeTracker.get_ice_masks
— Methodget_ice_masks(
falsecolor_image,
morph_residue,
landmask,
@@ -213,10 +213,10 @@
possible_ice_threshold,
k,
factor
-)
Get the ice masks from the falsecolor image and morphological residue given a particular tiling configuration.
Arguments
falsecolor_image
: The falsecolor image.morph_residue
: The morphological residue image.landmask
: The landmask.tiles
: The tiles.binarize::Bool=true
: Whether to binarize the tiling.band_7_threshold=5
: The threshold for band 7.band_2_threshold=230
: The threshold for band 2.band_1_threshold=240
: The threshold for band 1.band_7_threshold_relaxed=10
: The relaxed threshold for band 7.band_1_threshold_relaxed=190
: The relaxed threshold for band 1.possible_ice_threshold=75
: The threshold for possible ice.k=3
: The number of clusters to use for k-means segmentation.factor=255
: normalization factor to convert images to uint8.
Returns
- A named tuple
(icemask, bin)
where:icemask
: The ice mask.bin
: The binarized tiling.
sourceIceFloeTracker.get_matches
— Methodget_matches(matched_pairs)
Return a dataframe with the properties and goodness ratios of the matched pairs (right-hand matches) in matched_pairs
. Used in iterations 1:end
.
sourceIceFloeTracker.get_max_label
— Methodget_max_label(d::Dict{Int, Int})
Get the label k
in dictionary d
for which d[k] is maximal.
sourceIceFloeTracker.get_optimal_tile_size
— Methodget_optimal_tile_size(l0::Int, dims::Tuple{Int,Int}) -> Int
Calculate the optimal tile size in the range [l0-1, l0+1] for the given size l0
and image dimensions dims
.
Description
This function computes the optimal tile size for tiling an area with given dimensions. It ensures that the initial tile size l0
is at least 2 and not larger than any of the given dimensions. The function evaluates candidate tile sizes and selects the one that minimizes the area missed by its corresponding tiling. In case of a tie, it prefers the larger tile size.
Example
julia> get_optimal_tile_size(3, (10, 7))
-2
sourceIceFloeTracker.get_rgb_channels
— Methodget_rgb_channels(img)
Get the RBC (Red, Blue, and Green) channels of an image.
Arguments
img
: The input image.
Returns
An m x n x 3 array the Red, Blue, and Green channels of the input image.
sourceIceFloeTracker.get_tile_dims
— Methodget_tile_dims(tile)
Calculate the dimensions of a tile.
Arguments
tile::Tuple{UnitRange{Int},UnitRange{Int}}
: A tuple representing the tile dimensions.
Returns
Tuple{Int,Int}
: A tuple representing the width and height of the tile.
Examples
```julia julia> gettiledims((1:3, 1:4)) (4, 3)
sourceIceFloeTracker.get_tile_meta
— Methodget_tile_meta(tile)
Extracts metadata from a given tile.
Arguments
tile
: A collection of tuples, where each tuple represents a coordinate pair.
Returns
- A tuple
(a, b, c, d)
where:a
: The first element of the first tuple in tile
.b
: The last element of the first tuple in tile
.c
: The first element of the last tuple in tile
.d
: The last element of the last tuple in tile
.
sourceIceFloeTracker.get_tiles
— Methodget_tiles(array, side_length)
Generate a collection of tiles from an array.
Unlike TileIterator
, the function adjusts the bottom and right edges of the tile matrix if they are smaller than half the tile size side_length
.
sourceIceFloeTracker.get_tiles
— Methodget_tiles(array, t::Tuple{Int,Int})
Generate a collection of tiles from an array.
The function adjusts the bottom and right edges of the tile matrix if they are smaller than half the tile sizes in t
.
sourceIceFloeTracker.get_trajectory_heads
— Methodget_trajectory_heads(pairs)
Return the last row (most recent member) of each group (trajectory) in pairs
as a dataframe.
This is used for getting the initial floe properties for the next day in search for new pairs.
sourceIceFloeTracker.get_unmatched
— Methodget_unmatched(props, matched)
Return the floes in props
that are not in matched
.
sourceIceFloeTracker.getbboxcolumns
— Methodgetbboxcolumns(props::DataFrame)
Return the column names of the bounding box columns in props
.
sourceIceFloeTracker.getbestmatchdata
— Methodgetbestmatchdata(idx, r, props_day1, matching_floes)
Collect the data for the best match between the r
th floe in props_day1
and the idx
th floe in matching_floes
. Return a tuple of the floe properties for day 1 and day 2 and the ratios.
sourceIceFloeTracker.getcentroid
— Methodgetcentroid(props_day::DataFrame, r)
Get the coordinates of the r
th floe in props_day
.
sourceIceFloeTracker.getcentroidcolumns
— Methodgetcentroidcolumns(props::DataFrame)
Returns the column names of the centroid columns in props
.
sourceIceFloeTracker.getcollisions
— Methodgetcollisions(matchedpairs)
Get nonunique rows in matchedpairs
.
sourceIceFloeTracker.getcollisionslocs
— Methodgetcollisionslocs(df)
Return a vector of tuples of the row and the index of the row in df
that has a collision with another row in df
.
sourceIceFloeTracker.getfit
— Methodgetfit(dims::Tuple{Int,Int}, side_length::Int)::Tuple{Int,Int}
Calculate how many tiles of a given side length fit into the given dimensions.
Arguments
dims::Tuple{Int,Int}
: A tuple representing the dimensions (width, height).side_length::Int
: The side length of the tile.
Returns
Tuple{Int,Int}
: A tuple representing the number of tiles that fit along each dimension.
Examples
``` julia> getfit((10, 20), 5) (2, 4)
julia> getfit((15, 25), 5) (3, 5)
sourceIceFloeTracker.getfloemasks
— Methodgetfloemasks(props::DataFrame, floeimg::BitMatrix)
Return a vector of cropped floe masks from floeimg
using the bounding box data in props
.
sourceIceFloeTracker.getidxmostminimumeverything
— Methodgetidxmostminimumeverything(ratiosdf)
Return the index of the row in ratiosdf
with the most minima across its columns. If there are multiple columns with the same minimum value, return the index of the first column with the minimum value. If ratiosdf
is empty, return NaN
.
sourceIceFloeTracker.getidxofrow
— Methodgetidxofrow(rw0, df)
Return the indices of the rows in df
that are equal to rw0
.
sourceIceFloeTracker.getpropsday1day2
— Methodgetpropsday1day2(properties, dayidx::Int64)
Return the floe properties for day dayidx
and day dayidx+1
.
sourceIceFloeTracker.getsizecomparability
— Methodgetsizecomparability(s1, s2)
Check if the size of two floes s1
and s2
are comparable. The size is defined as the product of the floe dimensions.
Arguments
s1
: size of floe 1s2
: size of floe 2
sourceIceFloeTracker.grad
— Methoddx, dy = grad(A::Matrix{<:Number})
Make gradient vector field for the set of points with coordinates in the rows of the matrix A
with x-coordinates down column 1 and y-coordinates down column 2. Return a tuple with dx
and dy
in that order.
sourceIceFloeTracker.grad
— Methoddx, dy = grad(x::Vector{<:Number}, y::Vector{<:Number})
Make gradient vector field for the set of points with coordinates in vectors x
and y
. Return a tuple with dx
and dy
in that order.
sourceIceFloeTracker.hbreak!
— Methodhbreak!(img::AbstractArray{Bool})
Inplace version of hbreak
. See also hbreak
.
sourceIceFloeTracker.hbreak
— Methodhbreak(img::AbstractArray{Bool})
Remove H-connected pixels in the binary image img
. See also hbreak!
for an inplace version of this function.
Examples
```jldoctest; setup = :(using IceFloeTracker)
julia> h1 = trues(3,3); h1[[1 3], 2] .= false; h1 3×3 BitMatrix: 1 0 1 1 1 1 1 0 1
julia> h2 = trues(3,3); h2[2, [1 3]] .= false; h2 3×3 BitMatrix: 1 1 1 0 1 0 1 1 1
julia> hbreak!(h1); h1 # modify h1 inplace 3×3 BitMatrix: 1 0 1 1 0 1 1 0 1
julia> hbreak(h2) 3×3 BitMatrix: 1 1 1 0 0 0 1 1 1
sourceIceFloeTracker.histeq
— Methodhisteq(img)
-histeq(img; nbins=64)
Histogram equalization of img
using nbins
bins.
sourceIceFloeTracker.imadjust
— Methodimadjust(img; low, high)
Adjust the contrast of an image using linear stretching. The image is normalized to [0, 1] and then stretched to the range [low, high].
Arguments
img
: The input image.low
: The lower bound of the stretched image. Default is 0.01.high
: The upper bound of the stretched image. Default is 0.99.
Returns
The contrast-adjusted image in the range [0, 255].
sourceIceFloeTracker.imbrighten
— Methodimbrighten(img, brighten_mask, bright_factor)
Brighten the image using a mask and a brightening factor.
Arguments
img
: The input image.brighten_mask
: A mask indicating the pixels to brighten.bright_factor
: The factor by which to brighten the pixels.
Returns
- The brightened image.
sourceIceFloeTracker.imextendedmin
— Functionimextendedmin(img)
Mimics MATLAB's imextendedmin function that computes the extended-minima transform, which is the regional minima of the H-minima transform. Regional minima are connected components of pixels with a constant intensity value. This function returns a transformed bitmatrix.
Arguments
img
: image objecth
: suppress minima below this depth thresholdconn
: neighborhood connectivity; in 2D 1 = 4-neighborhood and 2 = 8-neighborhood
sourceIceFloeTracker.imgradientmag
— Methodimgradientmag(img)
Compute the gradient magnitude of an image using the Sobel operator.
sourceIceFloeTracker.imhist
— Functionimhist(img, imgtype::AbstractString="uint8")
Compute the histogram of an image where each possible value is represented in the histogram. The function returns a tuple with the bins and counts of each bin.
Example
```jldoctest; setup = :(using IceFloeTracker) julia> img = [ 4 4 4 4 4 3 4 5 4 3 3 5 5 5 3 3 4 5 4 3 4 4 4 4 4 ]
julia> bins, heights = imhist(img);
julia> [bins[heights .> 0] heights[heights .>0]] # display only non-zero bins and heights 3×2 Matrix{Int64}: 3 6 4 14 5 5
sourceIceFloeTracker.impose_minima
— Methodimpose_minima(I::AbstractArray{T}, BW::AbstractArray{Bool}) where {T<:Integer}
Use morphological reconstruction to enforce minima on the input image I
at the positions where the binary mask BW
is non-zero.
It supports both integer and grayscale images using different implementations for each.
sourceIceFloeTracker.imregionalmin
— Functionimregionalmin(img, conn=2)
Compute the regional minima of the image img
using the connectivity conn
.
Returns a bitmatrix of the same size as img
with the regional minima.
Arguments
img
: Image objectconn
: Neighborhood connectivity; in 2D, 1 = 4-neighborhood and 2 = 8-neighborhood
sourceIceFloeTracker.imsharpen
— Functionimsharpen(truecolor_image, landmask_no_dilate, lambda, kappa, niters, nbins, rblocks, cblocks, clip, smoothing_param, intensity)
Sharpen truecolor_image
.
Arguments
truecolor_image
: input image in truecolorlandmask_no_dilate
: landmask for region of interestlambda
: speed of diffusion (0–0.25)kappa
: conduction coefficient for diffusion (25–100)niters
: number of iterations of diffusionnbins
: number of bins during histogram equalizationrblocks
: number of row blocks to divide input image during equalizationcblocks
: number of column blocks to divide input image during equalizationclip
: Thresholds for clipping histogram bins (0–1); values closer to one minimize contrast enhancement, values closer to zero maximize contrast enhancementsmoothing_param
: pixel radius for gaussian blurring (1–10)intensity
: amount of sharpening to perform
sourceIceFloeTracker.imsharpen_gray
— Methodimsharpen_gray(imgsharpened, landmask)
Apply landmask and return Gray type image in colorview for normalization.
sourceIceFloeTracker.isfloegoodmatch
— Methodisfloegoodmatch(conditions, mct, area_mismatch, corr)
Return true
if the floes are a good match as per the set thresholds. Return false
otherwise.
Arguments
conditions
: tuple of booleans for evaluating the conditionsmct
: tuple of thresholds for the match correlation testarea_mismatch
and corr
: values returned by match_corr
sourceIceFloeTracker.isincountourlist
— MethodCheck P
is in countour list if so return the index of the contour that contains P
, otherwise return false.
sourceIceFloeTracker.ismember
— Methodismember(df1,df2)
Return a boolean array indicating whether each row in df1
is a member of df2
.
sourceIceFloeTracker.kmeans_segmentation
— Functionkmeans_segmentation(gray_image, ice_labels;)
Apply k-means segmentation to a gray image to isolate a cluster group representing sea ice. Returns a binary image with ice segmented from background.
Arguments
gray_image
: output image from ice-water-discrimination.jl
or gray ice floe leads image in segmentation_f.jl
ice_labels
: vector if pixel coordinates output from find_ice_labels.jl
sourceIceFloeTracker.loadimg
— Methodloadimg(; dir::String, fname::String)
Load an image from dir
with filename fname
into a matrix of Float64
values. Returns the loaded image.
sourceIceFloeTracker.long_tracker
— Methodlong_tracker(props, condition_thresholds, mc_thresholds)
Track ice floes over multiple days.
Trajectories are built in two steps:
- Get pairs of floes in day 1 and day 2. Any unmatched floes, in both day 1 and day 2, become the "heads" of their respective trajectories.
- For each subsequent day, find pairs of floes for the current trajectory heads. Again, any unmatched floe in the new prop table starts a new trajectory.
Arguments
props::Vector{DataFrame}
: A vector of DataFrames, each containing ice floe properties for a single day. Each DataFrame must have the following columns:- "area"
- "min_row"
- "min_col"
- "max_row"
- "max_col"
- "row_centroid"
- "col_centroid"
- "convex_area"
- "majoraxislength"
- "minoraxislength"
- "orientation"
- "perimeter"
- "mask": 2D array of booleans
- "passtime": A timestamp for the floe
- "psi": the psi-s curve for the floe
- "uuid": a universally unique identifier for each segmented floe
condition_thresholds
: 3-tuple of thresholds (each a named tuple) for deciding whether to match floe i
from day k
to floe j from day k+1
mc_thresholds
: thresholds for area mismatch and psi-s shape correlation
Returns
A DataFrame with the above columns, plus two extra columns, "area_mismatch" and "corr", which are the area mismatch and correlation between a floe and the one that follows it in the trajectory. Trajectories are identified by a unique identifier, "uuid".
sourceIceFloeTracker.make_filename
— Methodmake_filename()
Makes default filename with timestamp.
sourceIceFloeTracker.make_hbreak_dict
— Methodmake_hbreak_dict()
Build dict with the two versions of an H-connected 3x3 neighboorhood.
h1 = [1 0 1 1 1 1 1 0 1]
h2 = [1 1 1 0 1 0 1 1 1]
sourceIceFloeTracker.make_lut
— Methodmake_lut(lutfunc::Function)
Generate lookup table (lut) for 3x3 neighborhoods according to lutfunc
.
sourceIceFloeTracker.make_psi_s
— Methodmake_psi_s(XY::Matrix{<:Number};rangeout::Bool=true,
-unwrap::Bool=true)
Alternate method of make_psi_s
accepting input vectors x
and y
as a 2-column matrix [x y]
in order to facillitate workflow (output from resample_boundary
).
sourceIceFloeTracker.make_psi_s
— Methodmake_psi_s(x::Vector{<:Number},
+)
Get the ice masks from the falsecolor image and morphological residue given a particular tiling configuration.
Arguments
falsecolor_image
: The falsecolor image.morph_residue
: The morphological residue image.landmask
: The landmask.tiles
: The tiles.binarize::Bool=true
: Whether to binarize the tiling.band_7_threshold=5
: The threshold for band 7.band_2_threshold=230
: The threshold for band 2.band_1_threshold=240
: The threshold for band 1.band_7_threshold_relaxed=10
: The relaxed threshold for band 7.band_1_threshold_relaxed=190
: The relaxed threshold for band 1.possible_ice_threshold=75
: The threshold for possible ice.k=3
: The number of clusters to use for k-means segmentation.factor=255
: normalization factor to convert images to uint8.
Returns
- A named tuple
(icemask, bin)
where:icemask
: The ice mask.bin
: The binarized tiling.
sourceIceFloeTracker.get_matches
— Methodget_matches(matched_pairs)
Return a dataframe with the properties and goodness ratios of the matched pairs (right-hand matches) in matched_pairs
. Used in iterations 1:end
.
sourceIceFloeTracker.get_max_label
— Methodget_max_label(d::Dict{Int, Int})
Get the label k
in dictionary d
for which d[k] is maximal.
sourceIceFloeTracker.get_optimal_tile_size
— Methodget_optimal_tile_size(l0::Int, dims::Tuple{Int,Int}) -> Int
Calculate the optimal tile size in the range [l0-1, l0+1] for the given size l0
and image dimensions dims
.
Description
This function computes the optimal tile size for tiling an area with given dimensions. It ensures that the initial tile size l0
is at least 2 and not larger than any of the given dimensions. The function evaluates candidate tile sizes and selects the one that minimizes the area missed by its corresponding tiling. In case of a tie, it prefers the larger tile size.
Example
julia> get_optimal_tile_size(3, (10, 7))
+2
sourceIceFloeTracker.get_rgb_channels
— Methodget_rgb_channels(img)
Get the RBC (Red, Blue, and Green) channels of an image.
Arguments
img
: The input image.
Returns
An m x n x 3 array the Red, Blue, and Green channels of the input image.
sourceIceFloeTracker.get_tile_dims
— Methodget_tile_dims(tile)
Calculate the dimensions of a tile.
Arguments
tile::Tuple{UnitRange{Int},UnitRange{Int}}
: A tuple representing the tile dimensions.
Returns
Tuple{Int,Int}
: A tuple representing the width and height of the tile.
Examples
```julia julia> gettiledims((1:3, 1:4)) (4, 3)
sourceIceFloeTracker.get_tile_meta
— Methodget_tile_meta(tile)
Extracts metadata from a given tile.
Arguments
tile
: A collection of tuples, where each tuple represents a coordinate pair.
Returns
- A tuple
(a, b, c, d)
where:a
: The first element of the first tuple in tile
.b
: The last element of the first tuple in tile
.c
: The first element of the last tuple in tile
.d
: The last element of the last tuple in tile
.
sourceIceFloeTracker.get_tiles
— Methodget_tiles(array, side_length)
Generate a collection of tiles from an array.
Unlike TileIterator
, the function adjusts the bottom and right edges of the tile matrix if they are smaller than half the tile size side_length
.
sourceIceFloeTracker.get_tiles
— Methodget_tiles(array, t::Tuple{Int,Int})
Generate a collection of tiles from an array.
The function adjusts the bottom and right edges of the tile matrix if they are smaller than half the tile sizes in t
.
sourceIceFloeTracker.get_trajectory_heads
— Methodget_trajectory_heads(pairs)
Return the last row (most recent member) of each group (trajectory) in pairs
as a dataframe.
This is used for getting the initial floe properties for the next day in search for new pairs.
sourceIceFloeTracker.get_unmatched
— Methodget_unmatched(props, matched)
Return the floes in props
that are not in matched
.
sourceIceFloeTracker.getbboxcolumns
— Methodgetbboxcolumns(props::DataFrame)
Return the column names of the bounding box columns in props
.
sourceIceFloeTracker.getbestmatchdata
— Methodgetbestmatchdata(idx, r, props_day1, matching_floes)
Collect the data for the best match between the r
th floe in props_day1
and the idx
th floe in matching_floes
. Return a tuple of the floe properties for day 1 and day 2 and the ratios.
sourceIceFloeTracker.getcentroid
— Methodgetcentroid(props_day::DataFrame, r)
Get the coordinates of the r
th floe in props_day
.
sourceIceFloeTracker.getcentroidcolumns
— Methodgetcentroidcolumns(props::DataFrame)
Returns the column names of the centroid columns in props
.
sourceIceFloeTracker.getcollisions
— Methodgetcollisions(matchedpairs)
Get nonunique rows in matchedpairs
.
sourceIceFloeTracker.getcollisionslocs
— Methodgetcollisionslocs(df)
Return a vector of tuples of the row and the index of the row in df
that has a collision with another row in df
.
sourceIceFloeTracker.getfit
— Methodgetfit(dims::Tuple{Int,Int}, side_length::Int)::Tuple{Int,Int}
Calculate how many tiles of a given side length fit into the given dimensions.
Arguments
dims::Tuple{Int,Int}
: A tuple representing the dimensions (width, height).side_length::Int
: The side length of the tile.
Returns
Tuple{Int,Int}
: A tuple representing the number of tiles that fit along each dimension.
Examples
``` julia> getfit((10, 20), 5) (2, 4)
julia> getfit((15, 25), 5) (3, 5)
sourceIceFloeTracker.getfloemasks
— Methodgetfloemasks(props::DataFrame, floeimg::BitMatrix)
Return a vector of cropped floe masks from floeimg
using the bounding box data in props
.
sourceIceFloeTracker.getidxmostminimumeverything
— Methodgetidxmostminimumeverything(ratiosdf)
Return the index of the row in ratiosdf
with the most minima across its columns. If there are multiple columns with the same minimum value, return the index of the first column with the minimum value. If ratiosdf
is empty, return NaN
.
sourceIceFloeTracker.getidxofrow
— Methodgetidxofrow(rw0, df)
Return the indices of the rows in df
that are equal to rw0
.
sourceIceFloeTracker.getpropsday1day2
— Methodgetpropsday1day2(properties, dayidx::Int64)
Return the floe properties for day dayidx
and day dayidx+1
.
sourceIceFloeTracker.getsizecomparability
— Methodgetsizecomparability(s1, s2)
Check if the size of two floes s1
and s2
are comparable. The size is defined as the product of the floe dimensions.
Arguments
s1
: size of floe 1s2
: size of floe 2
sourceIceFloeTracker.grad
— Methoddx, dy = grad(A::Matrix{<:Number})
Make gradient vector field for the set of points with coordinates in the rows of the matrix A
with x-coordinates down column 1 and y-coordinates down column 2. Return a tuple with dx
and dy
in that order.
sourceIceFloeTracker.grad
— Methoddx, dy = grad(x::Vector{<:Number}, y::Vector{<:Number})
Make gradient vector field for the set of points with coordinates in vectors x
and y
. Return a tuple with dx
and dy
in that order.
sourceIceFloeTracker.hbreak!
— Methodhbreak!(img::AbstractArray{Bool})
Inplace version of hbreak
. See also hbreak
.
sourceIceFloeTracker.hbreak
— Methodhbreak(img::AbstractArray{Bool})
Remove H-connected pixels in the binary image img
. See also hbreak!
for an inplace version of this function.
Examples
```jldoctest; setup = :(using IceFloeTracker)
julia> h1 = trues(3,3); h1[[1 3], 2] .= false; h1 3×3 BitMatrix: 1 0 1 1 1 1 1 0 1
julia> h2 = trues(3,3); h2[2, [1 3]] .= false; h2 3×3 BitMatrix: 1 1 1 0 1 0 1 1 1
julia> hbreak!(h1); h1 # modify h1 inplace 3×3 BitMatrix: 1 0 1 1 0 1 1 0 1
julia> hbreak(h2) 3×3 BitMatrix: 1 1 1 0 0 0 1 1 1
sourceIceFloeTracker.histeq
— Methodhisteq(img)
+histeq(img; nbins=64)
Histogram equalization of img
using nbins
bins.
sourceIceFloeTracker.imadjust
— Methodimadjust(img; low, high)
Adjust the contrast of an image using linear stretching. The image is normalized to [0, 1] and then stretched to the range [low, high].
Arguments
img
: The input image.low
: The lower bound of the stretched image. Default is 0.01.high
: The upper bound of the stretched image. Default is 0.99.
Returns
The contrast-adjusted image in the range [0, 255].
sourceIceFloeTracker.imbrighten
— Methodimbrighten(img, brighten_mask, bright_factor)
Brighten the image using a mask and a brightening factor.
Arguments
img
: The input image.brighten_mask
: A mask indicating the pixels to brighten.bright_factor
: The factor by which to brighten the pixels.
Returns
- The brightened image.
sourceIceFloeTracker.imextendedmin
— Functionimextendedmin(img)
Mimics MATLAB's imextendedmin function that computes the extended-minima transform, which is the regional minima of the H-minima transform. Regional minima are connected components of pixels with a constant intensity value. This function returns a transformed bitmatrix.
Arguments
img
: image objecth
: suppress minima below this depth thresholdconn
: neighborhood connectivity; in 2D 1 = 4-neighborhood and 2 = 8-neighborhood
sourceIceFloeTracker.imgradientmag
— Methodimgradientmag(img)
Compute the gradient magnitude of an image using the Sobel operator.
sourceIceFloeTracker.imhist
— Functionimhist(img, imgtype::AbstractString="uint8")
Compute the histogram of an image where each possible value is represented in the histogram. The function returns a tuple with the bins and counts of each bin.
Example
```jldoctest; setup = :(using IceFloeTracker) julia> img = [ 4 4 4 4 4 3 4 5 4 3 3 5 5 5 3 3 4 5 4 3 4 4 4 4 4 ]
julia> bins, heights = imhist(img);
julia> [bins[heights .> 0] heights[heights .>0]] # display only non-zero bins and heights 3×2 Matrix{Int64}: 3 6 4 14 5 5
sourceIceFloeTracker.impose_minima
— Methodimpose_minima(I::AbstractArray{T}, BW::AbstractArray{Bool}) where {T<:Integer}
Use morphological reconstruction to enforce minima on the input image I
at the positions where the binary mask BW
is non-zero.
It supports both integer and grayscale images using different implementations for each.
sourceIceFloeTracker.imregionalmin
— Functionimregionalmin(img, conn=2)
Compute the regional minima of the image img
using the connectivity conn
.
Returns a bitmatrix of the same size as img
with the regional minima.
Arguments
img
: Image objectconn
: Neighborhood connectivity; in 2D, 1 = 4-neighborhood and 2 = 8-neighborhood
sourceIceFloeTracker.imsharpen
— Functionimsharpen(truecolor_image, landmask_no_dilate, lambda, kappa, niters, nbins, rblocks, cblocks, clip, smoothing_param, intensity)
Sharpen truecolor_image
.
Arguments
truecolor_image
: input image in truecolorlandmask_no_dilate
: landmask for region of interestlambda
: speed of diffusion (0–0.25)kappa
: conduction coefficient for diffusion (25–100)niters
: number of iterations of diffusionnbins
: number of bins during histogram equalizationrblocks
: number of row blocks to divide input image during equalizationcblocks
: number of column blocks to divide input image during equalizationclip
: Thresholds for clipping histogram bins (0–1); values closer to one minimize contrast enhancement, values closer to zero maximize contrast enhancementsmoothing_param
: pixel radius for gaussian blurring (1–10)intensity
: amount of sharpening to perform
sourceIceFloeTracker.imsharpen_gray
— Methodimsharpen_gray(imgsharpened, landmask)
Apply landmask and return Gray type image in colorview for normalization.
sourceIceFloeTracker.isfloegoodmatch
— Methodisfloegoodmatch(conditions, mct, area_mismatch, corr)
Return true
if the floes are a good match as per the set thresholds. Return false
otherwise.
Arguments
conditions
: tuple of booleans for evaluating the conditionsmct
: tuple of thresholds for the match correlation testarea_mismatch
and corr
: values returned by match_corr
sourceIceFloeTracker.isincountourlist
— MethodCheck P
is in countour list if so return the index of the contour that contains P
, otherwise return false.
sourceIceFloeTracker.ismember
— Methodismember(df1,df2)
Return a boolean array indicating whether each row in df1
is a member of df2
.
sourceIceFloeTracker.kmeans_segmentation
— Functionkmeans_segmentation(gray_image, ice_labels;)
Apply k-means segmentation to a gray image to isolate a cluster group representing sea ice. Returns a binary image with ice segmented from background.
Arguments
gray_image
: output image from ice-water-discrimination.jl
or gray ice floe leads image in segmentation_f.jl
ice_labels
: vector if pixel coordinates output from find_ice_labels.jl
sourceIceFloeTracker.loadimg
— Methodloadimg(; dir::String, fname::String)
Load an image from dir
with filename fname
into a matrix of Float64
values. Returns the loaded image.
sourceIceFloeTracker.long_tracker
— Methodlong_tracker(props, condition_thresholds, mc_thresholds)
Track ice floes over multiple days.
Trajectories are built in two steps:
- Get pairs of floes in day 1 and day 2. Any unmatched floes, in both day 1 and day 2, become the "heads" of their respective trajectories.
- For each subsequent day, find pairs of floes for the current trajectory heads. Again, any unmatched floe in the new prop table starts a new trajectory.
Arguments
props::Vector{DataFrame}
: A vector of DataFrames, each containing ice floe properties for a single day. Each DataFrame must have the following columns:- "area"
- "min_row"
- "min_col"
- "max_row"
- "max_col"
- "row_centroid"
- "col_centroid"
- "convex_area"
- "majoraxislength"
- "minoraxislength"
- "orientation"
- "perimeter"
- "mask": 2D array of booleans
- "passtime": A timestamp for the floe
- "psi": the psi-s curve for the floe
- "uuid": a universally unique identifier for each segmented floe
condition_thresholds
: 3-tuple of thresholds (each a named tuple) for deciding whether to match floe i
from day k
to floe j from day k+1
mc_thresholds
: thresholds for area mismatch and psi-s shape correlation
Returns
A DataFrame with the above columns, plus two extra columns, "area_mismatch" and "corr", which are the area mismatch and correlation between a floe and the one that follows it in the trajectory. Trajectories are identified by a unique identifier, "uuid".
sourceIceFloeTracker.make_filename
— Methodmake_filename()
Makes default filename with timestamp.
sourceIceFloeTracker.make_hbreak_dict
— Methodmake_hbreak_dict()
Build dict with the two versions of an H-connected 3x3 neighboorhood.
h1 = [1 0 1 1 1 1 1 0 1]
h2 = [1 1 1 0 1 0 1 1 1]
sourceIceFloeTracker.make_lut
— Methodmake_lut(lutfunc::Function)
Generate lookup table (lut) for 3x3 neighborhoods according to lutfunc
.
sourceIceFloeTracker.make_psi_s
— Methodmake_psi_s(XY::Matrix{<:Number};rangeout::Bool=true,
+unwrap::Bool=true)
Alternate method of make_psi_s
accepting input vectors x
and y
as a 2-column matrix [x y]
in order to facillitate workflow (output from resample_boundary
).
sourceIceFloeTracker.make_psi_s
— Methodmake_psi_s(x::Vector{<:Number},
y::Vector{<:Number};
rangeout::Bool=true,
unwrap::Bool=true)::Tuple{Vector{Float64}, Vector{Float64}}
Builds the ψ-s curve defined by vectors x
and y
.
Returns a tuple of vectors with the phases ψ
in the first component and the traversed arclength in the second component.
Following the convention in [1], the wrapped ψ-s curve has values in [0, 2π) by default; use rangeout
to control this behavior.
See also bwtraceboundary
, resample_boundary
Arguments
x
: Vector of x-coordinatesy
: corresponding vector of y-coordinatesrangeout
: true
(default) for phase values in [0, 2π); false
for phase values in (-π, π].unwrap
: set to true
to get "unwrapped" phases (default).
Reference
[1] McConnell, Ross, et al. "psi-s correlation and dynamic time warping: two methods for tracking ice floes in SAR images." IEEE Transactions on Geoscience and Remote sensing 29.6 (1991): 1004-1012.
Example
The example below builds a cardioid and obtains its ψ-s curve.
julia> t = range(0,2pi,201);
@@ -252,7 +252,7 @@
7.99877 9.35147
7.99926 9.39336
- julia> plot(s, psi) # inspect psi-s curve -- should be a straight line from (0, 0) to (8, 3π)
sourceIceFloeTracker.makecontourstartatP
— MethodMake contour start at point P by permuting the elements in contour.
sourceIceFloeTracker.makeemptydffrom
— Methodmakeemptydffrom(df::DataFrame)
Return an object with an empty dataframe with the same column names as df
and an empty dataframe with column names area
, majoraxis
, minoraxis
, convex_area
, area_mismatch
, and corr
for similarity ratios.
sourceIceFloeTracker.makeemptyratiosdf
— Methodmakeemptyratiosdf()
Return an empty dataframe with column names area
, majoraxis
, minoraxis
, convex_area
, area_mismatch
, and corr
for similarity ratios.
sourceIceFloeTracker.matchcorr
— Methodmatchcorr(
+ julia> plot(s, psi) # inspect psi-s curve -- should be a straight line from (0, 0) to (8, 3π)
sourceIceFloeTracker.makecontourstartatP
— MethodMake contour start at point P by permuting the elements in contour.
sourceIceFloeTracker.makeemptydffrom
— Methodmakeemptydffrom(df::DataFrame)
Return an object with an empty dataframe with the same column names as df
and an empty dataframe with column names area
, majoraxis
, minoraxis
, convex_area
, area_mismatch
, and corr
for similarity ratios.
sourceIceFloeTracker.makeemptyratiosdf
— Methodmakeemptyratiosdf()
Return an empty dataframe with column names area
, majoraxis
, minoraxis
, convex_area
, area_mismatch
, and corr
for similarity ratios.
sourceIceFloeTracker.matchcorr
— Methodmatchcorr(
f1::T,
f2::T,
Δt::F,
@@ -262,18 +262,18 @@
comp::F=0.25,
mm::F=0.22
)
-where {T<:AbstractArray{Bool,2},S<:Int64,F<:Float64}
Compute the mismatch mm
and psi-s-correlation c
for floes with masks f1
and f2
.
The criteria for floes to be considered equivalent is as follows: - c
greater than mm
- _mm
is less than mm
A pair of NaN
is returned for cases for which one of their mask dimension is too small or their sizes are not comparable.
Arguments
f1
: mask of floe 1f2
: mask of floe 2Δt
: time difference between floesmxrot
: maximum rotation (in degrees) allowed between floesn (default: 10)psi
: psi-s-correlation threshold (default: 0.95)sz
: size threshold (default: 16)comp
: size comparability threshold (default: 0.25)mm
: mismatch threshold (default: 0.22)
sourceIceFloeTracker.mean
— Methodmean(x,y)
Compute the mean of x
and y
.
sourceIceFloeTracker.mismatch
— Functionmismatch(fixed::AbstractArray,
+where {T<:AbstractArray{Bool,2},S<:Int64,F<:Float64}
Compute the mismatch mm
and psi-s-correlation c
for floes with masks f1
and f2
.
The criteria for floes to be considered equivalent is as follows: - c
greater than mm
- _mm
is less than mm
A pair of NaN
is returned for cases for which one of their mask dimension is too small or their sizes are not comparable.
Arguments
f1
: mask of floe 1f2
: mask of floe 2Δt
: time difference between floesmxrot
: maximum rotation (in degrees) allowed between floesn (default: 10)psi
: psi-s-correlation threshold (default: 0.95)sz
: size threshold (default: 16)comp
: size comparability threshold (default: 0.25)mm
: mismatch threshold (default: 0.22)
sourceIceFloeTracker.mean
— Methodmean(x,y)
Compute the mean of x
and y
.
sourceIceFloeTracker.mismatch
— Functionmismatch(fixed::AbstractArray,
moving::AbstractArray,
mxshift::Tuple{Int64, Int64}=(10,10),
mxrot::Float64=pi/4;
kwargs...
- )
Estimate a rigid transformation (translation + rotation) that minimizes the 'mismatch' of aligning moving
with fixed
using the QuadDIRECT algorithm.
Returns a pair with the mismatch score mm
and the associated registration angle rot
.
Arguments
fixed
,moving
: images to align via a rigid transformationmxshift
: maximum allowed translation in units of array indices (default set to (10,10)
)mxrot
: maximum allowed rotation in radians (default set to π/4
)thresh
: minimum sum-of-squared-intensity overlap between the images (default is 10% of the sum-of-squared-intensity of fixed
)kwargs
: other arguments such as tol
, ftol
, and fvalue
(see QuadDIRECT.analyze for details)
```
sourceIceFloeTracker.modenan
— Methodmodenan(x::AbstractVector{<:Float64})
Return the mode of x
or NaN
if x
is empty.
sourceIceFloeTracker.morph_fill
— Methodmorph_fill(bw::T)::T where {T<:AbstractArray{Bool}}
Fill holes in binary image bw
by setting 0-valued pixels to 1 if they are surrounded by 1-valued pixels.
Examples
```jldoctest; setup = :(using IceFloeTracker) julia> bw = Bool[ 0 0 0 0 0 0 1 1 1 0 0 1 0 1 0 0 1 1 1 0 0 0 0 0 0 ];
julia> morph_fill(bw) 5×5 Matrix{Bool}: 0 0 0 0 0 0 1 1 1 0 0 1 1 1 0 0 1 1 1 0 0 0 0 0 0
sourceIceFloeTracker.move
— Methodmove from current pixel to the next in given direction
sourceIceFloeTracker.norm
— Methodnorm(v)
Get the euclidean norm of the vector v
.
sourceIceFloeTracker.normalize_image
— Methodnormalize_image(image_sharpened, image_sharpened_gray, landmask, struct_elem;)
Adjusts sharpened land-masked image to highlight ice floe features.
Does reconstruction and landmasking to image_sharpened
.
Arguments
image_sharpened
: sharpened image (output of imsharpen
)image_sharpened_gray
: grayscale, landmasked sharpened image (output of imsharpen_gray(image_sharpened)
)landmask
: landmask for region of intereststruct_elem
: structuring element for dilation
sourceIceFloeTracker.normalizeangle
— Functionnormalizeangle(revised,t=180)
Normalize angle to be between -180 and 180 degrees.
sourceIceFloeTracker.padnhood
— Methodpadnhood(img, I, nhood)
Pad the matrix img[nhood]
with zeros according to the position of I
within the edgesimg
.
Returns img[nhood]
if I
is not an edge index.
sourceIceFloeTracker.pairfloes
— Methodpairfloes(
+ )
Estimate a rigid transformation (translation + rotation) that minimizes the 'mismatch' of aligning moving
with fixed
using the QuadDIRECT algorithm.
Returns a pair with the mismatch score mm
and the associated registration angle rot
.
Arguments
fixed
,moving
: images to align via a rigid transformationmxshift
: maximum allowed translation in units of array indices (default set to (10,10)
)mxrot
: maximum allowed rotation in radians (default set to π/4
)thresh
: minimum sum-of-squared-intensity overlap between the images (default is 10% of the sum-of-squared-intensity of fixed
)kwargs
: other arguments such as tol
, ftol
, and fvalue
(see QuadDIRECT.analyze for details)
```
sourceIceFloeTracker.modenan
— Methodmodenan(x::AbstractVector{<:Float64})
Return the mode of x
or NaN
if x
is empty.
sourceIceFloeTracker.morph_fill
— Methodmorph_fill(bw::T)::T where {T<:AbstractArray{Bool}}
Fill holes in binary image bw
by setting 0-valued pixels to 1 if they are surrounded by 1-valued pixels.
Examples
```jldoctest; setup = :(using IceFloeTracker) julia> bw = Bool[ 0 0 0 0 0 0 1 1 1 0 0 1 0 1 0 0 1 1 1 0 0 0 0 0 0 ];
julia> morph_fill(bw) 5×5 Matrix{Bool}: 0 0 0 0 0 0 1 1 1 0 0 1 1 1 0 0 1 1 1 0 0 0 0 0 0
sourceIceFloeTracker.move
— Methodmove from current pixel to the next in given direction
sourceIceFloeTracker.norm
— Methodnorm(v)
Get the euclidean norm of the vector v
.
sourceIceFloeTracker.normalize_image
— Methodnormalize_image(image_sharpened, image_sharpened_gray, landmask, struct_elem;)
Adjusts sharpened land-masked image to highlight ice floe features.
Does reconstruction and landmasking to image_sharpened
.
Arguments
image_sharpened
: sharpened image (output of imsharpen
)image_sharpened_gray
: grayscale, landmasked sharpened image (output of imsharpen_gray(image_sharpened)
)landmask
: landmask for region of intereststruct_elem
: structuring element for dilation
sourceIceFloeTracker.normalizeangle
— Functionnormalizeangle(revised,t=180)
Normalize angle to be between -180 and 180 degrees.
sourceIceFloeTracker.padnhood
— Methodpadnhood(img, I, nhood)
Pad the matrix img[nhood]
with zeros according to the position of I
within the edgesimg
.
Returns img[nhood]
if I
is not an edge index.
sourceIceFloeTracker.pairfloes
— Methodpairfloes(
segmented_imgs::Vector{BitMatrix},
props::Vector{DataFrame},
passtimes::Vector{DateTime},
latlonrefimage::AbstractString,
condition_thresholds,
-mc_thresholds,
)
Pair floes in props[k]
to floes in props[k+1]
for k=1:length(props)-1
.
The main steps of the algorithm are as follows:
- Crop floes from
segmented_imgs
using bounding box data in props
. Floes in the edges are removed. - For each floekr in
props[k]
, compare to floek+1s in props[k+1]
by computing similarity ratios, set of conditions
, and drift distance dist
. If the conditions are met, compute the area mismatch mm
and psi-s correlation c
for this pair of floes. Pair these two floes if mm
and c
satisfy the thresholds in mc_thresholds
. - If there are collisions (i.e. floe
s
in props[k+1]
is paired to more than one floe in props[k]
), then the floe in props[k]
with the best match is paired to floe s
in props[k+1]
. - Drop paired floes from
props[k]
and props[k+1]
and repeat steps 2 and 3 until there are no more floes to match in props[k]
. - Repeat steps 2-4 for
k=2:length(props)-1
.
Arguments
segmented_imgs
: array of images with segmented floesprops
: array of dataframes containing floe propertiespasstimes
: array of DateTime
objects containing the time of the image in which the floes were capturedlatlonrefimage
: path to geotiff reference image for getting latitude and longitude of floe centroidscondition_thresholds
: 3-tuple of thresholds (each a named tuple) for deciding whether to match floe i
from day k
to floe j from day k+1
mc_thresholds
: thresholds for area mismatch and psi-s shape correlation
Returns a dataframe containing the following columns:
ID
: unique ID for each floe pairing.passtime
: time of the image in which the floes were captured.area
: area of the floe in sq. kilometersconvex_area
: area of the convex hull of the floe in sq. kilometersmajor_axis_length
: length of the major axis of the floe in kilometersminor_axis_length
: length of the minor axis of the floe in kilometersorientation
: angle between the major axis and the x-axis in radiansperimeter
: perimeter of the floe in kilometerslatitude
: latitude of the floe centroidlongitude
: longitude of the floe centroidx
: x-coordinate of the floe centroidy
: y-coordinate of the floe centroidarea_mismatch
: area mismatch between the two floes in rowi and rowi+1 after registrationcorr
: psi-s shape correlation between the two floes in rowi and rowi+1
sourceIceFloeTracker.reconstruct
— Functionreconstruct(img, se, type, invert)
Perform closing/opening by reconstruction on img
.
Arguments
img::AbstractArray
: The input image.se::AbstractArray
: The structuring element.type::String
: The type of morphological operation to perform. Must be either "dilation" (close by reconstruction) or "erosion"
(open by reconstruction).invert::Bool=true
: Invert marker and mask before reconstruction.
sourceIceFloeTracker.regionprops
— Functionregionprops(label_img, ; properties, connectivity)
A wrapper of the regionprops
function from the skimage python library.
See its full documentation at https://scikit-image.org/docs/stable/api/skimage.measure.html#skimage.measure.regionprops.
Arguments
label_img
: Image with the labeled objects of interestintensity_img
: (Optional) Used for generating extra_properties
, integer/float array from which (presumably) label_img
was generatedextra_properties
: (Optional) not yet implemented. It will be set to nothing
See also regionprops_table
Examples
julia> Random.seed!(123);
+mc_thresholds,
)
Pair floes in props[k]
to floes in props[k+1]
for k=1:length(props)-1
.
The main steps of the algorithm are as follows:
- Crop floes from
segmented_imgs
using bounding box data in props
. Floes in the edges are removed. - For each floekr in
props[k]
, compare to floek+1s in props[k+1]
by computing similarity ratios, set of conditions
, and drift distance dist
. If the conditions are met, compute the area mismatch mm
and psi-s correlation c
for this pair of floes. Pair these two floes if mm
and c
satisfy the thresholds in mc_thresholds
. - If there are collisions (i.e. floe
s
in props[k+1]
is paired to more than one floe in props[k]
), then the floe in props[k]
with the best match is paired to floe s
in props[k+1]
. - Drop paired floes from
props[k]
and props[k+1]
and repeat steps 2 and 3 until there are no more floes to match in props[k]
. - Repeat steps 2-4 for
k=2:length(props)-1
.
Arguments
segmented_imgs
: array of images with segmented floesprops
: array of dataframes containing floe propertiespasstimes
: array of DateTime
objects containing the time of the image in which the floes were capturedlatlonrefimage
: path to geotiff reference image for getting latitude and longitude of floe centroidscondition_thresholds
: 3-tuple of thresholds (each a named tuple) for deciding whether to match floe i
from day k
to floe j from day k+1
mc_thresholds
: thresholds for area mismatch and psi-s shape correlation
Returns a dataframe containing the following columns:
ID
: unique ID for each floe pairing.passtime
: time of the image in which the floes were captured.area
: area of the floe in sq. kilometersconvex_area
: area of the convex hull of the floe in sq. kilometersmajor_axis_length
: length of the major axis of the floe in kilometersminor_axis_length
: length of the minor axis of the floe in kilometersorientation
: angle between the major axis and the x-axis in radiansperimeter
: perimeter of the floe in kilometerslatitude
: latitude of the floe centroidlongitude
: longitude of the floe centroidx
: x-coordinate of the floe centroidy
: y-coordinate of the floe centroidarea_mismatch
: area mismatch between the two floes in rowi and rowi+1 after registrationcorr
: psi-s shape correlation between the two floes in rowi and rowi+1
sourceIceFloeTracker.reconstruct
— Functionreconstruct(img, se, type, invert)
Perform closing/opening by reconstruction on img
.
Arguments
img::AbstractArray
: The input image.se::AbstractArray
: The structuring element.type::String
: The type of morphological operation to perform. Must be either "dilation" (close by reconstruction) or "erosion"
(open by reconstruction).invert::Bool=true
: Invert marker and mask before reconstruction.
sourceIceFloeTracker.regionprops
— Functionregionprops(label_img, ; properties, connectivity)
A wrapper of the regionprops
function from the skimage python library.
See its full documentation at https://scikit-image.org/docs/stable/api/skimage.measure.html#skimage.measure.regionprops.
Arguments
label_img
: Image with the labeled objects of interestintensity_img
: (Optional) Used for generating extra_properties
, integer/float array from which (presumably) label_img
was generatedextra_properties
: (Optional) not yet implemented. It will be set to nothing
See also regionprops_table
Examples
julia> Random.seed!(123);
julia> bw_img = rand([0, 1], 5, 10)
5×10 Matrix{Int64}:
@@ -299,7 +299,7 @@
13 11.621320343559642
1 0.0
1 0.0
-7 4.621320343559642
sourceIceFloeTracker.regionprops_table
— Functionregionprops_table(label_img, intensity_img; properties, connectivity, extra_properties)
A wrapper of the regionprops_table
function from the skimage python library.
See its full documentation at https://scikit-image.org/docs/stable/api/skimage.measure.html#regionprops-table.
Arguments
label_img
: Image with the labeled objects of interestintensity_img
: (Optional) Used for generating extra_properties
, integer/float array from which (presumably) label_img
was generated properties
: List (Vector
or Tuple
) of properties to be generated for each connected component in label_img
extra_properties
: (Optional) not yet implemented. It will be set to nothing
Notes
- Zero indexing has been corrected for the
bbox
and centroid
properties bbox
data (max_col
and max_row
) are inclusivecentroid
data are rounded to the nearest integer
See also regionprops
Examples
julia> using IceFloeTracker, Random
+7 4.621320343559642
sourceIceFloeTracker.regionprops_table
— Functionregionprops_table(label_img, intensity_img; properties, connectivity, extra_properties)
A wrapper of the regionprops_table
function from the skimage python library.
See its full documentation at https://scikit-image.org/docs/stable/api/skimage.measure.html#regionprops-table.
Arguments
label_img
: Image with the labeled objects of interestintensity_img
: (Optional) Used for generating extra_properties
, integer/float array from which (presumably) label_img
was generated properties
: List (Vector
or Tuple
) of properties to be generated for each connected component in label_img
extra_properties
: (Optional) not yet implemented. It will be set to nothing
Notes
- Zero indexing has been corrected for the
bbox
and centroid
properties bbox
data (max_col
and max_row
) are inclusivecentroid
data are rounded to the nearest integer
See also regionprops
Examples
julia> using IceFloeTracker, Random
julia> Random.seed!(123);
@@ -332,17 +332,17 @@
1 │ 13 11.6213
2 │ 1 0.0
3 │ 1 0.0
- 4 │ 7 4.62132
sourceIceFloeTracker.regularize_fill_holes
— Methodregularize_fill_holes(img, local_maxima_mask, factor, segment_mask, L0mask)
Regularize img
by: 1. increasing the maxima of img
by a factor of factor
2. filtering img
at positions where either segment_mask
or L0mask
are true 3. filling holes
Arguments
img
: The morphological residue image.local_maxima_mask
: The local maxima mask.factor
: The factor to apply to the local maxima mask.segment_mask
: The segment mask – intersection of bw1 and bw2 in first tiled workflow of master.m
.L0mask
: zero-labeled pixels from watershed.
sourceIceFloeTracker.regularize_sharpening
— Methodregularize_sharpening(img, L0mask, radius, amount, local_maxima_mask, factor, segment_mask, se)
Regularize img
via sharpening, filtering, reconstruction, and maxima elevating.
Arguments
img
: The input image.L0mask
: zero-labeled pixels from watershed.radius
: The radius of the unsharp mask.amount
: The amount of unsharp mask.local_maxima_mask
: The local maxima mask.factor
: The factor to apply to the local maxima mask.segment_mask
: The segment mask – intersection of bw1 and bw2 in first tiled workflow of master.m
.
sourceIceFloeTracker.remove_padding
— Methodremove_padding(paddedimg, border_spec)
Removes padding from the boundary of padded image paddedimg
according to the border specification border_spec
type. Returns the cropped image.
Arguments
paddedimg
: Pre-padded image.border_spec
: Type representing the style of padding (such as Pad
or Fill
) with which paddedimg
is assumend to be pre-padded. Example: Pad((1,2), (3,4))
specifies 1 row on the top, 2 columns on the left, 3 rows on the bottom, and 4 columns on the right boundary.
See also add_padding
sourceIceFloeTracker.renamecols!
— Methodrenamecols!(props::DataFrame, oldnames::Vector{T}, newnames::Vector{T}) where T<:Union{Symbol,String}
Rename the oldnames
columns in props
to newnames
.
sourceIceFloeTracker.resample_boundary
— Functionresample_boundary(bd_points::Vector{<:CartesianIndex}, reduc_factor::Int64=2, bd::String="natural")
Get a uniform set of resampled boundary points from bd_points
using cubic splines with specified boundary conditions
The resampled set of points is obtained using parametric interpolation of the points in bd_points
. It is assumed that the separation between a pair of adjacent points is 1.
Arguments
bd_points
: Sequetial set of boundary points for the object of interestreduc_factor
: Factor by which to reduce the number of points in bd_points
(2 by default)
-bd
: Boundary condition, either 'natural' (default) or 'periodic'
See also bwtraceboundary
Example
```jldoctest; setup = :(using IceFloeTracker) julia> A = zeros(Int, 13, 16); A[2:6, 2:6] .= 1; A[4:8, 7:10] .= 1; A[10:12,13:15] .= 1; A[10:12,3:6] .= 1; julia> A 13×16 Matrix{Int64}: 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 0 0 0 0 0 0 1 1 1 0 0 0 1 1 1 1 0 0 0 0 0 0 1 1 1 0 0 0 1 1 1 1 0 0 0 0 0 0 1 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
julia> boundary = bwtraceboundary(A);
julia> boundary[3] 9-element Vector{CartesianIndex}: CartesianIndex(10, 13) CartesianIndex(11, 13) CartesianIndex(12, 13) CartesianIndex(12, 14) CartesianIndex(12, 15) CartesianIndex(11, 15) CartesianIndex(10, 15) CartesianIndex(10, 14) CartesianIndex(10, 13)
julia> resample_boundary(boundary[3]) 4×2 Matrix{Float64}: 10.0 13.0 12.0357 13.5859 10.5859 15.0357 10.0 13.0
sourceIceFloeTracker.reset_id!
— Functionreset_id!(df, col)
Reset the distinct values in the column col
of df
to be consecutive integers starting from 1.
sourceIceFloeTracker.rgb2gray
— Methodrgb2gray(rgbchannels::Array{Float64, 3})
Convert an array of RGB channel data to grayscale in the range [0, 255].
Identical to MATLAB rgb2gray
(https://www.mathworks.com/help/matlab/ref/rgb2gray.html).
sourceIceFloeTracker.rgb2gray
— Methodrgb2gray(img::Matrix{RGB{Float64}})
Convert an RGB image to grayscale in the range [0, 255].
sourceIceFloeTracker.roundtoint!
— Methodroundtoint!(props::DataFrame, colnames::Vector{T}) where T<:Union{Symbol,String}
Round the colnames
columns in props
to Int
.
sourceIceFloeTracker.segmentation_A
— Methodsegmentation_A(segmented_ice_cloudmasked; min_opening_area)
Apply k-means segmentation to a gray image to isolate a cluster group representing sea ice. Returns an image segmented and processed as well as an intermediate files needed for downstream functions.
Arguments
segmented_ice_cloudmask
: bitmatrix with open water/clouds = 0, ice = 1, output from segmented_ice_cloudmasking()
min_opening_area
: minimum size of pixels to use during morphological openingfill_range
: range of values dictating the size of holes to fill
sourceIceFloeTracker.segmentation_B
— Functionsegmentation_B(sharpened_image, cloudmask, segmented_a_ice_mask, struct_elem; fill_range, isolation_threshold, alpha_level, adjusted_ice_threshold)
Performs image processing and morphological filtering with intermediate files from normalization.jl and segmentation_A to further isolate ice floes, returning a mask of potential ice.
Arguments
sharpened_image
: non-cloudmasked but sharpened image, output from normalization.jl
cloudmask
: bitmatrix cloudmask for region of interestsegmented_a_ice_mask
: binary cloudmasked ice mask from segmentation_a_direct.jl
struct_elem
: structuring element for dilationfill_range
: range of values dictating the size of holes to fillisolation_threshold
: threshold used to isolated pixels from sharpened_image
; between 0-1alpha_level
: alpha threshold used to adjust contrastgamma_factor
: amount of gamma adjustmentadjusted_ice_threshold
: threshold used to set ice equal to one after gamma adjustment
sourceIceFloeTracker.segmentation_F
— Methodsegmentation_F(
+ 4 │ 7 4.62132
sourceIceFloeTracker.regularize_fill_holes
— Methodregularize_fill_holes(img, local_maxima_mask, factor, segment_mask, L0mask)
Regularize img
by: 1. increasing the maxima of img
by a factor of factor
2. filtering img
at positions where either segment_mask
or L0mask
are true 3. filling holes
Arguments
img
: The morphological residue image.local_maxima_mask
: The local maxima mask.factor
: The factor to apply to the local maxima mask.segment_mask
: The segment mask – intersection of bw1 and bw2 in first tiled workflow of master.m
.L0mask
: zero-labeled pixels from watershed.
sourceIceFloeTracker.regularize_sharpening
— Methodregularize_sharpening(img, L0mask, radius, amount, local_maxima_mask, factor, segment_mask, se)
Regularize img
via sharpening, filtering, reconstruction, and maxima elevating.
Arguments
img
: The input image.L0mask
: zero-labeled pixels from watershed.radius
: The radius of the unsharp mask.amount
: The amount of unsharp mask.local_maxima_mask
: The local maxima mask.factor
: The factor to apply to the local maxima mask.segment_mask
: The segment mask – intersection of bw1 and bw2 in first tiled workflow of master.m
.
sourceIceFloeTracker.remove_padding
— Methodremove_padding(paddedimg, border_spec)
Removes padding from the boundary of padded image paddedimg
according to the border specification border_spec
type. Returns the cropped image.
Arguments
paddedimg
: Pre-padded image.border_spec
: Type representing the style of padding (such as Pad
or Fill
) with which paddedimg
is assumend to be pre-padded. Example: Pad((1,2), (3,4))
specifies 1 row on the top, 2 columns on the left, 3 rows on the bottom, and 4 columns on the right boundary.
See also add_padding
sourceIceFloeTracker.renamecols!
— Methodrenamecols!(props::DataFrame, oldnames::Vector{T}, newnames::Vector{T}) where T<:Union{Symbol,String}
Rename the oldnames
columns in props
to newnames
.
sourceIceFloeTracker.resample_boundary
— Functionresample_boundary(bd_points::Vector{<:CartesianIndex}, reduc_factor::Int64=2, bd::String="natural")
Get a uniform set of resampled boundary points from bd_points
using cubic splines with specified boundary conditions
The resampled set of points is obtained using parametric interpolation of the points in bd_points
. It is assumed that the separation between a pair of adjacent points is 1.
Arguments
bd_points
: Sequetial set of boundary points for the object of interestreduc_factor
: Factor by which to reduce the number of points in bd_points
(2 by default)
-bd
: Boundary condition, either 'natural' (default) or 'periodic'
See also bwtraceboundary
Example
```jldoctest; setup = :(using IceFloeTracker) julia> A = zeros(Int, 13, 16); A[2:6, 2:6] .= 1; A[4:8, 7:10] .= 1; A[10:12,13:15] .= 1; A[10:12,3:6] .= 1; julia> A 13×16 Matrix{Int64}: 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 0 0 0 0 0 0 1 1 1 0 0 0 1 1 1 1 0 0 0 0 0 0 1 1 1 0 0 0 1 1 1 1 0 0 0 0 0 0 1 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
julia> boundary = bwtraceboundary(A);
julia> boundary[3] 9-element Vector{CartesianIndex}: CartesianIndex(10, 13) CartesianIndex(11, 13) CartesianIndex(12, 13) CartesianIndex(12, 14) CartesianIndex(12, 15) CartesianIndex(11, 15) CartesianIndex(10, 15) CartesianIndex(10, 14) CartesianIndex(10, 13)
julia> resample_boundary(boundary[3]) 4×2 Matrix{Float64}: 10.0 13.0 12.0357 13.5859 10.5859 15.0357 10.0 13.0
sourceIceFloeTracker.reset_id!
— Functionreset_id!(df, col)
Reset the distinct values in the column col
of df
to be consecutive integers starting from 1.
sourceIceFloeTracker.rgb2gray
— Methodrgb2gray(rgbchannels::Array{Float64, 3})
Convert an array of RGB channel data to grayscale in the range [0, 255].
Identical to MATLAB rgb2gray
(https://www.mathworks.com/help/matlab/ref/rgb2gray.html).
sourceIceFloeTracker.rgb2gray
— Methodrgb2gray(img::Matrix{RGB{Float64}})
Convert an RGB image to grayscale in the range [0, 255].
sourceIceFloeTracker.roundtoint!
— Methodroundtoint!(props::DataFrame, colnames::Vector{T}) where T<:Union{Symbol,String}
Round the colnames
columns in props
to Int
.
sourceIceFloeTracker.segmentation_A
— Methodsegmentation_A(segmented_ice_cloudmasked; min_opening_area)
Apply k-means segmentation to a gray image to isolate a cluster group representing sea ice. Returns an image segmented and processed as well as an intermediate files needed for downstream functions.
Arguments
segmented_ice_cloudmask
: bitmatrix with open water/clouds = 0, ice = 1, output from segmented_ice_cloudmasking()
min_opening_area
: minimum size of pixels to use during morphological openingfill_range
: range of values dictating the size of holes to fill
sourceIceFloeTracker.segmentation_B
— Functionsegmentation_B(sharpened_image, cloudmask, segmented_a_ice_mask, struct_elem; fill_range, isolation_threshold, alpha_level, adjusted_ice_threshold)
Performs image processing and morphological filtering with intermediate files from normalization.jl and segmentation_A to further isolate ice floes, returning a mask of potential ice.
Arguments
sharpened_image
: non-cloudmasked but sharpened image, output from normalization.jl
cloudmask
: bitmatrix cloudmask for region of interestsegmented_a_ice_mask
: binary cloudmasked ice mask from segmentation_a_direct.jl
struct_elem
: structuring element for dilationfill_range
: range of values dictating the size of holes to fillisolation_threshold
: threshold used to isolated pixels from sharpened_image
; between 0-1alpha_level
: alpha threshold used to adjust contrastgamma_factor
: amount of gamma adjustmentadjusted_ice_threshold
: threshold used to set ice equal to one after gamma adjustment
sourceIceFloeTracker.segmentation_F
— Methodsegmentation_F(
segmentation_B_not_ice_mask::Matrix{Gray{Float64}},
segmentation_B_ice_intersect::BitMatrix,
segmentation_B_watershed_intersect::BitMatrix,
ice_labels::Vector{Int64},
cloudmask::BitMatrix,
landmask::BitMatrix;
-min_area_opening::Int64=20
)
Cleans up past segmentation images with morphological operations, and applies the results of prior watershed segmentation, returning the final cleaned image for tracking with ice floes segmented and isolated.
Arguments
segmentation_B_not_ice_mask
: gray image output from segmentation_b.jl
segmentation_B_ice_intersect
: binary mask output from segmentation_b.jl
segmentation_B_watershed_intersect
: ice pixels, output from segmentation_b.jl
ice_labels
: vector of pixel coordinates output from find_ice_labels.jl
cloudmask.jl
: bitmatrix cloudmask for region of interestlandmask.jl
: bitmatrix landmask for region of interestmin_area_opening
: threshold used for area opening; pixel groups greater than threshold are retained
sourceIceFloeTracker.segmented_ice_cloudmasking
— Methodsegmented_ice_cloudmasking(gray_image, cloudmask, ice_labels;)
Apply cloudmask to a bitmatrix of segmented ice after kmeans clustering. Returns a bitmatrix with open water/clouds = 0, ice = 1).
Arguments
gray_image
: output image from ice-water-discrimination.jl
or gray ice floe leads image in segmentation_f.jl
cloudmask
: bitmatrix cloudmask for region of interestice_labels
: vector if pixel coordinates output from find_ice_labels.jl
sourceIceFloeTracker.sort_floes_by_area!
— Methodsort_floes_by_area!(props)
Sort floes in props
by area in descending order.
sourceIceFloeTracker.timestamp
— Methodtimestamp(fname)
Attach timestamp to fname
.
sourceIceFloeTracker.trackercond1
— Functiontrackercond1(p1, p2, delta_time, t1=(dt = (30, 100, 1300), dist=(15, 30, 120)))
Return true
if the floe at p1
and the floe at p2
are within a certain distance of each other and the displacement time is within a certain range. Return false
otherwise.
Arguments
p1
: coordinates of floe 1p2
: coordinates of floe 2delta_time
: time elapsed from image day 1 to image day 2t
: tuple of thresholds for elapsed time and distance
sourceIceFloeTracker.trackercond2
— Functiontrackercond2(area1, ratios, t2=(area=1200, arearatio=0.28, majaxisratio=0.10, minaxisratio=0.12, convex_area=0.14))
Set of conditions for "big" floes. Return true
if the area of the floe is greater than t2.area
and the similarity ratios are less than the corresponding thresholds in t2
. Return false
otherwise.
sourceIceFloeTracker.trackercond3
— Functiontrackercond3(area1, ratios, t3=(area=1200, arearatio=0.18, majaxisratio=0.07, minaxisratio=0.08, convex_area=0.09))
Set of conditions for "small" floes. Return true
if the area of the floe is less than t3.area
and the similarity ratios are less than the corresponding thresholds in t3
. Return false
otherwise
sourceIceFloeTracker.unsharp_mask
— Methodunsharp_mask(image_gray, smoothing_param, intensity, clampmax)
Apply unsharp masking on (equalized) grayscale ([0, clampmax
]) image to enhance its sharpness.
Arguments
image_gray
: The input grayscale image, typically already equalized.smoothing_param::Int
: The pixel radius for Gaussian blurring (typically between 1 and 10).intensity
: The amount of sharpening to apply. Higher values result in more pronounced sharpening.clampmax
: upper limit of intensity values in the returned image.`
Returns
The sharpened grayscale image with values clipped between 0 and clapmax
.
sourceIceFloeTracker.unsharp_mask
— Method(Deprecated)
-unsharp_mask(image_gray, smoothing_param, intensity)
Apply unsharp masking on (equalized) grayscale image to enhance its sharpness.
Does not perform clamping after the smoothing step. Kept for legacy tests of IceFloeTracker.jl.
Arguments
image_gray
: The input grayscale image, typically already equalized.smoothing_param::Int
: The pixel radius for Gaussian blurring (typically between 1 and 10).intensity
: The amount of sharpening to apply. Higher values result in more pronounced sharpening.
Returns
The sharpened grayscale image with values clipped between 0 and clapmax
.
sourceIceFloeTracker.update!
— Methodupdate!(match_total::MatchedPairs, matched_pairs::MatchedPairs)
Update match_total
with the data from matched_pairs
.
sourceIceFloeTracker.watershed_ice_floes
— Methodwatershed_ice_floes(intermediate_segmentation_image;)
Performs image processing and watershed segmentation with intermediate files from segmentation_b.jl to further isolate ice floes, returning a binary segmentation mask indicating potential sparse boundaries of ice floes.
Arguments
-intermediate_segmentation_image
: binary cloudmasked and landmasked intermediate file from segmentation B, either SegB.not_ice_bit
or SegB.ice_intersect
sourceIceFloeTracker.watershed_product
— Methodwatershed_product(watershed_B_ice_intersect, watershed_B_not_ice;)
Intersects the outputs of watershed segmentation on intermediate files from segmentation B, indicating potential sparse boundaries of ice floes.
Arguments
watershed_B_ice_intersect
: binary segmentation mask from watershed_ice_floes
watershed_B_not_ice
: binary segmentation mask from watershed_ice_floes
sourceIceFloeTracker.@persist
— Macro@persist img fname
+min_area_opening::Int64=20
)
Cleans up past segmentation images with morphological operations, and applies the results of prior watershed segmentation, returning the final cleaned image for tracking with ice floes segmented and isolated.
Arguments
segmentation_B_not_ice_mask
: gray image output from segmentation_b.jl
segmentation_B_ice_intersect
: binary mask output from segmentation_b.jl
segmentation_B_watershed_intersect
: ice pixels, output from segmentation_b.jl
ice_labels
: vector of pixel coordinates output from find_ice_labels.jl
cloudmask.jl
: bitmatrix cloudmask for region of interestlandmask.jl
: bitmatrix landmask for region of interestmin_area_opening
: threshold used for area opening; pixel groups greater than threshold are retained
sourceIceFloeTracker.segmented_ice_cloudmasking
— Methodsegmented_ice_cloudmasking(gray_image, cloudmask, ice_labels;)
Apply cloudmask to a bitmatrix of segmented ice after kmeans clustering. Returns a bitmatrix with open water/clouds = 0, ice = 1).
Arguments
gray_image
: output image from ice-water-discrimination.jl
or gray ice floe leads image in segmentation_f.jl
cloudmask
: bitmatrix cloudmask for region of interestice_labels
: vector if pixel coordinates output from find_ice_labels.jl
sourceIceFloeTracker.sort_floes_by_area!
— Methodsort_floes_by_area!(props)
Sort floes in props
by area in descending order.
sourceIceFloeTracker.timestamp
— Methodtimestamp(fname)
Attach timestamp to fname
.
sourceIceFloeTracker.trackercond1
— Functiontrackercond1(p1, p2, delta_time, t1=(dt = (30, 100, 1300), dist=(15, 30, 120)))
Return true
if the floe at p1
and the floe at p2
are within a certain distance of each other and the displacement time is within a certain range. Return false
otherwise.
Arguments
p1
: coordinates of floe 1p2
: coordinates of floe 2delta_time
: time elapsed from image day 1 to image day 2t
: tuple of thresholds for elapsed time and distance
sourceIceFloeTracker.trackercond2
— Functiontrackercond2(area1, ratios, t2=(area=1200, arearatio=0.28, majaxisratio=0.10, minaxisratio=0.12, convex_area=0.14))
Set of conditions for "big" floes. Return true
if the area of the floe is greater than t2.area
and the similarity ratios are less than the corresponding thresholds in t2
. Return false
otherwise.
sourceIceFloeTracker.trackercond3
— Functiontrackercond3(area1, ratios, t3=(area=1200, arearatio=0.18, majaxisratio=0.07, minaxisratio=0.08, convex_area=0.09))
Set of conditions for "small" floes. Return true
if the area of the floe is less than t3.area
and the similarity ratios are less than the corresponding thresholds in t3
. Return false
otherwise
sourceIceFloeTracker.unsharp_mask
— Methodunsharp_mask(image_gray, smoothing_param, intensity, clampmax)
Apply unsharp masking on (equalized) grayscale ([0, clampmax
]) image to enhance its sharpness.
Arguments
image_gray
: The input grayscale image, typically already equalized.smoothing_param::Int
: The pixel radius for Gaussian blurring (typically between 1 and 10).intensity
: The amount of sharpening to apply. Higher values result in more pronounced sharpening.clampmax
: upper limit of intensity values in the returned image.`
Returns
The sharpened grayscale image with values clipped between 0 and clapmax
.
sourceIceFloeTracker.unsharp_mask
— Method(Deprecated)
+unsharp_mask(image_gray, smoothing_param, intensity)
Apply unsharp masking on (equalized) grayscale image to enhance its sharpness.
Does not perform clamping after the smoothing step. Kept for legacy tests of IceFloeTracker.jl.
Arguments
image_gray
: The input grayscale image, typically already equalized.smoothing_param::Int
: The pixel radius for Gaussian blurring (typically between 1 and 10).intensity
: The amount of sharpening to apply. Higher values result in more pronounced sharpening.
Returns
The sharpened grayscale image with values clipped between 0 and clapmax
.
sourceIceFloeTracker.update!
— Methodupdate!(match_total::MatchedPairs, matched_pairs::MatchedPairs)
Update match_total
with the data from matched_pairs
.
sourceIceFloeTracker.watershed_ice_floes
— Methodwatershed_ice_floes(intermediate_segmentation_image;)
Performs image processing and watershed segmentation with intermediate files from segmentation_b.jl to further isolate ice floes, returning a binary segmentation mask indicating potential sparse boundaries of ice floes.
Arguments
-intermediate_segmentation_image
: binary cloudmasked and landmasked intermediate file from segmentation B, either SegB.not_ice_bit
or SegB.ice_intersect
sourceIceFloeTracker.watershed_product
— Methodwatershed_product(watershed_B_ice_intersect, watershed_B_not_ice;)
Intersects the outputs of watershed segmentation on intermediate files from segmentation B, indicating potential sparse boundaries of ice floes.
Arguments
watershed_B_ice_intersect
: binary segmentation mask from watershed_ice_floes
watershed_B_not_ice
: binary segmentation mask from watershed_ice_floes
sourceIceFloeTracker.@persist
— Macro@persist img fname
@persist(img,fname)
@persist img
@persist(img)
@persist img fname ts
-@persist(img, fname, ts)
Given a reference to an image object img
, the macro persists (saves to a file) img
to the current working directory using fname
as filename. Returns img
.
Arguments
img
: Symbol expression representing an image object loaded in memory.fname
: Optional filename for the persisted image.ts
: Optional boolean to attach timestamp to fname
.
sourceIceFloeTracker.MatchedPairs
— TypeContainer for matched pairs of floes. props1
and props2
are dataframes with the same column names as the input dataframes. ratios
is a dataframe with column names area
, majoraxis
, minoraxis
, convex_area
, area_mismatch
, and corr
for similarity ratios. dist
is a vector of (pixel) distances between paired floes.
sourceIceFloeTracker.MatchedPairs
— MethodMatchedPairs(df)
Return an object of type MatchedPairs
with an empty dataframe with the same column names as df
, an empty dataframe with column names area
, majoraxis
, minoraxis
, convex_area
, area_mismatch
, and corr
for similarity ratios, and an empty vector for distances.
sourceIceFloeTracker.Tracked
— TypeContainer for final matched pairs of floes. data
is a vector of MatchedPairs
objects.
sourceIndex
IceFloeTracker.MatchedPairs
IceFloeTracker.MatchedPairs
IceFloeTracker.Tracked
Base.isequal
Base.sort!
IceFloeTracker._adjust_histogram
IceFloeTracker._bin9todec
IceFloeTracker._branch_filter
IceFloeTracker._generate_se!
IceFloeTracker._operator_lut
IceFloeTracker._swap_last_values!
IceFloeTracker.absdiffmeanratio
IceFloeTracker.add_padding
IceFloeTracker.add_passtimes!
IceFloeTracker.addfloemasks!
IceFloeTracker.addlatlon!
IceFloeTracker.addmatch!
IceFloeTracker.adduuid!
IceFloeTracker.appendrows!
IceFloeTracker.apply_cloudmask
IceFloeTracker.apply_landmask
IceFloeTracker.atan2
IceFloeTracker.binarize_landmask
IceFloeTracker.branch
IceFloeTracker.branch_candidates_func
IceFloeTracker.bridge
IceFloeTracker.bump_tile
IceFloeTracker.bwareamaxfilt
IceFloeTracker.bwareamaxfilt!
IceFloeTracker.bwdist
IceFloeTracker.bwperim
IceFloeTracker.bwtraceboundary
IceFloeTracker.callmatchcorr
IceFloeTracker.check_fname
IceFloeTracker.clockwise
IceFloeTracker.compute_ratios
IceFloeTracker.compute_ratios_conditions
IceFloeTracker.conditional_histeq
IceFloeTracker.conditional_histeq
IceFloeTracker.connected_background_count
IceFloeTracker.consolidate_matched_pairs
IceFloeTracker.convertcentroid!
IceFloeTracker.converttounits!
IceFloeTracker.corr
IceFloeTracker.corr
IceFloeTracker.counterclockwise
IceFloeTracker.create_cloudmask
IceFloeTracker.create_landmask
IceFloeTracker.cropfloe
IceFloeTracker.cropfloe
IceFloeTracker.cropfloe
IceFloeTracker.crosscorr
IceFloeTracker.deleteallbut!
IceFloeTracker.detect_move!
IceFloeTracker.discriminate_ice_water
IceFloeTracker.dist
IceFloeTracker.filt_except_label
IceFloeTracker.filt_except_label!
IceFloeTracker.find_floe_matches
IceFloeTracker.find_ice_labels
IceFloeTracker.find_reflectance_peaks
IceFloeTracker.fixzeroindexing!
IceFloeTracker.fname_ext_splice
IceFloeTracker.fname_ext_split
IceFloeTracker.from_to
IceFloeTracker.get_area_missed
IceFloeTracker.get_areas
IceFloeTracker.get_brighten_mask
IceFloeTracker.get_dt
IceFloeTracker.get_final
IceFloeTracker.get_ice_masks
IceFloeTracker.get_matches
IceFloeTracker.get_max_label
IceFloeTracker.get_optimal_tile_size
IceFloeTracker.get_rgb_channels
IceFloeTracker.get_tile_dims
IceFloeTracker.get_tile_meta
IceFloeTracker.get_tiles
IceFloeTracker.get_tiles
IceFloeTracker.get_trajectory_heads
IceFloeTracker.get_unmatched
IceFloeTracker.getbboxcolumns
IceFloeTracker.getbestmatchdata
IceFloeTracker.getcentroid
IceFloeTracker.getcentroidcolumns
IceFloeTracker.getcollisions
IceFloeTracker.getcollisionslocs
IceFloeTracker.getfit
IceFloeTracker.getfloemasks
IceFloeTracker.getidxmostminimumeverything
IceFloeTracker.getidxofrow
IceFloeTracker.getpropsday1day2
IceFloeTracker.getsizecomparability
IceFloeTracker.grad
IceFloeTracker.grad
IceFloeTracker.hbreak
IceFloeTracker.hbreak!
IceFloeTracker.histeq
IceFloeTracker.imadjust
IceFloeTracker.imbrighten
IceFloeTracker.imextendedmin
IceFloeTracker.imgradientmag
IceFloeTracker.imhist
IceFloeTracker.impose_minima
IceFloeTracker.imregionalmin
IceFloeTracker.imsharpen
IceFloeTracker.imsharpen_gray
IceFloeTracker.isfloegoodmatch
IceFloeTracker.isincountourlist
IceFloeTracker.ismember
IceFloeTracker.kmeans_segmentation
IceFloeTracker.loadimg
IceFloeTracker.long_tracker
IceFloeTracker.make_filename
IceFloeTracker.make_hbreak_dict
IceFloeTracker.make_lut
IceFloeTracker.make_psi_s
IceFloeTracker.make_psi_s
IceFloeTracker.makecontourstartatP
IceFloeTracker.makeemptydffrom
IceFloeTracker.makeemptyratiosdf
IceFloeTracker.matchcorr
IceFloeTracker.mean
IceFloeTracker.mismatch
IceFloeTracker.modenan
IceFloeTracker.morph_fill
IceFloeTracker.move
IceFloeTracker.norm
IceFloeTracker.normalize_image
IceFloeTracker.normalizeangle
IceFloeTracker.padnhood
IceFloeTracker.pairfloes
IceFloeTracker.reconstruct
IceFloeTracker.regionprops
IceFloeTracker.regionprops_table
IceFloeTracker.regularize_fill_holes
IceFloeTracker.regularize_sharpening
IceFloeTracker.remove_padding
IceFloeTracker.renamecols!
IceFloeTracker.resample_boundary
IceFloeTracker.reset_id!
IceFloeTracker.rgb2gray
IceFloeTracker.rgb2gray
IceFloeTracker.roundtoint!
IceFloeTracker.segmentation_A
IceFloeTracker.segmentation_B
IceFloeTracker.segmentation_F
IceFloeTracker.segmented_ice_cloudmasking
IceFloeTracker.sort_floes_by_area!
IceFloeTracker.timestamp
IceFloeTracker.trackercond1
IceFloeTracker.trackercond2
IceFloeTracker.trackercond3
IceFloeTracker.unsharp_mask
IceFloeTracker.unsharp_mask
IceFloeTracker.update!
IceFloeTracker.watershed_ice_floes
IceFloeTracker.watershed_product
IceFloeTracker.@persist
Settings
This document was generated with Documenter.jl version 0.27.21 on Thursday 19 December 2024. Using Julia version 1.7.3.
+@persist(img, fname, ts)
Given a reference to an image object img
, the macro persists (saves to a file) img
to the current working directory using fname
as filename. Returns img
.
Arguments
img
: Symbol expression representing an image object loaded in memory.fname
: Optional filename for the persisted image.ts
: Optional boolean to attach timestamp tofname
.
IceFloeTracker.MatchedPairs
— TypeContainer for matched pairs of floes. props1
and props2
are dataframes with the same column names as the input dataframes. ratios
is a dataframe with column names area
, majoraxis
, minoraxis
, convex_area
, area_mismatch
, and corr
for similarity ratios. dist
is a vector of (pixel) distances between paired floes.
IceFloeTracker.MatchedPairs
— MethodMatchedPairs(df)
Return an object of type MatchedPairs
with an empty dataframe with the same column names as df
, an empty dataframe with column names area
, majoraxis
, minoraxis
, convex_area
, area_mismatch
, and corr
for similarity ratios, and an empty vector for distances.
IceFloeTracker.Tracked
— TypeContainer for final matched pairs of floes. data
is a vector of MatchedPairs
objects.
Index
IceFloeTracker.MatchedPairs
IceFloeTracker.MatchedPairs
IceFloeTracker.Tracked
Base.isequal
Base.sort!
IceFloeTracker._adjust_histogram
IceFloeTracker._bin9todec
IceFloeTracker._branch_filter
IceFloeTracker._generate_se!
IceFloeTracker._operator_lut
IceFloeTracker._swap_last_values!
IceFloeTracker.absdiffmeanratio
IceFloeTracker.add_padding
IceFloeTracker.add_passtimes!
IceFloeTracker.addfloemasks!
IceFloeTracker.addlatlon!
IceFloeTracker.addmatch!
IceFloeTracker.adduuid!
IceFloeTracker.adduuid!
IceFloeTracker.appendrows!
IceFloeTracker.apply_cloudmask
IceFloeTracker.apply_landmask
IceFloeTracker.atan2
IceFloeTracker.binarize_landmask
IceFloeTracker.branch
IceFloeTracker.branch_candidates_func
IceFloeTracker.bridge
IceFloeTracker.bump_tile
IceFloeTracker.bwareamaxfilt
IceFloeTracker.bwareamaxfilt!
IceFloeTracker.bwdist
IceFloeTracker.bwperim
IceFloeTracker.bwtraceboundary
IceFloeTracker.callmatchcorr
IceFloeTracker.check_fname
IceFloeTracker.clockwise
IceFloeTracker.compute_ratios
IceFloeTracker.compute_ratios_conditions
IceFloeTracker.conditional_histeq
IceFloeTracker.conditional_histeq
IceFloeTracker.connected_background_count
IceFloeTracker.consolidate_matched_pairs
IceFloeTracker.convertcentroid!
IceFloeTracker.converttounits!
IceFloeTracker.corr
IceFloeTracker.corr
IceFloeTracker.counterclockwise
IceFloeTracker.create_cloudmask
IceFloeTracker.create_landmask
IceFloeTracker.cropfloe
IceFloeTracker.cropfloe
IceFloeTracker.cropfloe
IceFloeTracker.crosscorr
IceFloeTracker.deleteallbut!
IceFloeTracker.detect_move!
IceFloeTracker.discriminate_ice_water
IceFloeTracker.dist
IceFloeTracker.filt_except_label
IceFloeTracker.filt_except_label!
IceFloeTracker.find_floe_matches
IceFloeTracker.find_ice_labels
IceFloeTracker.find_reflectance_peaks
IceFloeTracker.fixzeroindexing!
IceFloeTracker.fname_ext_splice
IceFloeTracker.fname_ext_split
IceFloeTracker.from_to
IceFloeTracker.get_area_missed
IceFloeTracker.get_areas
IceFloeTracker.get_brighten_mask
IceFloeTracker.get_dt
IceFloeTracker.get_final
IceFloeTracker.get_ice_masks
IceFloeTracker.get_matches
IceFloeTracker.get_max_label
IceFloeTracker.get_optimal_tile_size
IceFloeTracker.get_rgb_channels
IceFloeTracker.get_tile_dims
IceFloeTracker.get_tile_meta
IceFloeTracker.get_tiles
IceFloeTracker.get_tiles
IceFloeTracker.get_trajectory_heads
IceFloeTracker.get_unmatched
IceFloeTracker.getbboxcolumns
IceFloeTracker.getbestmatchdata
IceFloeTracker.getcentroid
IceFloeTracker.getcentroidcolumns
IceFloeTracker.getcollisions
IceFloeTracker.getcollisionslocs
IceFloeTracker.getfit
IceFloeTracker.getfloemasks
IceFloeTracker.getidxmostminimumeverything
IceFloeTracker.getidxofrow
IceFloeTracker.getpropsday1day2
IceFloeTracker.getsizecomparability
IceFloeTracker.grad
IceFloeTracker.grad
IceFloeTracker.hbreak
IceFloeTracker.hbreak!
IceFloeTracker.histeq
IceFloeTracker.imadjust
IceFloeTracker.imbrighten
IceFloeTracker.imextendedmin
IceFloeTracker.imgradientmag
IceFloeTracker.imhist
IceFloeTracker.impose_minima
IceFloeTracker.imregionalmin
IceFloeTracker.imsharpen
IceFloeTracker.imsharpen_gray
IceFloeTracker.isfloegoodmatch
IceFloeTracker.isincountourlist
IceFloeTracker.ismember
IceFloeTracker.kmeans_segmentation
IceFloeTracker.loadimg
IceFloeTracker.long_tracker
IceFloeTracker.make_filename
IceFloeTracker.make_hbreak_dict
IceFloeTracker.make_lut
IceFloeTracker.make_psi_s
IceFloeTracker.make_psi_s
IceFloeTracker.makecontourstartatP
IceFloeTracker.makeemptydffrom
IceFloeTracker.makeemptyratiosdf
IceFloeTracker.matchcorr
IceFloeTracker.mean
IceFloeTracker.mismatch
IceFloeTracker.modenan
IceFloeTracker.morph_fill
IceFloeTracker.move
IceFloeTracker.norm
IceFloeTracker.normalize_image
IceFloeTracker.normalizeangle
IceFloeTracker.padnhood
IceFloeTracker.pairfloes
IceFloeTracker.reconstruct
IceFloeTracker.regionprops
IceFloeTracker.regionprops_table
IceFloeTracker.regularize_fill_holes
IceFloeTracker.regularize_sharpening
IceFloeTracker.remove_padding
IceFloeTracker.renamecols!
IceFloeTracker.resample_boundary
IceFloeTracker.reset_id!
IceFloeTracker.rgb2gray
IceFloeTracker.rgb2gray
IceFloeTracker.roundtoint!
IceFloeTracker.segmentation_A
IceFloeTracker.segmentation_B
IceFloeTracker.segmentation_F
IceFloeTracker.segmented_ice_cloudmasking
IceFloeTracker.sort_floes_by_area!
IceFloeTracker.timestamp
IceFloeTracker.trackercond1
IceFloeTracker.trackercond2
IceFloeTracker.trackercond3
IceFloeTracker.unsharp_mask
IceFloeTracker.unsharp_mask
IceFloeTracker.update!
IceFloeTracker.watershed_ice_floes
IceFloeTracker.watershed_product
IceFloeTracker.@persist