Index

Base.findmaxMethod
findmax(C::ClimGrid; skipnan::Bool=false)

Return the maximum element of ClimGrid C and its index. If there are multiple maximal elements, then the first one will be returned. If any data element is NaN, this element is returned. The result is in line with max. As climate data can often have NaN values (irregular polygons, missing weather data), the option to skip NaNs is provided as a keyword argument.

source
Base.findminMethod
findmin(C::ClimGrid; skipnan::Bool=false)

Return the minimum element of ClimGrid C and its index. If there are multiple minimal elements, then the first one will be returned. If any data element is NaN, this element is returned. The result is in line with min. As climate data can often have NaN values (irregular polygons, missing weather data), the option to skip NaNs is provided as a keyword argument.

source
Base.mergeMethod
merge(A::ClimGrid, B::ClimGrid)

Combines two ClimGrid. Based on the AxisArrays method.

source
Base.writeMethod
write(C::ClimGrid, filename::String)

Write to disk ClimGrid C to netCDF file. Still experimental, some attributes are hardcoded.

source
ClimateTools.RXdayMethod
RXday(C::ClimGrid; days=5)

Annual maximum consecutive x days precipitation amount (defaults to days=5).

Let RRkj be the precipitation amount for the 5-day interval ending k, period j. Then maximum 5-day values for period j are:

Rx5dayj = max (RRkj)

source
ClimateTools.annualmaxMethod
annualmax(C::ClimGrid)

Annual maximum of array data.

Let data[i,j] be daily time serie on day i in year j. Extract the highest value for year j.

source
ClimateTools.annualmeanMethod
annualmean(C::ClimGrid)

Annual mean of array data.

Let data[i,j] be daily time serie on day i in year j. Calculate the mean value for year j.

source
ClimateTools.annualminMethod
annualmin(C::ClimGrid)

Annual minimum of array data.

Let data[i,j] be daily time serie on day i in year j. Extract the lowest value for year j.

source
ClimateTools.annualsumMethod
annualsum(C::ClimGrid)

Annual sum of array data.

Let data[i,j] be daily time serie on day i in year j. Sums daily values for year j.

source
ClimateTools.applymaskMethod
applymask(A::AbstractArray{N, n}, mask::AbstractArray{N, n})

Applies a mask on the array A. Return an AbstractArray{N, n}.

source
ClimateTools.approx_surfacepressureMethod

approxsurfacepressure(sealevelpressure::ClimGrid, orography::ClimGrid, daily_temperature::ClimGrid)

Returns the approximated surface pressure (sp) (Pa) using sea level pressure (psl) (Pa), orography (orog) (m), and daily mean temperature (tas) (K).

$sp = psl * 10^{x}$

where $x = \frac{-orog}{18400 * tas / 273.15}$

source
ClimateTools.biascorrect_extremesMethod
biascorrect_extremes(obs::ClimGrid, ref::ClimGrid, fut::ClimGrid; detrend=false, window::Int=15, rankn::Int=50, thresnan::Float64=0.1, keep_original::Bool=false, interp=Linear(), extrap=Flat(), gev_params::DataFrame, frac=0.25, power=1.0)

Combine the empirial Quantile-Quantile mapping (see qqmap) and Generalized Pareto Distribution bias-correction methods.

The tail of the distribution is corrected with a paramatric GPD, using provided parameters μ, σ and ξ contained in gevparams DataFrame. The DataFrame gevparams has column :lat, :lon, :mu, :sigma and :xi, representing the GEV parameters and the spatial location of the parameters. For each grid point, the function extracts the closest GEV parameters.

Options specific to GPD

gevparams::DataFrame represents the (external) GEV parameters to correct each grid-points.

frac=0.25 is the fraction where the cutoff happens between QQM and GPD, as defined by (maximum(x) - minimum(x))*frac (e.g. For a maximum value of 150mm and a minimum value of 30mm, the linear transition will be between 30mm and 60mm).

power=1.0 is the shape of the transition.

Options specific to QQM

detrend::Bool = false (default). A 4th order polynomial is adjusted to the time series and the residuals are corrected.

window::Int = 15 (default). The size of the window used to extract the statistical characteristics around a given julian day.

rankn::Int = 50 (default). The number of bins used for the quantile estimations. The quantiles uses by default 50 bins between 0.01 and 0.99. The bahavior between the bins is controlled by the interp keyword argument. The behaviour of the quantile-quantile estimation outside the 0.01 and 0.99 range is controlled by the extrap keyword argument.

