SARProcessing

Sentinel1BurstInformation

returns structure of Sentinel1BurstInformation from metadata in .xml Sentinel1BurstInformation contain information from Sentinel1DopplerCentroid and Sentinel1AzimuthFmRate

" Sentinel1BurstInformation

Constructors for the Sentinel1BurstInformation structure.

It takes a dictionary containing the full sentinel-1 swath metadata and extracts the Sentinel1BurstInformation specific data for a single burst as a structure.

Input

• meta_dict[dict]: a dictionary of the metadata.burst_number
• burst_number[Int]: Integer value of burst number.

Output

• Sentinel1BurstInformation[structure of Sentinel1BurstInformation]

Sentinel1GeolocationGrid

returns structure of Sentinel1GeolocationGrid from metadata in .xml

" Sentinel1GeolocationGrid

Constructors for the Sentinel1GeolocationGrid structure.

It takes a dictionary containing the full sentinel-1 swath metadata and extracts the Sentinel1GeolocationGrid as a structure. Input in the Sentinel1GeolocationGrid file:

• lines: Reference image MDS line to which this geolocation grid point applies.
• samples,
• latitude: Geodetic latitude of grid point [degrees].
• longitude: Geodetic longitude of grid point [degrees].
• azimuth_time: Zero Doppler azimuth time to which grid point applies [UTC].
• slant_range_time_seconds: Two-way slant range time to grid point.
• elevation_angle: Elevation angle to grid point [degrees].
• incidence_angle: Incidence angle to grid point [degrees].
• height: Height of the grid point above sea level [m].

Example

# accesing the geolocation data from the full metadata.
xmlPath = "s1a-iw1-slc-vh-20220220t144146-20220220t144211-041998-050092-001.xml"
Metadata1 = Sentinel1MetaData(xmlPath)
geolocation = Metadata1.geolocation

# accessing the geolocation directly from the xml.
meta_dict = read_xml_as_dict(xmlPath)
geolocation = Sentinel1GeolocationGrid(meta_dict)

Input

• meta_dict[dict]: a dictionary of the metadata.

Output

• Sentinel1GeolocationGrid[structure of Sentinel1GeolocationGrid]

Sentinel1Header

returns structure of Sentinel1Header from metadata in .xml

" Sentinel1Header

Constructors for the Sentinel1Header structure.

It takes a dictionary containing the full sentinel-1 swath metadata and extracts the Sentinel1Header as a structure. Input in the header file:

• missionId: Mission identifier for this data set.
• productType: Product type for this data set.
• polarisation: Polarisation for this data set.
• mission_data_take_id: Mission data take identifier.
• swath: Swath identifier for this data set. This element identifies the swath that applies to all data contained within this data set. The swath identifier "EW" is used for products in which the 5 EW swaths have been merged. Likewise, "IW" is used for products in which the 3 IW swaths have been merged.
• mode: Sensor mode for this data set.
• start_time: Zero Doppler start time of the output image [UTC].
• stop_time: Zero Doppler stop time of the output image [UTC].
• absolute_orbit_number: Absolute orbit number at data set start time.
• image_number: Image number. For WV products the image number is used to distinguish between vignettes. For SM, IW and EW modes the image number is still used but refers instead to each swath and polarisation combination (known as the 'channel') of the data.

Input

• meta_dict[dict]: a dictionary of the metadata.

Output

• Sentinel1Header[structure of Sentinel1Header]

Sentinel1ImageInformation

returns structure of Sentinel1ImageInformation from metadata in .xml

" Sentinel1ImageInformation

Constructor for the Sentinel1ImageInformation structure.

It takes a dictionary containing the full sentinel-1 swath metadata and extracts the Sentinel1ImageInformation as a structure. Input in the Sentinel1ImageInformation file:

• range_pixel_spacing: Pixel spacing between range samples [m].
• azimuth_frequency: Azimuth line frequency of the output image [Hz]. This is the inverse of the azimuth_timeInterval.
• slant_range_time_seconds: Two-way slant range time to first sample.
• incidence_angle_mid_swath: Incidence angle at mid swath [degrees].
• azimuth_pixel_spacing: Nominal pixel spacing between range lines [m].
• number_of_samples: Total number of samples in the output image (image width).

