API documentation

This section of the documentation details the auto-generated API from Read the Docs.

class RADWave.altiwave.waveAnalysis(cycloneCSV=None, altimeterURL=None, bbox=None, stime=None, etime=None, satNames=None)[source]

Bases: object

RadWave is a package to perform wave analysis using AODN Radar Altimeter Data. Radar altimeters are mounted on satellites and measure a footprint of the ocean directly under the satellite (approx. 5-7 km wide). When the water surface is calm and flat, the reflection of the radar pulse back to the altimeter is almost instantaneous. In contrast, when waves are present the pulse is first reflected at the crest of the wave, then progressively reflected as the pulse reaches the wave trough.

Method:
Calibrated altimeter dataset provide wave height and wind speed through a series of equations described in [RibalYoung2019]
[RibalYoung2019]Ribal, A. & Young, I. R. - 33 Years of Globally Calibrated Wave Height and Wind Speed Data Based on Altimeter Observations. Sci Data 6, 77, DOI: 10.1038/s41597-019-0083-9, 2019.
Data:
Altimeter data are read from netcdf files that are provided via a list of URLs. Here we propose a series of examples using the Australian Ocean Data Network (AODN) (https://portal.aodn.org.auIMOS) but other web portals could be easely added.
Parameters:
  • cycloneCSV (str) – path as a csv file containing cyclone tracks with in the header the following names lon, lat & datetime [default: None]
  • altimeterURL (str) – list of NetCDF URLs downloaded from the wave data portal containing the radar altimeter data [default: None]
  • bbox (list) – bounding box specifying the geographical extent of the uploaded dataset following the convention [lon min,lon max,lat min,lat max] [default: None]
  • stime (list) – starting time of wave climate analysis following the convention [year, month, day] [default: None]
  • etime (list) – ending time of wave climate analysis following the convention [year, month, day] [default: None]
  • satNames (list) – list of satellites to use for the analysis - AODN portal provide the record from 10 satellites for altimeter data [default: None]

Note

This remote sensing technique does not provide individual wave statistics, but rather returns the average value over satellite footprint, up to 7 km wide.

close2Track(radius=2.0, dtmax=6)[source]

From all cyclone tracks, this function finds the closest processed altimeter geographical locations that have been recorded in the database. In addition to their coordinates, the altimeter dataset has to be recorded during a user defined time lapse close enough to the cyclone path time at each position.

The recorded dataframe cyclone_data contains the following variables:

  • altimeter significant wave height (m) ‘wH’
  • altimeter wave period (s) ‘period’
  • altimeter wave energy flux (kW/m) ‘power’
  • altimeter wave average energy density (J/m2) ‘energy’
  • altimeter wave group velocity (m/s) ‘speed’
  • distance between altimeter coordinates and cyclone path (km) ‘dist’
  • cyclone date (datetime) ‘cdate’
  • difference in time between recorded cyclone date and altimeter data for specific position (hours) ‘hours’
  • cyclone latitude position ‘clat’
  • cyclone longitude position ‘clon’
  • altimeter data latitude position ‘lat’
  • altimeter data longitude position ‘lon’
Parameters:
  • radius (float) – maximum radius distance in degree between cyclone position and altimeter data coordinates [default: ‘2.’]
  • dtmax (float) – maximum difference in time between recorded cyclone date and picked altimeter data (hours) [default: ‘6.’]

Note

The class waveAnalysis() saves a pandas dataframe called cyclone_data that stores the closest points to each cyclone path coordinates for further analysis.

Returns:pandas dataframe containing cyclone track information
Return type:cycdata (dataframe)
computeSeasonalCharacteristics(series='wh', time=None, lonlat=None, fsave=None, plot=True)[source]

This function allows the user to analyse the seasonal characteristics of each parameter obtained from the altimeter dataset. For a specified time interval and geographical extent, it computes for a given wave variable the monthly seasonality. Obtained monthly averaged values are stored and returned with a Pandas dataframe. User has the option to plot the computed wave paraneter characteristics as a heatmap, a box plot and a standard deviation graph.

For the wave height series, a Seasonal Mann-Kendall test is also performed to determine monotonic trends in computed dataset using the package from Hussain & Mahmud (2019).

Hussain & Mahmud, 2019: pyMannKendall: a python package for non parametric Mann Kendall family of trend tests - JOSS, 4(39), 1556.

Note

This function uses the seaborn and matplotlib libraries to plot the heatmap and the seasonal boxplot.

Parameters:
  • series (str) – name of the series to plot based on timeseries dataframe, choices are: ‘wh’, ‘period’, ‘power’, ‘energy’ and ‘speed’ [default: ‘wh’]
  • time (list) – extent of years to plot for plotting time series from 1995 to 2010 user will set ‘time’ to [1995,2010] [default: None]
  • lonlat (list) – specifying the geographical extent of the season characteristics computation following the convention [lon min,lon max,lat min,lat max] [default: None]
  • fsave (str) – saved image name without extension that will be written as a PNG file [default: None]
  • plot (bool) – flag specifying if plots have to been done [plot: True]
Returns:

pandas dataframe containing for chosen series variable the seasonality parameter for specified time interval

Return type:

dfseason (dataframe)

generateTimeSeries(days=30)[source]

Time series of wave characteristics are obtained via both the significant wave height and wind speed parameters from the altimeter dataset.

Time series compute the following instantaneous and monthly variables:

  • significant wave height (m) ‘wh’ & ‘wh_rolling’
  • wave period (s) ‘period’ & ‘period_rolling’
  • wave energy flux (kW/m) ‘power’ & ‘power_rolling’
  • wave average energy density (J/m2) ‘energy’ & ‘energy_rolling’
  • wave group velocity (m/s) ‘speed’ & ‘speed_rolling’
Parameters:days (int) – moving average window for time series calculation [default: ‘30’]

Note

The class waveAnalysis() saves a pandas dataframe called timeseries that stores the wave time series for further analysis.

Returns:pandas dataframe containing time series of wave characteristics
Return type:timeseries (dataframe)
meanEnergy(H, grav=9.80665, seadensity=1025.0)[source]

In a sea state, the average(mean) energy density per unit area of gravity waves on the water surface is proportional to the wave height squared, according to linear wave theory [Goda 2000].

Method:

Goda Y., 2000: Random Seas and Design of Maritime Structures - World Scientific, ISBN 978-981-02-3256-6.

Holthuijsen L.H., 2007: Waves in oceanic and coastal waters - Cambridge: Cambridge University Press, ISBN 978-0-521-86028-4.

Parameters:
  • H (numpy array) – significant wave height in metres
  • grav (float) – acceleration of gravity (m/s2) [default: 9.80665]
  • seadensity (float) – sea density (kg/m3) [default: 1025.]
Returns:

average(mean) energy density per unit area of gravity waves (J/m2)

Return type:

E (numpy array)

plotCycloneAltiPoint(showinfo=False, extent=None, addcity=None, markersize=100, zoom=4, fsize=(12, 10))[source]

This function plots a series of figures of the geographical coordinates for processed altimeter data close to each cyclone path position.

Parameters:
  • showinfo (bool) – show considered cyclone path coordiantes and associated recorded time and associate wave parameters [default: False]
  • extent (int) – geographical extent of the figure following the convention [lon min,lon max,lat min,lat max] [default: None]
  • addcity (list) – defined a specific location using the following convention [‘City Name’, longitude, latitude] [default: None]
  • markersize (int) – size of the markers used to specify cyclone tracks [default: 100]
  • zoom (int) – given target zoom for the target domain background (must be a value >= 1) [default: 4]
  • fsize – size of the image [default: (12,10)]

Note

This function relies on cartopy and matplotlib libraries. We use cartopy’s ability to draw map tiles which are downloaded on demand from the Stamen tile server.

plotCycloneTracks(title='Cyclone data tracks', markersize=100, zoom=4, extent=None, fsize=(12, 10), fsave=None)[source]

This function plots and saves in a figure a specific cyclone track loaded during the initialisation phase of the waveAnalysis class.

Example of cyclone track for Australia can be found on the Australian Bureau of Meteorology (BOM).

Parameters:
  • title (str) – title for the figure [default: “Cyclone data tracks”]
  • markersize (int) – size of the markers used to specify cyclone tracks [default: 100]
  • zoom (int) – given target zoom for the target domain background (must be a value >= 1) [default: 4]
  • extent (int) – geographical extent of the figure following the convention [lon min,lon max,lat min,lat max] [default: None]
  • fsize – size of the image [default: (12,10)]
  • fsave (str) – saved image name without extension that will be written as a PNG file [default: None]

Note

This function relies on cartopy and matplotlib libraries. We use cartopy’s ability to draw map tiles which are downloaded on demand from the Stamen tile server.

Warning

The cyclone tracks contains the following variables lon, lat & datetime.

plotTimeSeries(time='all', series='H', fsize=(12, 5), fsave=None)[source]

This function plots and saves in a figure a time series for a specific parameter from the processed altimeter data.

The choices of time series to plot are:

  • ‘H’: for wave height,
  • ‘T’: for wave period,
  • ‘P’: for wave power,
  • ‘E’: for wave energy, and
  • ‘Cg’: for wave group velocity

It also provides the following information:

  • Maximum parameter value
  • Mean parameter value
  • Median parameter value
  • 95th percentile parameter value
Parameters:
  • time (list) – extent of years to plot for plotting time series from 1995 to 2010 user will set ‘time’ to [1995,2010], to plot the entire record use the keyword ‘all’ [default: ‘all’]
  • series (str) – name of the series to plot choices are: ‘H’, ‘T’, ‘P’, ‘E’ and ‘Cg’ [default: ‘H’]
  • fsize – size of the image [default: (12,5)]
  • fsave (str) – saved image name without extension that will be written as a PNG file [default: None]
processAltimeterData(altimeter_pick='all', saveCSV='altimeterData.csv')[source]

From the list of OPeNDAP data URL’s this function extracts the altimeter data information.

The function can take some times to execute depending on the number of files to load and the size of the dataset to query.

Note

This function relies mostly on Pandas (library) and writes the processed dataset to file that can be later used to access more efficiently altimeter information.

Warning

Because the data is accessed via the THREDDS catalog, there can be a latency when querying the NetCDF files.

Parameters:
  • altimeter_pick (list) – list of satellites to use for the analysis as an example AODN portal provide the record from 10 satellites for altimeter data [default: ‘all’]
  • saveCSV (str) – filename used to save processed altimeter data obtained from the OPeNDAP web service [default: ‘altimeterData.csv’]
readAltimeterData(saveCSV=None)[source]

In case where the processAltimeterData function has already been executed, one can load directly the processed data from the created CSV file.

Parameters:saveCSV (str) – filename used to save processed altimeter data obtained from the OPeNDAP web service [default: None]
visualiseData(title='Altimeter data coordinates', extent=None, addcity=None, markersize=100, zoom=4, fsize=(12, 10), fsave=None)[source]

This function plots and saves in a figure the geographical coordinates of processed altimeter data. In case where a cyclone track has also been loaded during the initialisation phase of the waveAnalysis class, the cyclone path is also added to the figure.

Parameters:
  • title (str) – title for the figure [default: “Altimeter data coordinates”]
  • extent (int) – geographical extent of the figure following the convention [lon min,lon max,lat min,lat max] [default: None]
  • addcity (list) – defined a specific location using the following convention [‘City Name’, longitude, latitude] [default: None]
  • markersize (int) – size of the markers used to specify cyclone tracks [default: 100]
  • zoom (int) – given target zoom for the target domain background (must be a value >= 1) [default: 4]
  • fsize – size of the image [default: (12,10)]
  • fsave (str) – saved image name without extension that will be written as a PNG file [default: None]

Note

This function relies on cartopy and matplotlib libraries. We use cartopy’s ability to draw map tiles which are downloaded on demand from the Stamen tile server.

waveEnergyFlux(H, T, grav=9.80665, seadensity=1025.0)[source]

The rate at which energy is carried by waves is determined using total wave energy and wave energy speed (kW/m)

Parameters:
  • H (numpy array) – significant wave height in metres
  • T (numpy array) – wave period in seconds
  • grav (float) – acceleration of gravity (m/s2) [default: 9.80665]
  • seadensity (float) – sea density (kg/m3) [default: 1025.]
Returns:

wave energy flux (kW/m)

Return type:

P (numpy array)

waveGroupVelocity(T, grav=9.80665)[source]

As the waves propagate, their energy is transported and the energy transport velocity is the group velocity. From Airy wave theory, the group velocity (cg) is calculated in (m/s).

Parameters:
  • T (numpy array) – wave period in seconds
  • grav (float) – acceleration of gravity (m/s2) [default: 9.80665]
Returns:

group velocity (m/s)

Return type:

cg (numpy array)

wavePeriod(H, U, grav=9.80665)[source]

This function computes the wave period from the following altimeter data: wave height and surface wind speed.

Method:

Altimeters do not directly measure wave period, an essential component to adequately characterise wave climate. To calculate wave period from altimeter backscatter we use the algorithm proposed by Remya et al. (2010) based on wind speed (U) and significant wave height (H). This method has a high accuracy of estimating wave period (T) in both wind and swell-generated seas (root mean square error = 0.76 s).

Remya et al., 2010: Altimeter-derived ocean wave period using genetic algorithm - IEEE Geoscience and Remote Sensing Letters 8(2), 354–358.

Parameters:
  • H (numpy array) – significant wave height in metres
  • U (numpy array) – surface wind speed in metres per second
  • grav (float) – acceleration of gravity (m/s2) [default: 9.80665]
Returns:

wave period following Remya et al. (2010) algorithm (second)

Return type:

period (numpy array)