Page History

...

Computes the geostrophic wind from geopotential fields defined on pressure levels.  For a given z geopotential field the computation of the geostrophic wind components is based on the following formulas:

...

The derivatives are computed with a second order finite-difference approximation. The resulting fieldset contains two fields for each input field: the zonal derivative followed by the meridional derivative. The output fields are bitmapped on the poles (they contain missing values there). Please note that this function is only implemented for regular latitude-longitude grids.

...

Please be aware that in versions before Metview 5.5.3 the resulting fields of gradient() appear in the wrong order.  The correct order should look like this:

...

.

...

but instead of it the following order is produced:

[gradx_f1, gradx_f2,..., grady_f1, grady_f2, ....]

For the efficient retrieval of multiple GRIB keys from a fieldset. A single call to grib_get can replace multiple calls to the other grib_get_* functions and is hence more efficient. The keys are provided as a list for the second argument; by default they will be retrieved as strings, but their type can be specified by adding a modifier to their names, following the convention used by grib_ls where the key name is followed by a colon and then one or two characters which specify the type (:

• • s=string

...

• • l=long

...

• • d=double

...

• • la=long array

...

• • da=double array
  • n=native type New in Metview version 5.14.0

For ). For example, the key 'centre' can be retrieved as a string with 'centre' or 'centre:s', or as a number with 'centre:l'. Each GRIB key has a ‘native type’, e.g. long or string. If the type is specified as “n” then the type that is returned. The native type for the key ‘centre’ is str, so ‘centre:n’ will return a str.

The result is always a list of lists; by default, or if the optional third argument is 'field', the result will be grouped by field, containing one list per field, each of these lists containing one element per key; if the optional third parameter is 'key', the result will be grouped by key, containing one list per key, each of these lists containing one element per field. Example - the following lines of Macro code on a particular 6-field fieldset:

...

f = grib_set(f, ["date", 20150601,       # integer
                 "time", 0600,           # integer
                 "stepType", "avg",      # string
                 "startStep", 0 ,        # integer
                 "endStep", 31,          # integer
                 "unitOfTimeRange", "D", # string
                 "longitudeOfLastGridPointInDegrees", 100.5])   100.5  #  double

...

This function computes the average of each a field in a fieldset over an area. The area of the individual cells can be is approximated with the following formula:

...

This formula is only used for (regular and reduced) latitude-longitude and Gaussian grids. For all other grid types integrate() simply returns the mathematical average of the values (just like the average() function does).

...

Interpolate a fieldset at a given point. The method used is bilinear interpolation. If a list is given, it must contain two numbers - latitude and longitude. If two numbers are given, the first is the latitude, the second the longitude. The field must be a gridded field. If the fieldset has only one field, a number is returned; otherwise a list is returned. Where it is not possible to generate a sensible value due to lack of valid data in the fieldset, nil is returned. Note that a similar function, nearest_gridpoint() , also exists.

geopoints interpolate ( fieldset,geopoints )

Generates a set of geopoints from a field. The first field of the input fieldset is used. The field is interpolated for each position of the geopoints given as a second parameter. The method used is bilinear interpolation. The output geopoints take their date, time and level from the fieldset. Where it is not possible to generate a sensible value due to lack of valid data in the fieldset, the internal geopoints missing value is used (this value can be checked for with the built-in variable geo_missing_value or removed with the function remove_missing_values). Note that a similar function, nearest_gridpoint() , also exists.

...

• • Takes the grid values in the first fieldset and uses them as index in the second fieldset. E.g. a grid value of n in the first fieldset, retrieves the corresponding grid point value of the (n-1)th field of the second fieldset (indexing is 0 based). The output fieldset is built from these values and has as many fields as the second first fieldset.
  • Takes the grid values in the first fieldset and uses them as index in the list - real numbers are truncated, not rounded. E.g. a grid value of n in the first fieldset, retrieves the (n-1)th list element (indexing is 0 based). The output fieldset is built from these values and has as many fields as the first fieldset.

...

For each field of the input fieldset, this function creates a field containing grid point values of 0 or 1 according to whether they are outside or inside a defined geographical area.

The list parameter must contain exactly four numbers representing a geographical area. These numbers should be in the order north, west, south and east (negative values for western and southern coordinates). Non-rectangular masks, and even convex masks can be created by using the operators and , or and not . To create the

If "missing" is specified as the third argument it will change the behaviour so that points outside the area will become missing values and points inside the area retain their original value. This option is new in Metview version 5.13.0.

Non-rectangular masks, and even convex masks can be created by using the operators and , or and not . To create the following mask :

First decompose into basic rectangles :

...

Merge several fieldsets. The same as the operator &. The output is a fieldset with as many fields as the total number of fields in all merged fieldsets. Merging with the value nil does nothing, and is used to initialise when building a fieldset from nothing.

...

...

fieldset   ml_to_hl (mfld: fieldset, z: fieldset, zs: fieldset, hlist: list, reflev: string, method: string, [fs_surf: fieldset])

Interpolates a fieldset on on model levels   (i.e. on hybrid or eta levels used by the IFS) onto onto height levels (in m) above sea or ground level.   At gridpoints where interpolation is not possible missing value is returned. This function has the following positional arguments:

• • mfld: the fieldset to be interpolated
  • z: the geopotential fieldset on model

...

• • levels (it must contain the same levels as

...

• • mfld

...

• • but the order of the levels can be different)
  • zs: the surface geopotential field (if the

...

• • reflev

...

• • argument is set to

...

• • "sea"

...

• • it should be set to

...

• • nil).
  • hlist: the list of target height levels (they can came in any given order)
  • reflev: specifies the reference level for the target heights. The possible values are

...

• • "sea"

...

• • and