Input

• meta_dict[dict]: a dictionary of the metadata.

Output

• Sentinel1ImageInformation[structure of Sentinel1ImageInformation]

Sentinel1MetaData

Metadata structure for the Sentinel-1 satellite for each burst in the swath.

General metadata info is kept in the following structures:

• Sentinel1Header
• Sentinel1ProductInformation
• Sentinel1ImageInformation
• Sentinel1SwathTiming
• Sentinel1GeolocationGrid

Sentinel1BurstInformation specific Info is kept in

• Vector{Sentinel1BurstInformation}

Example

slcMetadata = Sentinel1MetaData(meta_dict)

Input

• meta_dict: xml file.

Can be accessed as, e.g., slcMetadata.product.radar_frequency –> 5.40500045433435e9::Float64 slcMetadata.header.swath –> 1::Int slcMetadata.header.mode –> "IW"::String slcMetadata.header.polarisation –> "VH"::String

" Sentinel1MetaData

Constructors for the Sentinel1MetaData structure.

It takes a Sentinel-1 single swath metafile in .xml format and constructs the metadata structure using the individual sub-structures in the metadata. The individual sub-structures in the metadata are:

• Sentinel1Header
• Sentinel1ProductInformation
• Sentinel1ImageInformation
• Sentinel1SwathTiming
• Sentinel1BurstInformation
• Sentinel1GeolocationGrid

Input

• xmlFile[string]: path of swath specific metadata in xml.format.

Output

• Sentinel1MetaData[structure of Sentinel1MetaData]: Structure with all Sentinel-1 metadata for a swath.

Example

Getting the t0 for the 5th burst: metadata = Sentinel1MetaData(annotation.xml) metadata.bursts.azimuthFmRate[5].t0

Sentinel1Polynomial

Sentinel1ProductInformation

returns structure of product information

" Sentinel1ProductInformation

Constructors for the Sentinel1ProductInformation structure.

It takes a dictionary containing the full sentinel-1 swath metadata and extracts the Sentinel1ProductInformation as a structure. Sentinel1ProductInformation file:

• pass: Direction of the orbit (ascending, descending)
• timeliness_category: Timeliness category under which the product was produced, i.e. time frame from the data acquisition
• platform_heading: Platform heading relative to North [degrees].
• projection: Projection of the image, either slant range or ground range.
• range_sampling_rate: Range sample rate [Hz]
• radar_frequency: Radar (carrier) frequency [Hz]
• azimuth_steering_rate: Azimuth steering rate for IW and EW modes [degrees/s].

Input

• meta_dict[dict]: a dictionary of the metadata.

Output

• Sentinel1ProductInformation[structure of Sentinel1ProductInformation]

Sentinel1SwathTiming

returns structure of Sentinel1SwathTiming from metadata in .xml

" Sentinel1SwathTiming

Constructors for the Sentinel1SwathTiming structure.

It takes a dictionary containing the full sentinel-1 swath metadata and extracts the Sentinel1SwathTiming as a structure.

Input

• meta_dict[dict]: a dictionary of the metadata.

Output

• Sentinel1SwathTiming[structure of Sentinel1SwathTiming]

TandemxDEM <: DEM

A DEM implementation for Tandem-X DEM data.

Computes Doppler centroid frequency (fetac), as given by equation 13 in the document "Definition of the TOPS SLC deramping function for products generated by the S-1 IPF" by Miranda (2014): https://sentinel.esa.int/documents/247904/1653442/sentinel-1-tops-slcderamping

Computes Doppler FM rate (ka), as given by equation 11 in the document "Definition of the TOPS SLC deramping function for products generated by the S-1 IPF" by Miranda (2014): https://sentinel.esa.int/documents/247904/1653442/sentinel-1-tops-slcderamping

" despeckled pixel value

Weighting coefficient for lee filter