thresnan::Float64 = 0.1 (default). The fraction is missing values authorized for the estimation of the quantile-quantile mapping for a given julian days. If there is more than treshnan missing values, the output for this given julian days returns NaNs.

keep_original::Bool = false (default). If keep_original is set to true, the values are set to the original values in presence of too many NaNs.

interp = Interpolations.Linear() (default). When the data to be corrected lies between 2 quantile bins, the value of the transfer function is linearly interpolated between the 2 closest quantile estimation. The argument is from Interpolations.jl package.

extrap = Interpolations.Flat() (default). The bahavior of the quantile-quantile transfer function outside the 0.01-0.99 range. Setting it to Flat() ensures that there is no "inflation problem" with the bias correction. The argument is from Interpolation.jl package.

See also: qqmap

source
ClimateTools.customthresoverMethod
customthresover(C::ClimGrid, thres)

customthresover, annual number of days over a specified threshold.

Let TS[i,j] be a daily time serie value on day i in year j. Count the number of days where:

TS[i,j] > thres.

source
ClimateTools.customthresunderMethod
customthresunder(C::ClimGrid, thres)

customthresunder, annual number of days under a specified threshold.

Let TS[i,j] be a daily time serie value on day i in year j. Count the number of days where:

TS[i,j] < thres.

source
ClimateTools.diurnaltemperatureMethod

diurnaltemperature(temperatureminimum::ClimGrid, temperaturemaximum::ClimGrid, α::Float64)

Returns an estimation of the diurnal temperature (temperature between 7:00 (7am) and 17:00 (5pm)). The estimation is a linear combination of the daily minimum temperature (temperatureminimum) and daily maximum temperature (temperaturemaximum). The value of α has to be estimated seperatly from observations and depends on the location. The daily max and min must be in the same unit and in Celsius or Kelvin The diurnal temperature returned is in the same units as the daily minimum temperature and daily maximum temperature.

$Tdiu = α * Tmin + (1 - α) * Tmax$

source
ClimateTools.drought_dcMethod

drought_dc(pr::ClimGrid, tas::Climgrid)

Returns the drought index. The index is defined as a function of daily temperature and daily precipitation. This indice correspond to the DC indice in Wang et al. 2015.

Reference: Wang et al. 2015. Updated source code for calculating fire danger indices in the Canadian Forest Fire Weather Index System. Information Report NOR-X-424. Canadian Forst Service. 36pp.

source
ClimateTools.drought_dcMethod

drought_dc(prvec::Array{N,1}, tasvec::Array{N,1}, timevec)

Returns the drought index. The index is defined as a function of daily temperature and daily precipitation. This indice correspond to the DC indice in Wang et al. 2015.

Reference: Wang et al. 2015. Updated source code for calculating fire danger indices in the Canadian Forest Fire Weather Index System. Information Report NOR-X-424. Canadian Forst Service. 36pp.

source
ClimateTools.frostdaysMethod
frostdays(C::ClimGrid)

FD, Number of frost days: Annual count of days when TN (daily minimum temperature) < 0 Celsius.

Let TN[i,j] be daily minimum temperature on day i in year j. Count the number of days where:

TN[i,j] < 0 Celsius.

source
ClimateTools.get_dimnameMethod
get_dimname(::NCDatasets.Dataset, dim::String)

Scans dataset ds and extract the dimensions X, Y, Z and T. The dataset needs to be CF-compliant. See http://cfconventions.org/.

source
ClimateTools.getdim_latMethod
getdim_lat(ds::NCDatasets.Dataset)

Returns the name of the "latitude" dimension and the status related to a regular grid. The latitude dimension is usually "latitude", "lat", "y", "yc", "rlat".

source
ClimateTools.getdim_lonMethod
getdim_lon(ds::NCDatasets.Dataset)

Returns the name of the "longitude" dimension and the status related to a regular grid. The longitude dimension is usually "longitue", "lon", "x", "xc", "rlon".

source
ClimateTools.griddataMethod
C = griddata(A::ClimGrid, londest::AbstractArray{N, 1} where N, latdest::AbstractArray{N, 1} where N)

Interpolate ClimGrid A onto lat-lon grid defined by londest and latdest vector or array. If an array is provided, it is assumed that the grid is curvilinear (not a regular lon-lat grid) and the user needs to provide the dimension vector ("x" and "y") for such a grid.

source
ClimateTools.griddataMethod
C = griddata(A::ClimGrid, B::ClimGrid; min=[], max=[])