...

• • "ground"

...

• • * method: specifies the interpolation method. The possible values are

...

• • "linear"

...

• • and

...

• • "log".

Please note that geopotential is not archived operationally on model levels in MARS at ECMWF.  To compute geopotential on model levels use Metview's mvl_geopotential_on_ml() function. The following example shows how to use function ml_to_hl() together with mvl_geopotential_on_ml() :

• • fs_surf: (optional) specifies the field values on the surface. With this it is possible to interpolate to target heights between the surface and the bottom-most model level. If fs_surf is a number it defines a constant fieldset. Only available when ref_level is "ground". New in Metview version 5.14.0.

At gridpoints where the interpolation is not possible a missing value is returned. It can happen when the target height level is below the bottom-most model level or the surface (when fs_surf is used) or above the top-most level. Please note that model levels we are dealing with in ml_to_hl are "full-levels" and the bottom-most model level does match the surface but it is above it. If you need to interpolate to height levels close to the surface use fs_surf.

The following example shows how to use function ml_to_hl() together with mvl_geopotential_on_ml():

# retrieve the data on model levels - 
# surface geopotential (zs) is only available in the first forecast step!
common_retrieve_params = ( type : "fc", levtype : "ml", step : 12, grid : [1.5,1.5] )
t = retrieve param : "t", levelist : [1, 'to', 137], common_retrieve_params)
q = retrieve param : "q", levelist : [1, 'to', 137], common_retrieve_params)
lnsp = retrieve( param : "lnsp", levelist : 1, common

# retrieve the data on model levels - surface geopotential (zs) is only available in the first forecast step!

common_retrieve_params = ( type : "fc", levtype : "ml", step : 12, grid : [1.5,1.5] )
t = retrieve param : "t", levelist : [1, 'to', 137], common_retrieve_params)
q = retrieve param : "q", levelist : [1, 'to', 137], common_retrieve_params)
lnsp = retrieve( param : "lnsp", levelist : 1, common_retrieve_params)
zs = retrieve( param : "z", levelist : 1, type : "fc", levtype : "ml", step : 0, grid : [1.5,1.5])

#   
# compute geopotential on model levels
z = mvl_geopotential_on_ml(t, q, lnsp, zs)
   

# interpolate the t field onto a list of height levels above sea level
hlevs = [1000, 2000, 3000, 4000, 5000]
th = ml_to_hl (t, z, nil, hlevs, "sea", "linear")

...

Computes geopotential on model levels.

Parameter t should be a fieldset of temperature on model levels in ascending numeric order (e.g. 1-137), q a fieldset of specific humidity on model levels in ascending numeric order, lnsp a field of log of surface pressure on model level 1, zs a field of geopotential on model level 1 (available from MARS). All fields must be GRIDDED gridpoint data - no spherical harmonics, and they must all be on the same grid, with the same number of points. The function mvl_geopotential_on_ml() assumes that there are no other dimensions contained in the data, e.g. all fields should have the same date and time.

The return value is a fieldset of geopotential defined on the model levels .

The code below illustrates how to use this function:

# retrieves analysis data on model levels

r = (date: -1, time: 12, levtype: "ml", grid: [1.5,1.5])
t    = retrieve(r,levelist: [1,"to",137],param: "t")
q    = retrieve(r,levelist: [1,"to",137],param: "q")
zs   = retrieve(r,levelist: 1,param: "z")
lnsp = retrieve(r,levelist: 1,param: "lnsp")


# computes the geopotential


z_ml = mvl_geopotential_on_ml(t, q, lnsp, zs)

...

present in the input data sorted by ascending numeric level order.

The required levels and their ordering in t and q depend on the Metview version used:

• • from Metview version 5.14.0: t and q must contain the same levels in the same order but there is no restriction on the actual level ordering. The model level range must be contiguous and must include the bottom-most level. E.g. if the current vertical coordinate system has 137 model levels using only a subset of levels between e.g. 137-96 is allowed.
  • in previous Metview versions: t and q must contain the full model level range in ascending numeric order. E.g. if the current vertical coordinate system has 137 model levels t and q must contain all the levels ordered as 1,..., 137.

The code below illustrates how to use this function:

# retrieves analysis data on model levels