Add padding to array

Q for simon:

approx_line_of_sight(orbit_state::OrbitState,incidence_angle_mid::Real)

Output

• line_of_sight::Array{float}(3): Line of sight to mid swath

#TODO, interpolate geolocationGridPoint from metadata instead?

azimuth_time2row(azimuth_time::Real,metadata::MetaData)

Returns the row corresponding to a specific azimuth time.

Note: That the burst overlap is not considered in this function. The actual image row will thus differ.

binarize_array(image::Matrix{T}, threshold::Real= 0.0001) where T<:Real

Making a mask with boolean values of image

Examples

bool_array = binarize_array(rand(10,10), threshold= 0.0001)

The two paramter (TP)-CFAR object detection method described in The State-of-the-Art in Ship Detection in Synthetic Aperture Radar imagery, D.J. Crips, 2004, in section 5.2 Adaptive threshold algorithms

Finding CFAR for all pixels in an image.

column2range(column::Real,metadata::MetaData)

Returns the range corresponding to the image column

" The The constant false alarm rate with convolution and pooling (CP-CFAR) object detection method described in: Z. Cui, H. Quan, Z. Cao, S. Xu, C. Ding and J. Wu, "SAR Target CFAR Detection Via GPU Parallel Operation," in IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, vol. 11, no. 12, pp. 4884-4894, Dec. 2018, doi: 10.1109/JSTARS.2018.2879082.

Convolution function copied from Yosi Pramajaya. Credits goes to him. In his blogpost, he showed this implementation was faster than many others.. See https://towardsdatascience.com/understanding-convolution-by-implementing-in-julia-3ed744e2e933

dont want to have too many packages. I therefore wont use Convolution pkg.

Input

input::Matrix{Float64}. The input image,
filter::Matrix{Float64}. The  filter/kernel to convolve
stride::Int64 = 1. Stride of the convolution.
padding::String = "valid". If padding is used ["valid" or "same"]

Output

result::Matrix{Float64}. convolved image.

Example

#define a filter.
average_pool_filter = filters.meanFilter([2,2])
#perform convolution.
image = operations.conv2d(image, average_pool_filter,2, "same")

ecef2SAR_index(
    ecef_coordinate::Array{T,1},
    interpolator,
    range_pixel_spacing::Real,
    azimuth_frequency::Real,
    near_range::Real,
    image_duration_seconds::Real
    ) where T <: Real

Convert ECEF-coordinates [X,Y,Z] to SARindex (rowfromfirstburst, image_column)

ecef2geodetic(ecef_coordinate::Array{Real,1};
                    semi_major_axis=6378137., flattening=1/298.257223563,
                    tolerance_latitude = 1.e-12, tolerance_height = 1.e-5)

Convert ECEF-coordinates [X,Y,Z] to geodetic-coordinates [latitude(radians),longitude(radians),height] (WGS-84) radians