Interpolate ClimGrid A onto the lon-lat grid of ClimGrid B, where A and B are ClimGrid.

Min and max optional keyword are used to constraint the results of the interpolation. For example, interpolating bounded fields can lead to unrealilstic values, such as negative precipitation. In that case, one would use min=0.0 to convert negative precipitation to 0.0.

source
ClimateTools.icingdaysMethod
icingdays(C::ClimGrid)

ID, Number of summer days: Annual count of days when TX (daily maximum temperature) < 0 degree Celsius.

Let TX[i,j] be daily maximum temperature on day i in year j. Count the number of days where:

TX[i,j] < 0 Celsius.

source
ClimateTools.inpolyMethod
inpoly(p, poly::Matrix)

Determines if a point is inside a polygon.

  • p – point (x,y) or [x,y]
  • poly – polygon vertices [x1 x2 ... xn x1 y1 y2 ... yn y1]

Returns true if point has an odd winding number. This should label points as exterior which are inside outcrops. See test for a test.

Author: Github "Mauro3" / "Mauro"

source
ClimateTools.inpolygridMethod
inpolygrid(lon, lat, poly::AbstractArray{N,2} where N)

Used to test a grid of points. Returns a mask of ones and NaNs of the same size as lon and lat.

source
ClimateTools.loadMethod
load(file::String, vari::String; poly = Array{Float64}([]), start_date::Tuple, end_date::Tuple, data_units::String = "")

Returns a ClimGrid type with the data in file of variable vari inside the polygon poly. Metadata is built-in the ClimGrid type, from the netCDF attributes.

Inside the ClimgGrid type, the data is stored into an AxisArray data type, with time, longitude/x and latitude/y dimensions.

The polygon provided should be in the -180, +180 longitude format. If the polygon crosses the International Date Line, the polygon should be splitted in multiple parts (i.e. multi-polygons).

Options for data_units are for precipitation : "mm", which converts the usual "kg m-2 s-1" unit found in netCDF files. For temperature : "Celsius", which converts the usual "Kelvin" unit.

Temporal subsetting can be done by providing start_date and end-date Tuples of length 1 (year), length 3 (year, month, day) or 6 (hour, minute, second).

Note: load uses CF conventions. If you are unable to read the netCDF file with load, the user will need to read it with low-level functions available in NetCDF.jl package or NCDatasets.jl or re-create standartized netCDF files.

source
ClimateTools.loadMethod
load(files::Array{String,1}, vari::String; poly = ([]), start_date::Date = Date(-4000), end_date::Date = Date(-4000), data_units::String = "")

Loads and merge the files contained in the arrar files.

source
ClimateTools.load2DMethod
load2D(file::String, vari::String; poly=[], data_units::String="")

Returns a 2D array. Should be used for fixed data, such as orography.

source
ClimateTools.meantemperatureMethod

meantemperature(temperatureminimum::ClimGrid, temperaturemaximum::ClimGrid)

Returns the daily mean temperature calculated from the maximum and minimum temperature. Daily maximum and minimum temperature must be in the same units. The mean temperature returned is in the same units as the daily minimum temperature and daily maximum temperature.

$Tmean = \frac{Tmax + Tmin}{2}$

source
ClimateTools.meshgridMethod
X, Y = meshgrid{T}(vx::AbstractVector{T}, vy::AbstractVector{T})

This function creates a 2-D mesh-grid in a format consistent with Matlab's function meshgrid(). XV and YV are vectors.

source
ClimateTools.ndgridMethod
X, Y = ndgrid(XV, YV)

This function creates a 2-D mesh-grid in a format consistent with Matlab's function ndgrid(). XV and YV are vectors.

source
ClimateTools.polyfitMethod
polyfit(C::ClimGrid)

Returns an array of the polynomials functions of each grid points contained in ClimGrid C.

source
ClimateTools.polyvalMethod
polyval(C::ClimGrid, polynomial::Array{Poly{Float64}})

Returns a ClimGrid containing the values, as estimated from polynomial function polyn.
source
ClimateTools.prcp1Method
prcp1(C::ClimGrid)

Annual number with preciptation >= 1 mm. This function returns a ClimGrid. Input data should be in mm.

source
ClimateTools.qqmapMethod

qqmap(obs::ClimGrid, ref::ClimGrid, fut::ClimGrid; method="Additive", detrend=true, window::Int=15, rankn::Int=50, thresnan::Float64=0.1, keep_original::Bool=false, interp::Function = Linear(), extrap::Function = Flat())