r = (date: -1, time: 12, levtype: "ml",

...

fieldset mvl_ml2hPa(lnsp: fieldset, mfld: fieldset, plist: list)

Interpolates a fieldset currently on model levels onto pressure levels (in hPa). Locations where interpolation is not possible are returned as missing.

Parameter lnsp is a field of logarithm of surface pressure; mfld is the fieldset to be interpolated and should be on model levels; plist is a list of pressure levels in hPa - the result will be the mfld fieldset interpolated onto these levels. Neither mfld nor plist need to be sorted.

The following code shows a simple example:

# retrieve the data in model levels

common_retrieve_params = ( type : "fc", levtype : "ml", step : 12, grid: [1.5,1.5])
tmod t    = retrieve param : "t"(r,levelist: [1, '"to'", 91137],common_retrieve_paramsparam: "t")
lnsp q    = retrieve(param : "lnsp", r,levelist: [1,"to",137],param: "q")
zs   = retrieve(r,levelist: 1,common_retrieve_params)

# interpolate onto a list of pressure levels

plevels = [1000, 900, 850, 500, 300, 100, 10, 1, 0.1]

tpres = mvl_ml2hPa (lnsp, tmod, plevels)

...

param: "z")
lnsp = retrieve(r,levelist: 1,param: "lnsp")


# computes the geopotential


z_ml = mvl_geopotential_on_ml(t, q, lnsp, zs)

fieldset mvl_ml2hPa(lnsp: fieldset, mfld: fieldset, plist: list)

Interpolates a fieldset currently on model levels onto pressure levels (in hPa). Locations where interpolation is not possible are returned as missing.

Parameter lnsp is a field of logarithm of surface pressure; mfld is the fieldset to be interpolated and should be on model levels; plist is a list of pressure levels in hPa - the result will be the mfld fieldset interpolated onto these levels. Neither mfld nor plist need to be sorted.

The following code shows a simple example:

# retrieve the data in model levels

common_retrieve_params = ( type : "fc", levtype : "ml", step : 12, grid : [1.5,1.5] )

tmod = retrieve param : "t", levelist : [1, 'to', 91], common_retrieve_params)

lnsp = retrieve( param : "lnsp", levelist : 1, common_retrieve_params)

# interpolate onto a list of pressure levels

plevels = [1000, 900, 850, 500, 300, 100, 10, 1, 0.1]

tpres = mvl_ml2hPa (lnsp, tmod, plevels)

Returns the value of the nearest point to a given location (or locations) in each field of a fieldset. The field must be a gridded field. If a list is given, it must contain two numbers - latitude and longitude. If two numbers are given, the first is the latitude, the second the longitude. For batch processing of multiple locations, two vectors can be given, the first is a vector of latitudes, the second the longitudes; this can be much more efficient than multiple calls with a single location each. If the fieldset has only one field, a number (or vector) is returned; otherwise a list of numbers (or a list of vectors) is returned.

By default, when the nearest gridpoint value is a missing value or the location is out of the grid area, nil is returned in the case of a single coordinate, or vector_missing_value in the case of a vector. If an extra parameter 'valid' is added to the function call, then of the surrounding points, the nearest valid one is returned; nil will still be returned if all the surrounding points are missing.

Note that a similar function, interpolate(), also exists.

geopoints nearest_gridpoint ( fieldset,geopoints )

Generates a set of geopoints from a field. The first field of the input fieldset is used. The result is a set of geopoints whose values are those of the nearest gridpoints in the field to the geopoints given as a second parameter. Where it is not possible to generate a sensible value due to lack of valid data in the fieldset, the internal geopoints missing value is used (this value can be checked for with the built-in variable geo_missing_value or removed with the function remove_missing_values). Note that a similar function, interpolate() , also exists.

Anchornearest_gridpoint_infonearest_gridpoint_infolist nearest_gridpoint_info ( fieldset,list[,string] )
list nearest_gridpoint_info ( fieldset,number,number[,string] )


Returns the value

...

of the nearest point to a given location (or locations) in each field of a fieldset. The field must be a gridded field. If a list is given, it must contain two numbers - latitude and longitude. If two numbers are given, the first is the latitude, the second the longitude.

...

For batch processing of multiple locations, two vectors can be given, the first is a vector of latitudes, the second the longitudes; this can be much more efficient than multiple calls with a single location each. If the fieldset has only one field, a number (or vector) is returned; otherwise a list of numbers (or a list of vectors) is returned.

By default, when the nearest gridpoint value is a missing value or the location is out of the grid area, nil is returned in the case of a single coordinate, or vector_missing_value in the case of a vector. If an extra parameter 'valid' is added to the function call, then of the surrounding points, the nearest valid one is returned; nil will still be returned if all the surrounding points are missing.

The following example illustrates how to use the function.

fs = read (strGribFile)
listdef = nearest_gridpoint_info(fs, 51.46, -1.33)
loop ngp in listdef
     print ("Value : ", ngp.value)
     print ("Latitude : ", ngp.latitude)
     print ("Longitude : ", ngp.longitude)
end loop

...

Returns a copy of the input fieldset (first argument) with all of its missing values replaced with the number specified by the second argument. See also bitmap .

fieldset percentile(...)

Computes the specified percentiles for a given fieldset. This is a Metview icon function, for detailed documentation please see Percentile.

...

This function creates fields of pressure from the logarithm of the surface pressure (lnsp) and a list of model levels. Note that this function only works with gridded fields and assumes that the parameter for lnsp is 152. A newer, more flexible version of this function exists - see unipressure () .

• • The first argument is always a fieldset containing an lnsp field. If no other parameter is given, the list of levels will range from 1 to (number of vertical coordinates/2)-1 as coded in the GRIB header of the lnsp parameter.
  • The second argument specifies the levels at which the output fields must be generated. To generate a single level, pass a number. For more than one level, either pass a list of levels or a fieldset. If a fieldset is passed as the second parameter, the level information is extracted from each field of the fieldset.

Missing values in the lnsp field are retained in the output fieldset.

...

Similar to mask , except that a round mask is computed with a given radius around a geographical centre point. These can be given by either:

• • three numbers : latitude, longitude (in degrees), radius (in meters)
  • a list containing the above three numbers

The name of this function is derived from the fact that it creates a "round mask" or a "radius mask".

...

Computes the root mean square of a fieldset. A missing value in any field will result in a missing value in the corresponding place in the output fieldset. With n fields in the input fieldset, if x_i ^k is the i^th value of the k^th input field and y_i is the i^th value of the resulting field, the formula can be written :

Image Removed

Note that the following lines are equivalent :

y = rms(x)
y = sqrt(mean(x*x))

...

Computes the second zonal (from West to East) partial derivative of each field in the fieldset. The computations for a field f are based on the following formula:

\frac {\partial^2 f}{\partial x^2} = \frac{1}{R^2 \ cos^2\phi}\frac{\partial^2 f}{\partial \lambda^2} 

where

• • R is the radius of the Earth
  • φ is the latitude
  • λ is the longitude.

The derivatives are computed with a second order finite-difference approximation. The resulting fields contain missing values on the poles. Please note that this function is only implemented for regular latitude-longitude grids.

...

Computes the second meridional (from South to North) partial derivative of each field in the fieldset. The computations for a field f are based on the following formula:

\frac {\partial^2 f}{\partial y^2} = \frac{1}{R^2}\frac{\partial^2 f}{\partial \phi^2} 

where

• • R is the radius of the Earth
  • φ is the latitude.

The derivatives are computed with a second order finite-difference approximation. The resulting fields contain missing values on the poles. Please note that this function is only implemented for regular latitude-longitude grids.

...

Creates a new fieldset with all the fields' values replaced by those supplied. If supplied as a single vector, the values are set in all fields; if a list of vectors is supplied then there must be the same number of vectors as there are fields in the fieldset. The default behaviour is to produce an error if the input fieldset and vector have different numbers of values. If, however, a third parameter (set to the string 'resize') is passed to the function, the resulting fieldset will instead be resized to have the same number of values as the input vector - this can be useful when creating a new fieldset from a template. Missing values in the vector(s) are retained as missing values in the fieldset.

...

For each field in the input fieldset, this function creates a field where each grid point has the value of the sine of its latitude. For example, the following macro adds the coriolis parameter to each grid point of a field :

# Computes absolute vorticity from vorticity
omega = 2 * pi / 86400
coriolis = 2 * omega * sinlat(vort)
absvort = vort + coriolis

...

This function accepts a fieldset as input and returns it sorted according to keys and rules specified in the other input arguments.

Specified with only the fieldset as its single argument, sort() sorts in ascending order the fieldset according to the following MARS keys: date, time, step, number (ensemble member), levelist and param (integer ID).

The second argument allows you to modify the precedence of the sorting keys - e.g. if the second argument is a string "param", then the sorting is done according to the param key. If the second argument is a list you specify a list of sorting keys - e.g. ["param", "date"] sorts on parameter and then date.

The third argument specifies a sorting direction. This can be a string (">" or "<") or a list ([">", "<", ">",...]). If it is a string, the sorting direction it specifies applies to all sorting keys specified in the second argument. If it is a list, then the second argument must also be a list with the same number of elements - the sorting directions apply to each sorting key specified.

...

Computes the standard deviation of a fieldset. A missing value in any field will result in a missing value in the corresponding place in the output fieldset. With n fields in the input fieldset, if xik is the ith value of the kth input field and yi is the ith value of the resulting field, the formula can be written :

Image Removed

Note that the following lines are equivalent :

y = stdev(x)
y = sqrt(mean(x*x)-mean(x)*mean(x))
y = sqrt(var(x))

...

Computes the standard deviation over a weighted area. The area, if specified, is a list of numbers representing North, West, South, East. If the area is not specified, the whole field will be used in the calculation. The result is a number for a single field, or a list for a multi-field fieldset.

...

Computes the sum of a fieldset. The result is a single-field fieldset. A missing value in any field will result in a missing value in the corresponding place in the output fieldset. With n fields in the input fieldset, if xik is the ith value of the kth input field and yi is the ith value of the resulting field, the formula can be written :

Image Removed

...

Returns the indexes of the four gridpoints surrounding the given location, ordered by increasing distance from the target point. If a list is given, it must contain two numbers - latitude and longitude. If two numbers are given, the first is the latitude, the second the longitude. The field must be a gridded field. If the fieldset has only one field, a single vector of indexes is returned; otherwise a list of vectors is returned. In the case where the field is a reduced Gaussian grid and the input location is at the North or South pole, beyond the most extreme row of points, there will be a 'circle' of surrounding points, and all of these indexes are returned.

For batch processing of multiple locations, two vectors can be given, the first is a vector of latitudes, the second the longitudes; this can be much more efficient than multiple calls with a single location each. If the fieldset has only one field, a single vector is returned; otherwise a list of vectors is returned.

By default, if any of the surrounding points are missing, the function will return nil. To prevent this, and to return all the points regardless, add the option 'all' as the last parameter of the function call.

...

For each field in the input fieldset, this function creates a field where each grid point has the value of the tangent of its latitude. The resulting fields contain missing values on the poles.

...

This function creates fields of thickness from the logarithm of the surface pressure (lnsp ) and a list of model levels. Note that this function only works with lat/long grids and assumes that the parameter for lnsp is 152. A newer, more flexible version of this function exists - see unithickness () .

• • The first argument is always a fieldset containing an lnsp field. If no other parameter is given, the list of levels will range from 1 to (number of vertical coordinates/2)-1 as coded in the GRIB header of the lnsp .
  • The second argument specifies the levels at which the output fields must be generated. To generate a single level, pass a number. For more than one level, either pass a list of levels or a fieldset. If a fieldset is passed as the second parameter, the level information is extracted from each field of the fieldset.

Missing values in the lnsp field are retained in the output fieldset.

...

This function creates fields of pressure from the logarithm of the surface pressure (lnsp) and a list of model levels. Unlike pressure() , this function works with all grid types known to Metview (not just lat/long); it also allows the user to override the parameter number for lnsp (default 152).

• • The first argument is always a fieldset containing an lnsp field. If no other parameter is given, then pressure is computed for all model levels that are described in the GRIB header of fieldset .
  • If number is given (always the last parameter) it is the lnsp parameter code (default is 152).
  • list should contain model levels for which pressure is to be computed. Note that also for a single model level one has to use a list (this is a signature difference compared to the old function pressure() ).
  • If fieldset is given as the second parameter then pressure is computed for those model levels found in the second fieldset.

Missing values in the lnsp field are retained in the output fieldset.

...

This function creates fields of thickness from the logarithm of the surface pressure (lnsp) and a list of model levels. Unlike thickness() , this function works with all grid types known to Metview (not just lat/long); it also allows the user to override the parameter number for lnsp (default 152).

• • The first argument is always a fieldset containing an lnsp field. If no other parameter is given, then thickness is computed for all model levels that are described in the GRIB header of fieldset .
  • If number is given (always the last parameter) it is the lnsp parameter code (default is 152).
  • list should contain model levels for which thickness is to be computed. Note that also for a single model level one has to use a list (this is a signature difference compared to the old function thickness() ).
  • If fieldset is given as the second parameter then thickness is computed for those model levels found in the second fieldset.

Missing values in the lnsp field are retained in the output fieldset.

...

This function performs a vertical integration for data on. pressure levels or ECMWF (hybrid) model levels. This function is similar to vertint() , except that univertint() also works with sparse levels whereas vertint() is restricted to continuous levels.

• Pressure levels: the function has to be called with one argument only (the input fieldset).
• Model levels: If just one fieldset is given, it must also contain the logarithm of the surface pressure (lnsp). If the function is called with two fieldsets, the first one is a fieldset containing an lnsp field, the second one is the multi-level fieldset. A number may optionally be given as the last parameter in order to specify the lnsp code used in the fieldset that contains the lnsp data; its default value is 152. If a list with two elements [top, bottom] is given as the third parameter, then the integration is performed between (and including) these layers.

A missing value in any field will result in a missing value in the corresponding place in the output fieldset.

The function computes:

\int_{bottom}^{top} f \frac{dp}{g}

where

• • f is the fieldset
  • p is the pressure
  • g is the acceleration of gravity (9.80665 m/s²).

...

Returns the valid dates (including the time components) of the given fields. If the fieldset has only one field, a date is returned; otherwise a list of dates is returned.

...

This function returns the grid point values as a vector. If the fieldset contains more than one field it returns a list of vectors. Each of these vectors contains as many elements as grid points in each field. Missing values are included in the results as vector_missing_value.

...

Computes the variance of a fieldset. A missing value in any field will result in a missing value in the corresponding place in the output fieldset. With n fields in the input fieldset, if xik is the ith value of the kth input field and yi is the ith value of the resulting field, the formula can be written :
Image Removed
Note that the following lines are equivalent :

y = var(x)
y = mean(x*x)-mean(x)*mean(x)

...

Computes the variance over a weighted area. The area, if specified, is a list of numbers representing North, West, South, East. If the area is not specified, the whole field will be used in the calculation. The result is a number for a single field, or a list for a multi-field fieldset.

...

fieldset vertint ( fieldset,fieldset )

This function performs a vertical integration of the input fieldset, which must contain a range of model levels for the same parameter. A missing value in any field will result in a missing value in the corresponding place in the output fieldset. If the function is called with the fieldset as its single argument, it must also contain the logarithm of the surface pressure (lnsp ). If the function is called with two parameters, the first one is a fieldset containing an lnsp field, the second one is the multi-level fieldset.

The function computes :

\int_{bottom}^{top} f \frac{dp}{g}

where

• • f is the fieldset
  • p is the pressure
  • g is the acceleration of gravity (9.80665 m/s²).

The following example computes the total amount of liquid water in the atmosphere by integrating the cloud liquid water content (clwc ) over all levels of the model

...

Note that a similar function, interpolate(), also exists.

geopoints nearest_gridpoint ( fieldset,geopoints )

Generates a set of geopoints from a field. The first field of the input fieldset is used. The result is a set of geopoints whose values are those of the nearest gridpoints in the field to the geopoints given as a second parameter. Where it is not possible to generate a sensible value due to lack of valid data in the fieldset, the internal geopoints missing value is used (this value can be checked for with the built-in variable geo_missing_value or removed with the function remove_missing_values). Note that a similar function, interpolate() , also exists.

Returns the value and location of the nearest point to a given location in each field of a fieldset. If a list is given, it must contain two numbers - latitude and longitude. If two numbers are given, the first is the latitude, the second the longitude. The field must be a lat-long field. The return value is a list of definitions, one for each field, and each containing the following members: value , latitude , longitude . Where it is not possible to generate a sensible value due to lack of valid data in the fieldset, nil is returned. If an extra parameter 'valid' is added to the function call, then of the surrounding points, the nearest valid one is returned; nil will still be returned if all the surrounding points are missing.

The following example illustrates how to use the function.

fs = read (strGribFile)
listdef = nearest_gridpoint_info(fs, 51.46, -1.33)
loop ngp in listdef
     print ("Value : ", ngp.value)
     print ("Latitude : ", ngp.latitude)
     print ("Longitude : ", ngp.longitude)
end loop

Returns a copy of the input fieldset (first argument) with all of its missing values replaced with the number specified by the second argument. See also bitmap .

fieldset percentile(...)

Computes the specified percentiles for a given fieldset. This is a Metview icon function, for detailed documentation please see Percentile.

This function creates fields of pressure from the logarithm of the surface pressure (lnsp) and a list of model levels. Note that this function only works with gridded fields and assumes that the parameter for lnsp is 152. A newer, more flexible version of this function exists - see unipressure () .

• • The first argument is always a fieldset containing an lnsp field. If no other parameter is given, the list of levels will range from 1 to (number of vertical coordinates/2)-1 as coded in the GRIB header of the lnsp parameter.
  • The second argument specifies the levels at which the output fields must be generated. To generate a single level, pass a number. For more than one level, either pass a list of levels or a fieldset. If a fieldset is passed as the second parameter, the level information is extracted from each field of the fieldset.

Missing values in the lnsp field are retained in the output fieldset.

Similar to mask , except that a round mask is computed with a given radius around a geographical centre point. These can be given by either:

• • three numbers : latitude, longitude (in degrees), radius (in meters)
  • a list containing the above three numbers

The name of this function is derived from the fact that it creates a "round mask" or a "radius mask".

Computes the root mean square of a fieldset. A missing value in any field will result in a missing value in the corresponding place in the output fieldset. With n fields in the input fieldset, if x_i ^k is the i^th value of the k^th input field and y_i is the i^th value of the resulting field, the formula can be written :

Image Added

Note that the following lines are equivalent :

y = rms(x)
y = sqrt(mean(x*x))

Computes the second zonal (from West to East) partial derivative of each field in the fieldset. The computations for a field f are based on the following formula:

\frac {\partial^2 f}{\partial x^2} = \frac{1}{R^2 \ cos^2\phi}\frac{\partial^2 f}{\partial \lambda^2} 

where

• • R is the radius of the Earth
  • φ is the latitude
  • λ is the longitude.

The derivatives are computed with a second order finite-difference approximation. The resulting fields contain missing values on the poles. Please note that this function is only implemented for regular latitude-longitude grids.

Computes the second meridional (from South to North) partial derivative of each field in the fieldset. The computations for a field f are based on the following formula:

\frac {\partial^2 f}{\partial y^2} = \frac{1}{R^2}\frac{\partial^2 f}{\partial \phi^2} 

where

• • R is the radius of the Earth
  • φ is the latitude.

The derivatives are computed with a second order finite-difference approximation. The resulting fields contain missing values on the poles. Please note that this function is only implemented for regular latitude-longitude grids.

Creates a new fieldset with all the fields' values replaced by those supplied. If supplied as a single vector, the values are set in all fields; if a list of vectors is supplied then there must be the same number of vectors as there are fields in the fieldset. The default behaviour is to produce an error if the input fieldset and vector have different numbers of values. If, however, a third parameter (set to the string 'resize') is passed to the function, the resulting fieldset will instead be resized to have the same number of values as the input vector - this can be useful when creating a new fieldset from a template. Missing values in the vector(s) are retained as missing values in the fieldset.

New in Metview version 5.13.0.

Computes the shear deformation of 2-dimensional vector fields. The computations for a vector field f=(fx,fy) are based on the following formula:

d(f) = \frac{1}{R \ cos\phi}\frac{\partial f_y}{\partial \lambda} + \frac{1}{R}\frac{\partial f_x}{\partial \phi} + \frac{f_x}{R}tan\phi

where

• • R is the radius of the Earth (m)
  • φ is the latitude
  • λ is the longitude.

The derivatives are computed with a second order finite-difference approximation. The resulting fields contain missing values on the poles.  Please note that this function is only implemented for regular latitude-longitude grids.

For each field in the input fieldset, this function creates a field where each grid point has the value of the sine of its latitude. For example, the following macro adds the coriolis parameter to each grid point of a field :

# Computes absolute vorticity from vorticity
omega = 2 * pi / 86400
coriolis = 2 * omega * sinlat(vort)
absvort = vort + coriolis

New in Metview version 5.14.0.

Computes the solar zenith angle for each gridpoint by using the following positional arguments:

• • fs: input fieldset
  • to_cosine: (optional) when this argument is specified as set to "to_cosine" the cosine of the solar zenith angle is returned

The result is the solar zenith angle in degrees (unless "to_cosine" is specified when the cosine of the solar zenith angle is returned). The computations are based on the following formula:

cos\theta_{s} = sin\phi\, sin\delta + cos\phi\, cos\delta\, cosh 

where

• • θ is the solar zenith angle
  • φ is the latitude
  • δ is the declination of the Sunδ
  • h is the hour angle in local solar time

The declination of the Sun is computed as:

\delta = - arcsin\left(0.39779 cos(0.98565\unicode{xB0} (N+10) + 1.914\unicode{xB0} sin(0.98565\unicode{xB0} (N-2))\right) 

where:

• • N is the day of the year beginning with N=0 at midnight Universal Time (UT) as January 1. It is a floating point number allowing for fractional days.

A missing value in any field in fs will result in a missing value in the corresponding grid point in the output fieldset.

The dates and times used in the computations are based on the "validityDate" and "validityTime" ecCodes keys. If these are not available for a given field the result will contain missing values for all the gridpoints for that field.

When "to_cosine" is not specified and the GRIB edition of the input field is 2 the ecCodes paramId in the output field is set to 260225 (shortName="solza"). For GRIB edition 1 this parameter is not defined.

When "to_cosine" is specified the ecCodes paramId in the output is set to 214001 (shortName="uvcossza").

This function accepts a fieldset as input and returns it sorted according to keys and rules specified in the other input arguments.

Specified with only the fieldset as its single argument, sort() sorts in ascending order the fieldset according to the following MARS keys: date, time, step, number (ensemble member), levelist and param (integer ID).

The second argument allows you to modify the precedence of the sorting keys - e.g. if the second argument is a string "param", then the sorting is done according to the param key. If the second argument is a list you specify a list of sorting keys - e.g. ["param", "date"] sorts on parameter and then date.

The third argument specifies a sorting direction. This can be a string (">" or "<") or a list ([">", "<", ">",...]). If it is a string, the sorting direction it specifies applies to all sorting keys specified in the second argument. If it is a list, then the second argument must also be a list with the same number of elements - the sorting directions apply to each sorting key specified.

New in Metview version 5.14.0.

Computes the wind speed from the u and v wind components.

The resulting values are speed values in the same units as the input fields. A missing value in either u or v will result in a missing value in the corresponding place in the output fieldset. The ecCodes paramId in the output is set as follows:

• • 10 (atmospheric wind speed)
  • 207 (10m wind speed)
  • 228249 (100m wind speed)
  • 228241 (200m wind speed)

In any other cases the ecCodes paramId is set to 10.

Computes the standard deviation of a fieldset. A missing value in any field will result in a missing value in the corresponding place in the output fieldset. With n fields in the input fieldset, if xik is the ith value of the kth input field and yi is the ith value of the resulting field, the formula can be written :

Image Added

Note that the following lines are equivalent :

y = stdev(x)
y = sqrt(mean(x*x)-mean(x)*mean(x))
y = sqrt(var(x))

Computes the standard deviation over a weighted area. The area, if specified, is a list of numbers representing North, West, South, East. If the area is not specified, the whole field will be used in the calculation. The result is a number for a single field, or a list for a multi-field fieldset.

New in Metview version 5.13.0.

Computes the stretch deformation of 2-dimensional vector fields. The computations for a vector field f=(fx,fy) are based on the following formula:

d(f) = \frac{1}{R \ cos\phi}\frac{\partial f_x}{\partial \lambda} - \frac{1}{R}\frac{\partial f_y}{\partial \phi} - \frac{f_y}{R}tan\phi

where

• • R is the radius of the Earth (m)
  • φ is the latitude
  • λ is the longitude.

The derivatives are computed with a second order finite-difference approximation. The resulting fields contain missing values on the poles.  Please note that this function is only implemented for regular latitude-longitude grids.

Computes the sum of a fieldset. The result is a single-field fieldset. A missing value in any field will result in a missing value in the corresponding place in the output fieldset. With n fields in the input fieldset, if xik is the ith value of the kth input field and yi is the ith value of the resulting field, the formula can be written :

Image Added

Returns the indexes of the four gridpoints surrounding the given location, ordered by increasing distance from the target point. If a list is given, it must contain two numbers - latitude and longitude. If two numbers are given, the first is the latitude, the second the longitude. The field must be a gridded field. If the fieldset has only one field, a single vector of indexes is returned; otherwise a list of vectors is returned. In the case where the field is a reduced Gaussian grid and the input location is at the North or South pole, beyond the most extreme row of points, there will be a 'circle' of surrounding points, and all of these indexes are returned.

For batch processing of multiple locations, two vectors can be given, the first is a vector of latitudes, the second the longitudes; this can be much more efficient than multiple calls with a single location each. If the fieldset has only one field, a single vector is returned; otherwise a list of vectors is returned.

By default, if any of the surrounding points are missing, the function will return nil. To prevent this, and to return all the points regardless, add the option 'all' as the last parameter of the function call.

For each field in the input fieldset, this function creates a field where each grid point has the value of the tangent of its latitude. The resulting fields contain missing values on the poles.

This function creates fields of thickness from the logarithm of the surface pressure (lnsp ) and a list of model levels. Note that this function only works with lat/long grids and assumes that the parameter for lnsp is 152. A newer, more flexible version of this function exists - see unithickness () .

• • The first argument is always a fieldset containing an lnsp field. If no other parameter is given, the list of levels will range from 1 to (number of vertical coordinates/2)-1 as coded in the GRIB header of the lnsp .
  • The second argument specifies the levels at which the output fields must be generated. To generate a single level, pass a number. For more than one level, either pass a list of levels or a fieldset. If a fieldset is passed as the second parameter, the level information is extracted from each field of the fieldset.

Missing values in the lnsp field are retained in the output fieldset.

This function creates fields of pressure from the logarithm of the surface pressure (lnsp) and a list of model levels. Unlike pressure() , this function works with all grid types known to Metview (not just lat/long); it also allows the user to override the parameter number for lnsp (default 152).

• • The first argument is always a fieldset containing an lnsp field. If no other parameter is given, then pressure is computed for all model levels that are described in the GRIB header of fieldset .
  • If number is given (always the last parameter) it is the lnsp parameter code (default is 152).
  • list should contain model levels for which pressure is to be computed. Note that also for a single model level one has to use a list (this is a signature difference compared to the old function pressure() ).
  • If fieldset is given as the second parameter then pressure is computed for those model levels found in the second fieldset.

Missing values in the lnsp field are retained in the output fieldset.

This function creates fields of thickness from the logarithm of the surface pressure (lnsp) and a list of model levels. Unlike thickness() , this function works with all grid types known to Metview (not just lat/long); it also allows the user to override the parameter number for lnsp (default 152).

• • The first argument is always a fieldset containing an lnsp field. If no other parameter is given, then thickness is computed for all model levels that are described in the GRIB header of fieldset .
  • If number is given (always the last parameter) it is the lnsp parameter code (default is 152).
  • list should contain model levels for which thickness is to be computed. Note that also for a single model level one has to use a list (this is a signature difference compared to the old function thickness() ).
  • If fieldset is given as the second parameter then thickness is computed for those model levels found in the second fieldset.

Missing values in the lnsp field are retained in the output fieldset.

This function performs a vertical integration for data on pressure levels or ECMWF (hybrid) model levels. This function is similar to vertint() , except that univertint() also works with sparse levels whereas vertint() is restricted to continuous levels.

• Pressure levels: the function has to be called with one argument only (the input fieldset).
• Model levels: If just one fieldset is given, it must also contain the logarithm of the surface pressure (lnsp). If the function is called with two fieldsets, the first one is a fieldset containing an lnsp field, the second one is the multi-level fieldset. A number may optionally be given as the last parameter in order to specify the lnsp code used in the fieldset that contains the lnsp data; its default value is 152. If a list with two elements [top, bottom] is given as the third parameter, then the integration is performed between (and including) these layers.

A missing value in any field will result in a missing value in the corresponding place in the output fieldset.

The function computes:

\int_{p_{top}}^{p_{bottom}} f \frac{dp}{g}

where

• • f is the fieldset
  • p is the pressure
  • g is the acceleration of gravity (9.80665 m/s²).

Returns the valid dates (including the time components) of the given fields. If the fieldset has only one field, a date is returned; otherwise a list of dates is returned.

This function returns the grid point values as a vector. If the fieldset contains more than one field it returns a list of vectors. Each of these vectors contains as many elements as grid points in each field. Missing values are included in the results as vector_missing_value.

# x is a fieldset of n fields
xgrid = values(x)
field1_values = xgrid[1]
gridpoint1 = field1_values[1]

# or equivalently
gridpoint1 = xgrid[1][1]

Computes the variance of a fieldset. A missing value in any field will result in a missing value in the corresponding place in the output fieldset. With n fields in the input fieldset, if xik is the ith value of the kth input field and yi is the ith value of the resulting field, the formula can be written :
Image Added
Note that the following lines are equivalent :

y = var(x)
y = mean(x*x)-mean(x)*mean(x)

Computes the variance over a weighted area. The area, if specified, is a list of numbers representing North, West, South, East. If the area is not specified, the whole field will be used in the calculation. The result is a number for a single field, or a list for a multi-field fieldset.

fieldset vertint ( fieldset,fieldset )

This function performs a vertical integration of the input fieldset, which must contain a range of model levels for the same parameter. A missing value in any field will result in a missing value in the corresponding place in the output fieldset. If the function is called with the fieldset as its single argument, it must also contain the logarithm of the surface pressure (lnsp ). If the function is called with two parameters, the first one is a fieldset containing an lnsp field, the second one is the multi-level fieldset.

The function computes:

\int_{p_{top}}^{p_{bottom}} f \frac{dp}{g}

where

• • f is the fieldset
  • p is the pressure
  • g is the acceleration of gravity (9.80665 m/s²).

The following example computes the total amount of liquid water in the atmosphere by integrating the cloud liquid water content (clwc ) over all levels of the model

# Retrieve clwc
clwc = retrieve(
   levtype : "ml",
   levelist : [1,"to",31],
   param : "clwc",
   date : -1,
   grid : [2.5,2.5]
)

# Retrieve lnsp
lnsp = retrieve(
   levtype : "ml",
   levelist : 1,
   param : "lnsp",
   date : -1,
   grid : [2.5,2.5]
)

# Integrate the field
x = vertint(lnsp,clwc)
plot(x)

Computes the vertical component of the curl differential operator for 2-dimensional vector fields. For wind fields (i.e. when the input fieldsets are u and v wind components) this computes the relative vorticity (ζ). The computations for a vector field f=(f_x ,f_y ) are based on the following formula :

 \zeta =\frac{1}{R \ cos\phi}\frac{\partial f_y}{\partial \lambda} - \frac{1}{R}\frac{\partial f_x}{\partial \phi} + \frac{f_x}{R}tan\phi

where

• • R is the radius of the Earth
  • φ is the latitude
  • λ is the longitude.

The derivatives are computed with a second order finite-difference approximation. The resulting fields contain missing values on the poles. If the input fields are horizontal wind components the GRIB paramId of the resulting field is set to 138 (relative vorticity). Please note that this function is only implemented for regular latitude-longitude grids.

number w_from_omega(omega: number,  t: number,  p: number)

vector w_from_omega(omega: vector,  t: vector,  p: vector)

fieldset w_from_omega(omega: fieldset, t: fieldset)

fieldset w_from_omega(omega: fieldset,  t: fieldset,  p: fieldset)

Computes the hydrostatic vertical velocity from  pressure velocity (omega) for a given temperature and pressure, where

• • omega: pressure velocity (Pa/s)
  • t: temperature (K)
  • p: pressure (Pa)

The result is the vertical velocity in m/s units. On error  nil is returned. The following rules are applied when omega is a fieldset:

• • if omega is a pressure level fieldset no pressure argument is needed
  • if omega is defined on ECMWF model levels (hybrid/eta) p is either a single LNSP (logarithm of surface pressure, identified by paramId=152) field or a fieldset defining the pressure on the levels of omega
  • for other level types p is a fieldset defining the pressure on the levels of omega

The computation is based on the following hydrostatic formula:

w = - \frac{\omega T R_{d}}{p g}

where

• Rd is the specific gas constant for dry air (287.058 J/(kg K)).
• g is the gravitational acceleration (9.81 m/s²)

This functions was introduced in version 5.10.0.

list of numbers xy_from_polar(magnitude: number,  dir: number)

list of vectors xy_from_polar(magnitude: vector,  dir: vector)

fieldset xy_from_polar(magnitue: fieldset, dir: fieldset)

Converts vector data from meteorological polar representation to xy representation. In polar representation the data is specified by two components:

• magnitude: represents the speed/magnitude
• dir: represents the direction of the vector in degrees. Angles measured from South in clockwise direction.

In the target xy representation the x axis points East while the y axis points North.

The type of the result depends on the type of the input data

• if the input is number the result is a list of two numbers
• if the input is vector the result is a list of two vectors, the first vector contains the x components while the second vector the y components
• if the input is fieldset the result is a fieldset where an x component field is immediately followed by the corresponding y component field.

This function was introduced in version 5.10.0

...

...

Computes the vertical component of the curl differential operator for 2-dimensional vector fields. For wind fields (i.e. when the input fieldsets are u and v wind components) this computes the relative vorticity (ζ). The computations for a vector field f=(f_x ,f_y ) are based on the following formula :

 \zeta =\frac{1}{R \ cos\phi}\frac{\partial f_y}{\partial \lambda} - \frac{1}{R}\frac{\partial f_x}{\partial \phi} + \frac{f_x}{R}tan\phi

where

• • R is the radius of the Earth
  • φ is the latitude
  • λ is the longitude.

...

.