compute is not available when accessing MARS via the Web API.
Mathematical computations can be done on retrieved GRIB fields using the compute verb. It uses keyword fieldset to store the result of the computation and formula which specifies the actual computation. The following conventions are used in this documentation:
X | upper case characters denote fieldsets |
Xk | is the k-th field of fieldset X |
Xijk | is the point (column i, line j) of field Xk |
a | lower case characters denote real scalars |
'...' | is a character string |
f(...) | is a function call |
formula = "..." | is the MARS formula, what is actually computed |
See also the fieldset syntax on how to access a subset of a fieldset.
Here is an example, where X and Y are the two input fieldsets and Z is the output fieldset:
compute, fieldset = Z, formula = "X/2+log(Y)*X" |
Any mathematical formula can be reduced to a series of elementary computations involving only one or two values. So the formula in the example can be reduced to:
[1] compute X/2 | [2] compute log(Y) |
[3] compute [2]*X | [4] compute [1]+[3] |
For each of the fields Xk and Yk of X and Y, a field Zk is computed using the formula and added to the fieldset Z.
Any operation involving a missing value will results in a missing value.
All compute operators available with the MARS compute verb are listed in the table below. The sign '-' can also be used to denote negation, in which case it is actually a function, not an operator. All the boolean operators return 1 if the test is true, 0 otherwise. MARS follows the mathematical order of precedence: 1. function calls, 2. power, 3. * and /, 4. + and -, 5. boolean operators. By using brackets, the order of computation can be changed.
MARS compute operators | |||
---|---|---|---|
Arithmetic operators | Boolean operators | ||
+ | addition | = | equal to |
- | subtraction | > | greater than |
* | multiplication | < | less than |
/ | division | >= | greater or equal to |
^ | power | <= | less than or equal to |
<> | not equal to |
Note: No check is done on the type of data combined in the formulae. It is the user's responsibility to make sure the computations have a physical meaning.
The GRIB header of the resulting fields is copied from the first field encountered. The resulting fields' GRIB headers are marked as being created with the compute facility by setting the model identification number to 255. Users can preserve the model identification number by setting the environment variable MARS_COMPUTE_FLAG to 0 before calling the mars command (see table 3).
The GRIB headers of resulting fields can be further modified using the GRIB tools, e.g. grib_set to read all fields from input file "speed" and write them as parameter "wind speed" to file out.grib:
grib_set -s paramId=10 speed out.grib
or
grib_set -s shortName=ws speed out.grib
Mixing fields and scalars
- If the two operands are two fields, the operation is carried out between each pair of corresponding grid points in the two fields. The result is a field.
- If one operand is a scalar and the other a field, the operation is carried out by applying the scalar to each grid point. The result is a field.
- If both operands are scalars, the result is a scalar.
Calling functions
The same logic applies to mathematical functions, which are listed in the table below. If the argument of a function is a field, the result is a field where each grid point is the result of the function at the corresponding grid point in the input field. If the argument is a scalar the output is also a scalar. Examples:
compute, fieldset = X, formula = "A * B" |
creates a fieldset X where Xijk = Aijk * Bijk. Note that * is scalar multiplication and not matrix multiplication.
compute, fieldset = X, formula = "A * 2" |
creates a fieldset X where Xijk = Aijk * 2.
compute, fieldset = X, formula = "(A>10) * B + (A<=10)*C" |
Creates a fieldset X, where
Xijk = Bijk | if | Aijk > 10 |
Xijk = Cijk | if | Aijk < = 10 |
MARS mathematical functions | ||||
---|---|---|---|---|
Name | No. of param. | Param. type | Result type | Description |
max | any | any | fieldset if one of the parameters is a fieldset, scalar otherwise | Maximum |
min | any | any | fieldset if one of the parameters is a fieldset, scalar otherwise | Minimum |
merge | any | fieldset | fieldset | Merges several fieldsets into a unique fieldset |
maxvalue | 1 | fieldset | scalar | Maximum grid point value of all fields |
minvalue | 1 | fieldset | scalar | Minimum grid point value of all fields |
sgn | 1 | any | same as parameter | 0 if the value is 0, +1 if the value is positive, -1 if the value is negative. |
int | 1 | any | same as parameter | Integer part |
exp | 1 | any | same as parameter | Exponential |
log | 1 | any | same as parameter | Natural logarithm |
log10 | 1 | any | same as parameter | Base 10 logarithm |
sin | 1 | any | same as parameter | Sine, input in radians |
cos | 1 | any | same as parameter | Cosine, input in radians |
tan | 1 | any | same as parameter | Tangent, input in radians |
asin | 1 | any | same as parameter | Arc-sine, return value(s) in radians |
acos | 1 | any | same as parameter | Arc-cosine, return value(s) in radians |
atan | 1 | any | same as parameter | Arc-tangent, return value(s) in radians |
atan2 | 2 | any | fieldset if one of the parameters is a fieldset, scalar otherwise | Arc-tangent of y/x, using the signs of both parameters to determine the quadrant of the result (in radians) |
abs | 1 | any | same as parameter | Absolute value |
sqrt | 1 | any | same as parameter | Square root |
count | 1 | fieldset | scalar | Number of fields in the fieldset |
sum | 1 | fieldset | 1 field | Sum of all the fields in the fieldset. |
mean | 1 | fieldset | 1 field | Mean of all the fields in the fieldset. |
rms | 1 | fieldset | 1 field | Root mean square of all the fields in the fieldset. |
stdev | 1 | fieldset | 1 field | Standard deviation of all the fields in the fieldset. |
var | 1 | fieldset | 1 field | Variance of all the fields in the fieldset. |
covar | 2 | fieldset, fieldset | 1 field | Covariance between the 2 fieldsets. |
distribution | 2 | fieldset, scalar | fieldset | Sorts a set of n (indicated by scalar) fields according to the parameter values, e.g. field 1 will have the lowest values, field 2 will have the second lowest values ..., field n, which will have the largest values |
bitmap | 2 | fieldset, scalar | fieldset | Creates missing values for those grid points equal to scalar |
bitmap | 2 | fieldset, field | fieldset | Creates missing values for those missing values of the second field |
bitmap | 2 | fieldset, fieldset | fieldset | Creates missing values on the first parameter for those grid points with missing values of the second fieldset |
nobitmap | 2 | fieldset, scalar | fieldset | Removes bitmap on the fieldset by changing missing values by the scalar provided |
repeat | 2 | field, scalar | fieldset | Repeats field as many times as indicated by scalar |
Bitmaps (Missing values)
Some fields may contain a bitmap, i.e., there is no data (missing values) for certain grid points. Those missing values are not considered in any computation, but copied to the resulting field. Only valid data values are computed. Functions bitmap and nobitmap enable users to handle bitmaps on fieldsets.
Examples:
compute, fieldset = X, formula = "bitmap(A,0)" |
creates a fieldset X, where
Xijk = Aijk | if | Aijk ≠ 0 |
Xijk = missing value | if | Aijk = 0 |
compute, fieldset = X, formula = "bitmap(A,B)" |
If B is a field with a bitmap, this request creates a fieldset X, where Xk is a copy of Ak with the bitmap of B.
compute, fieldset = X, formula = "nobitmap(A,1)" |
creates a fieldset X, where
Xijk = Aijk | if | Aijk ≠ missing value |
Xijk = 1 | if | Aijk = missing value |
Compute example
retrieve, type = analysis, time = 12, date = -1, grid = 2.5/2.5, param = u, fieldset = u retrieve, param = v, fieldset = v compute, formula = "sqrt(u*u + v*v)", target = "speed"
The GRIB headers of the resulting fields can be further modified by using e.g. grib_set, see above.
4 Comments
Renate Wilcke
Hi there,
I would like to use compute to retrieve daily maximum temperature (and accumulated precipitation) and am wondering if there is a function like "daymax" ("daysum") or so. I tried max(tmaxfield), hoping it would be only applied to daily steps, but obviously got only one timestep per year as output. Is there already a neat function for me to use?
Otherwise I'll retrieve subdaily data and let cdo handle that.
/Renate
Dominique Lucas
Hello Renate,
These functions do not exist in mars. The function max will return the daily max if your input fieldset only contains the grib fields for one day. Like this, you should obtain the same result as with cdo. Note that the precipitation fields in MARS are accumulated. Depending on the forecast range you look at, you may need to de-accumulate the archived precipitation fields to obtain daily values. You can achieve this with MARS by carefully selecting the subsets of your precipitation fieldset, before doing the de-accumulation.
Good luck.
Dominique
Renate Wilcke
Thank you, Dominique! I understand and was suspecting something like this. I have to download 55 years of UERRA data. I think, I'll be faster in doing the daily statistics as post-processing, than in retrieving/computing single days and then joining the daily files to years again. I have to do quite some post-processing anyway. Only was thinking about the trafic (retrieving/downloading hourly or daily data makes a "small" difference ).
Dominique Lucas
For large extractions like yours, it's certainly worthwhile to spend some time choosing an efficient solution. You can reduce the traffic back to your side by doing your post-processing at ECMWF.