Quantile-Quantile mapping bias correction.

For each julian day of the year (+/- window size), a transfer function is estimated through an empirical quantile-quantile mapping. The quantile-quantile transfer function between ref and obs is etimated on a julian day (and grid-point) basis with a moving window around the julian day. Hence, for a given julian day, the transfer function is then applied to the fut dataset for a given julian day.

Options

method::String = "Additive" (default) or "Multiplicative". Additive is used for most climate variables. Multiplicative is usually bounded variables such as precipitation and humidity.

detrend::Bool = true (default). A 4th order polynomial is adjusted to the time series and the residuals are corrected with the quantile-quantile mapping.

window::Int = 15 (default). The size of the window used to extract the statistical characteristics around a given julian day.

rankn::Int = 50 (default). The number of bins used for the quantile estimations. The quantiles uses by default 50 bins between 0.01 and 0.99. The bahavior between the bins is controlled by the interp keyword argument. The behaviour of the quantile-quantile estimation outside the 0.01 and 0.99 range is controlled by the extrap keyword argument.

thresnan::Float64 = 0.1 (default). The fraction is missing values authorized for the estimation of the quantile-quantile mapping for a given julian days. If there is more than treshnan missing values, the output for this given julian days returns NaNs.

keep_original::Bool = false (default). If keep_original is set to true, the values are set to the original values in presence of too many NaNs.

interp = Interpolations.Linear() (default). When the data to be corrected lies between 2 quantile bins, the value of the transfer function is linearly interpolated between the 2 closest quantile estimation. The argument is from Interpolations.jl package.

extrap = Interpolations.Flat() (default). The bahavior of the quantile-quantile transfer function outside the 0.01-0.99 range. Setting it to Flat() ensures that there is no "inflation problem" with the bias correction. The argument is from Interpolation.jl package.

source
ClimateTools.qqmapMethod

qqmap(obsvec::Array{N,1}, refvec::Array{N,1}, futvec::Array{N,1}, days, obsjul, refjul, futjul; method::String="Additive", detrend::Bool=true, window::Int64=15, rankn::Int64=50, thresnan::Float64=0.1, keeporiginal::Bool=false, interp=Linear(), extrap=Flat())

Quantile-Quantile mapping bias correction for single vector. This is a low level function used by qqmap(A::ClimGrid ..), but can work independently.

source
ClimateTools.resampleMethod
resample(C::ClimGrid, startmonth::Int64, endmonth::Int64)

Return a resampled subset of ClimGrid C based on months.

source
ClimateTools.resampleMethod
resample(C::ClimGrid, season::String)

Return a resampled subset of ClimGrid C for a given season. Season options are: "DJF" (December-February), "MAM" (March-May), "JJA" (June-August), "SON" (September-November)

source
ClimateTools.shapefile_coords_polyMethod
shapefile_coords_poly(poly::Shapefile.Polygon)

Return the polygons contained in shp.shapes[i]. It returns an array containing the polygons.

See also shapefile_coords, which returns vectors as opposed to array. Returned polygon is consistent with the data extraction of the load function.

source
ClimateTools.spatialsubsetMethod
spatialsubset(C::ClimGrid, poly::Array{N, 2})

Returns the spatial subset of ClimGrid C. The spatial subset is defined by the polygon poly, defined on a -180, +180 longitude reference.

source
ClimateTools.summerdaysMethod
summerdays(C::ClimGrid)

SD, Number of summer days: Annual count of days when TX (daily maximum temperature) > 25 degree Celsius.

Let TX[i,j] be daily maximum temperature on day i in year j. Count the number of days where:

TX[i,j] > 25 Celsius.

source
ClimateTools.temporalsubsetMethod
function temporalsubset(C::ClimGrid, startdate::Date, enddate::Date)

Returns the temporal subset of ClimGrid C. The temporal subset is defined by a start and end date.

source
ClimateTools.timeindexMethod
timeindex(timeVec, start_date, end_date, T)

Return the index of time vector specified by startdate and enddate. T is the DateTime type (see NCDatasets.jl documentation).

source
ClimateTools.tropicalnightsMethod
tropicalnights(C::ClimGrid)

TropicalNights, Number of tropical nights: Annual count of days when TN (daily maximum temperature) > 20 degree Celsius.

Let TN[i,j] be daily minimum temperature on day i in year j. Count the number of days where:

TN[i,j] > 20 Celsius.

source
ClimateTools.vaporpressureMethod