(Based on B.R. Bowring, "The accuracy of geodetic latitude and height equations", Survey Review, v28 #218, October 1985 pp.202-206).

ellipsoid_intersect(x_sat::Array{Real,1},normalised_line_of_sight::Array{Real,1};
                            semi_major_axis::Real=6378137.,flattening::Real=1/298.257223563)

Computes the intersection between the satellite line of sight and the earth ellipsoid in ECEF-coordinates

Arguments

• x_sat::Array{Real,1}: [X,Y,Z] position of the satellite in ECEF-coordinates.
• normalised_line_of_sight::Array{Real,1}: Normalised Line of sight

Output

• x_0::Array{Real,1}: intersection between line and ellipsoid in ECEF-coordinates.

Note:

Equations follows I. Cumming and F. Wong (2005) p. 558-559

geodetic2SAR_index(geodetic_coordinate::Array{T,1}, interpolator, metadata::MetaData) where T <: Real

Convert geodetic-coordinates [latitude(radians),longitude(radians),height] to SARindex (rowfromfirstburst, image_column)

geodetic2ecef(geodetic_coordinate::Array{Real,1}; semi_major_axis::Real=WGS84_SEMI_MAJOR_AXIS,
    flattening::Real=WGS84_FLATTENING)

Convert geodetic-coordinates [latitude(radians),longitude(radians),height] (WGS-84) to ECEF-coordinates [X,Y,Z]

get_burst_mid_times(image::T) where T <: SingleLookComplex

Returns a vector of the mid burst times for the burst in the full image

get_burst_numbers(image::Sentinel1SLC)

Returns list of burst included in the image subset

get_coordinate(dem::T,index::Tuple{Integer,Integer}) where T <: DEM

Get the coordinate for a certain index in the DEM. (All the DEM's are stored as a matrix of heights)

get_image_rows(meta_data::Sentinel1MetaData, row_from_first_burst)

Converts the row number representing a unique azimuth time, rowfromfirst_burst, to the number in the full image. Note that two results are returned when the row appears in two bursts

get_index(dem::T,coordinate::Tuple{Real,Real}) where T <: DEM

Get the index of the DEM corresponding to the coordinate. (All the DEM's are stored as a matrix of heights)

get_sentinel1_annotation_paths(safe_path::string)

Getting the paths for the annotation files for a SLC image using its .SAFE folder path.

Parameters

• safe_path::String: path of .SAFE folder for one image.

Returns

• annotationPaths::Vector: Vector of paths for annotation files in .SAFE folder

Extracting a subset from an image. The subset will be extraxted from the image row/column defined by coordinate and size subset_size.

Input

image: The image array
coordinate::Vector{Int64}. Center coordinate of subset, in image geometry.
subset_size::Vector{Int64}=[75,75]. Size of the subset.

Output

subset::Array{Float64, 3}. The three dimensional subset [rows,columns,dimensions.] with dimension=1 for an input 2D array.

get_window(image::Sentinel1SLC)

Returns the window of the complete Sentinel image covered the "image"

load_precise_orbit_sentinel1(path)

Loads a Sentinel 1 orbit file

Returns

• Vector{OrbitState}

load_tandemx_dem(tiffPath::String) -> TandemxDEM

Load a Tandem-X DEM tiff. The Tandem-X DEM's can be downloaded from https://download.geoservice.dlr.de/TDM90/. Note: Invalid values are replaced with missing

load_tiff(filepath::String, window=nothing; convertToDouble = true,flip = true)

Read a Sentinel 1 tiff file.

Examples

julia> filepath = "s1a-iw3-slc-vv-20220918t074921-20220918t074946-045056-056232-006.tiff"
julia> data = readSwathSLC(filepath, [(501,600),(501,650)]);
julia> typeof(data)
Matrix{ComplexF64}
julia> size(data)
(100,150)

Making a mask with 1/NaN values of image

Examples

array_with_nan = mask_array_nan(rand(10,10), 0.5)

" mean_filter

Creates a mean filter for a image

Examples

filter = meanFilter([3,3])
or
filter = meanFilter([3])

" Find center locations of objects in a binary image.

Example

object_centers = object_locations(binary_array)

orbit_state_interpolator(orbit_states::Vector{OrbitState}, image::SARImage,
    polynomial_degree::Integer=4, margin::Integer = 3 )

Create a polynomial interpolation function for orbit states valid in the time span from image start time to image end time.

Returns

Anonymous interpolation function. (Input: secondsfromt_start::Float64, Output: state::OrbitState)

phase_ramp(rows::Vector{T}, columns::Vector{T}, burst_number::Int64, v_s::Float64,
        k_psi::Float64, dc_coefficient::Vector{Float64},
        dc_tau_0::Float64, fm_coefficient::Vector{Float64},
        fm_tau_0::Float64, f_c::Float64, lines_per_burst::Int64,
        number_of_samples::Int64, delta_t_s::Float64,
        delta_tau_s::Float64, tau_0::Real, c=LIGHT_SPEED::Real) where T <: Real

Computes the phase ramp (phi) for the given burst number for input rows (lines) and columns (samples).

NOTES

reference: Equation numbers refer to the document "Definition of the TOPS SLC deramping function for products generated by the S-1 IPF" by Miranda (2014): https://sentinel.esa.int/documents/247904/1653442/sentinel-1-tops-slc_deramping

phase_ramp(rows::Vector{T}, columns::Vector{T},
        burst_number::Int64, mid_burst_speed::Float64, meta_data::Sentinel1MetaData) where T <: Integer

Extracts relevant parameters from metadata and calls phaseramp().

phaserampgrid(rows::AbstractRange, columns::AbstractRange, burstnumber::Int64, midburstspeed::Float64, metadata::Sentinel1MetaData)

Computes the phase ramp (phi) for the given burst number over the grid defined by rows::AbstractRange and columns::AbstractRange

range2column(range::Real,metadata::MetaData)

Returns the image column corresponding to the range

" readsentinel1burstinformation(metadict::AbstractDict)::Vector{Sentinel1BurstInformation}

Read all burst information from meta_dict.

row2azimuth_time(row_from_first_burst::Real,metadata::MetaData)

Returns the azimuth_time corresponding to a specific row (as counted from first burst ignoring burst overlap)

Note: That the burst overlap is not considered in this function.

sar2gray(data::AbstractArray; p_quantile = 0.85)

Maps the data to values between 0 and 1 and convert into a gray scaled image. The minimum data value is mapped to 0 and all values above the p_quantile is mapped to 1

sar_index2geodetic(row_from_first_burst,
    image_column,
    height,
    interpolator,
    metadata::MetaData)

Convert SARindex (rowfromfirstburst, image_column) to geodetic coordinates [latitude(radians),longitude(radians),height]

" search_directory(path, key)

Searching a directory for files with extension.

" sobelFilter(input::Matrix{Float64},stride::Int64 = 1,padding::String="same")::Matrix{Float64}

Creates a sobelFilter for a image using edgeVertical() and edgeHorizontal()

Examples

sobel_image = filters.sobelFilter(image)

solve_radar

Find the point that is range away from the satellite, orthogonal on the flight directions and "height" above the elipsiod using Newton_rhapsody method.

In the homogeneous areas, the ratio of the standard deviation to the mean is a good measure of speckle strength. For the filtered SAR images, this ratio is also frequently used to measure the amount of speckle reduction.

speckle_lee_filter(image::Matrix{T} where T<: Real ,filter_size::Vector{N} where N<:Integer =[3,3])::Matrix{Real}

Original speckle lee filter for SAR images, see 'Refined Filtering of Image Noise Using Local Statistics' (1980) The Lee filter is an adaptive filter in the sence that it used the local statistics to determine the amount of speckle filtering. In homogenous region, it will resemble a mean filter. In inhomogenous region, i.e., near cities or edges, it will do less filtering.

Inputs

• image::Matrix, Speckled input iamge. 2D matrix.
• filter_size::[integer,integer], Filter size of speckle filter

Outputs

• despeckled_image::Matrix, despeckled image

Examples

descpeckled_image = SARProcessing.speckle_lee_filter(speckle_image,[9,9]);
descpeckled_image = SARProcessing.speckle_lee_filter(speckle_image,[3,3]);

speckle_mean_filter(image::Matrix{T} where T<: Real ,filter_size::Vector{N} where N<:Integer =[3,3])::Matrix{Real}

Mean filter.

Also called BoxFilter for SAR speckle reduction. Sometimes even Average filter.

Examples

descpek_mean_3 = speckle_mean_filter(abseloute_image,[3,3])
descpek_mean_3 = speckle_mean_filter(abseloute_image,[9,9])

speckle_median_filter(image::Matrix{T} where T<: Real, filter_size::Vector{N} where N<:Integer =[3,3])::Matrix{Real}

Median filter.

The two paramter (TP)-CFAR object detection method described in The State-of-the-Art in Ship Detection in Synthetic Aperture Radar imagery, D.J. Crips, 2004, in section 5.2 Adaptive threshold algorithms

Finding CFAR for all pixels in an image.