Versions Compared

Key

  • This line was added.
  • This line was removed.
  • Formatting was changed.

...

In this part we will estimate the risk of the wind gust being higher than a certain threshold. We will compute the probability of the wind gust exceeding 22 m/s (about 80 kkm/mh) and generate the plot shown below:

...

Now duplicate the ensemble mean Macro and edit it. Find the code line computing the mean and replace it with this code block:

...

  • the value is 1 in each gridpoint where the condition meets (i.e. the value is larger than the threshold)
  • the value is 0 in all other gridpoints.

The second line simply derives the probability as the mean of these fields. We multiply the result by 100 to scale it into the 0-100 range for an easier interpretation.

Once you finished your Macro drag it into the bottom right map and customise , visualise your Geographical View icon and drag the Macro into the  plot. Customise it with the prob_shade Contouring icon. Also use a custom Text Plotting icon to define the title. As for the probabilities, you should see that there is some probability of high wind speeds.

...

Code Block
plot(dw[i],title,f,wgust_shade)

Having done so, run your Macro (this will take a minute or so) and try to identify the ENS members predicting high wind speeds in our area.

...

so that our map could show a larger (North Atlantic) area.

Next, drop define the contouring used for the "spaghetti" by dropping the cont_spag Contouring icon into the Macro. A code like this should be generated for you:

...

In this mcont() we turned contour labels off to keep the plot uncluttered and defined only a single contour value (for 560 gpm). 

Continue with reading in the GRIB file of the ENS forecasts used for the "spaghetti":

Code Block
 g = read("ens_spag.grib")

The "spaghetti" will be generated by plotting each perturbed forecasts member as a separate layer into the same map. To achieve this goal we need to write a loop like this:

...

Code Block
f=read(data: g,
	type: "pf", 
	number: number
 )	

Next, plot it into your view with your contour settings:

Code Block
plot(your_view,f,cont) 

By default, if no title definition is specified, Metview adds a title line for each field in the plot. Since we are about to plot 50 fields into the same map this would result in 50 titles in the plot! To avoid having too many titles we use a custom Text Plotting icon:

...

We finish the loop by plotting the field with out our contour settings and title:

...

The stamp plot only shows the perturbed ENS members but there is still space left to display additional fields, as well. Try to add the control forecast (from ENS) and the operational forecast to it.

Hints:

  • plot the control forecast into the 51st map (dw[51]). The control forecast is stored in the same file as the perturbed forecast members: fc_ens.grib.

    To access it (via the read() function) you need to omit the number parameter and set type to "cf".

    Read it in with this code:

    Code Block
    f = read(data: g, type: "cf", step: step)
  • plot the operational forecast into the 52nd map (dw[52]). The operational forecast is stored in fc_oper.grib.

    For the read() command you just need to specify the step.

    Read it in with this code:

    Code Block
    f =read(source: "fc_oper.grib", step: step)
Info

While setting up these extra plots it is a good idea to temporarily comment out the loop processing the perturbed forecast members.

...

The spaghetti plot only shows the perturbed ENS members. Try to add the control forecast (from ENS) and the operational forecast to it as well. You should use different isoline colours for them.

Hints:

  • plot the control forecast with a thick green isoline. The control forecast is stored in the same file as the perturbed forecast members: spag_ens.grib.

    To access it (via the read() function) you need to omit the number parameter and set type to "cf".

    Read it in with this code:

    Code Block
    f = read(data: g, type: "cf")
  • plot the operational  forecast with thick red isoline. The operational forecast is stored in spag_oper.grib.

    You just need to read

    Read it in with

    the read() command. 

    this code:

    Code Block
    f = read("spag_oper.grib")
Info

While setting up these extra layers it is a good idea to temporarily comment out the whole processing the perturbed forecasts members.