vaporpressure(specifichumidity::ClimGrid, sealevelpressure::ClimGrid, orography::ClimGrid, daily_temperature::ClimGrid)

Returns the vapor pressure (vp) (Pa) estimated with the specific humidity (q), the sea level pressure (psl) (Pa), the orography (orog) (m) and the daily mean temperature (tas) (K). An approximation of the surface pressure is first computed by using the sea level pressure, orography and the daily mean temperature (see approx_surfacepressure). Then, vapor pressure is calculated by:

$vp = \frac{q * sp}{q+0.622}$

source
ClimateTools.vaporpressureMethod

vaporpressure(surfacepressure::ClimGrid, specifichumidity::ClimGrid)

Returns the vapor pressure (vp) (Pa) based on the surface pressure (sp) (Pa) and the specific humidity (q).

$vp = \frac{q * sp}{q+0.622}$

source
ClimateTools.wbgtMethod

wbgt(diurnaltemperature::ClimGrid, vaporpressure::ClimGrid)

Returns the simplified wet-bulb global temperature (wbgt) (Celsius) calculated using the vapor pressure (Pa) of the day and the estimated mean diurnal temperature (Celsius; temperature between 7:00 (7am) and 17:00 (5pm)).

$wbgt = 0.567 * Tday + 0.00393 * vp + 3.94$

source
ClimateTools.yearmonthdayhourMethod
yearmonthdayhour(dt::AbstractCFDateTime) -> (Int64, Int64, Int64, Int64)

Simultaneously return the year, month, day and hour parts of dt. Author: Alexander-Barth (Github)

source
ClimateTools.ClimGridType
ClimGrid{A <: AxisArray}

In-memory representation of Climate Forecast netCDF files. struct ClimGrid

data::AxisArray # Data

longrid::AbstractArray{N,2} where N # the longitude grid

latgrid::AbstractArray{N,2} where N # the latitude grid

msk::Array{N, 2} where N # Data mask (NaNs and 1.0)

grid_mapping::Dict#{String, Any} # bindings for native grid

dimension_dict::Dict

model::String

frequency::String # Day, month, years

experiment::String # Historical, RCP4.5, RCP8.5, etc.

run::String

project::String # CORDEX, CMIP5, etc.

institute::String # UQAM, DMI, etc.

filename::String # Path of the original file

dataunits::String # Celsius, kelvin, etc.

latunits::String # latitude coordinate unit

lonunits::String # longitude coordinate unit

variable::String # Type of variable (i.e. can be the same as "typeofvar", but it is changed when calculating indices)

typeofvar::String # Variable type (e.g. tasmax, tasmin, pr)

typeofcal::String # Calendar type

timeattrib::Dict # Time attributes (e.g. days since ... )

varattribs::Dict # Variable attributes dictionary

globalattribs::Dict # Global attributes dictionary

end

source
ClimateTools.ClimGridMethod
ClimGrid(data; longrid=[], latgrid=[], msk=[], grid_mapping=Dict(), dimension_dict=Dict(), model="NA", frequency="NA", experiment="NA", run="NA", project="NA", institute="NA", filename="NA", dataunits="NA", latunits="NA", lonunits="NA", variable="NA", typeofvar="NA", typeofcal="NA", varattribs=Dict(), globalattribs=Dict())

Constructor of the ClimGrid function. Data is an AxisArray. Everything else is optional, but usually needed for further processing (mapping, interpolation, etc...). struct ClimGrid

data::AxisArray # Data

longrid::AbstractArray{N,2} where N # the longitude grid

latgrid::AbstractArray{N,2} where N # the latitude grid

msk::Array{N, 2} where N # Data mask (NaNs and 1.0)

grid_mapping::Dict#{String, Any} # bindings for native grid

dimension_dict::Dict

model::String

frequency::String # Day, month, years

experiment::String # Historical, RCP4.5, RCP8.5, etc.

run::String

project::String # CORDEX, CMIP5, etc.

institute::String # UQAM, DMI, etc.

filename::String # Path of the original file

dataunits::String # Celsius, kelvin, etc.

latunits::String # latitude coordinate unit

lonunits::String # longitude coordinate unit

variable::String # Type of variable (i.e. can be the same as "typeofvar", but it is changed when calculating indices)

typeofvar::String # Variable type (e.g. tasmax, tasmin, pr)

typeofcal::String # Calendar type

timeattrib::Dict # Time attributes (e.g. days since ... )

varattribs::Dict # Variable attributes dictionary

globalattribs::Dict # Global attributes dictionary

end

source