Page History

...

Note that when passing numeric dates such as 20150105 to other modules, such as the MARS Retrieval module, these do not need to be converted into date variables. However, MARS treats Date and Time as separate parameters, so a date variable would need to be split into these components.

Date arithmetic

When dealing with dates, the number 1 represents one day. So the expression d1 + 1 gives a date one day later than day 1. To compute the difference, in days, between two dates, it's simply:

...

Looping through dates

Three examples (no need to type these in), to get a feel for it:

for d = 2015-01-01 to 2015-03-01 do
    print(d)
     ... # each step is 1 day
end for

for d = 2015-01-01 to 2015-03-01 by 2 do
    print(d)
     ... # each step is 2 days
end for

for d = 2015-01-01 to 2015-03-01 by hour(6) do
    print(d)  # each step is 6 hours
end for

Computing mean precipitation rate over an area

As an exercise to put all of this together, write a new macro to do the following:

• compute the 'period' precipitation from precip.grib - this is what we already did earlier, so it's done! Just copy the code to your new macro and change the result variable name to precip_period
• loop through the fields in the original fieldset and for each:
  • get the date and time of the forecast step
  • combine these into a Metview date variable
  • add it to a list
• use syntax similar to the line of code used to compute the 'period' precipitation to find the differences between the times of adjacent fields (ok, we know it's 3 hours, but in theory it could be anything)
• extract a sub-area from the computed precipitation data
• compute the mean value over that area for each field (use the average() function)
• scale up from metres (as the data are stored) to mm by multiplying by 1000
• convert this into a rate, mm per hour, using the time differences computed earlier

Here are some hints to help.

You can get the date and time of a field as numbers like this:

field_valid_date = grib_get_string(precip_3h[i], 'validityDate')
field_valid_time = grib_get_string(precip_3h[i], 'validityTime')

then combine those numbers into proper date variables.

A list is built up like this:

dates = nil
for i = .... do
    d = .....  # construct a date variable
    dates = dates & [d]
end for

To extract a sub-area, use the following syntax (derived from using the GRIB Filter icon's Area  parameter):

 f = read(data: f, area:[N,W,S,E])

The average() function returns a list of values, a mean for each field. You can directly multiply a list variable by a number to obtain a new list where each element has been multiplied.

The final calculation requires converting the time intervals into hours (because if the time difference between two steps is 7 hours, then the rate of precip per hour is the mean precip value divided by 7).

 XXXX remove once tested!

precip = read("precip.grib")

n = count(precip) # the number of fields in the fieldset
precip_period = precip[2, n] - precip[1, n-1]

for i = 1 to count(precip) do
    d_valid_date = grib_get_string(precip_3h[i], 'validityDate')
    t = grib_get_string(precip_3h[i], 'validityTime') step is 6 hours
end for