#Metview Macro

# **************************** LICENSE START ***********************************
#
# Copyright 2012 ECMWF. This software is distributed under the terms
# of the Apache License version 2.0. In applying this license, ECMWF does not
# waive the privileges and immunities granted to it by virtue of its status as
# an Intergovernmental Organization or submit itself to any jurisdiction.
#
# ***************************** LICENSE END ************************************

# Read the data file

tuv_data = read("/home/graphics/cgi/metview//macro_tutorial/macro_tut1/TUV_Data")


# Filter out the temperature forecasts and analyses

tfc = read(
        logstats : "",
        type     : "fc",
        param    : 130,
        date     : 20120228,
        step     : 48,
        data     : tuv_data
        )

tan = read(
        logstats : "",
        type     : "an",
        param    : 130,
        date     : 20120301,
        data     : tuv_data
        )


# Calcuate 'forecast minus analysis'

fa_diff = tfc - tan

pos = mcont(
        contour_line_thickness       : 2,
        contour_line_colour          : "red",
        contour_highlight            : "off",
        contour_level_selection_type : "level_list",
        contour_max_level            : 10,
        contour_min_level            : 0.5,
        contour_level_list           : [0.5,1,2,4,10]
        )

neg = mcont(
        contour_line_thickness       : 2,
        contour_highlight            : "off",
        contour_level_selection_type : "level_list",
        contour_max_level            : -0.5,
        contour_min_level            : -10,
        contour_level_list           : [-10,-4,-2,-1,-0.5]
        )

acoast = mcoast(
        map_coastline_resolution        : "low",
        map_coastline_thickness         : 3,
        map_coastline_land_shade        : "on",
        map_coastline_land_shade_colour : "grey",
        map_coastline_sea_shade         : "on",
        map_coastline_sea_shade_colour  : "RGB(0.9,0.95,1)",
        map_grid_longitude_increment    : 10
        )

ps_atlantic = geoview(
        map_projection      : "polar_stereographic",
        map_area_definition : "corners",
        area                : [30,-25,50,65],
        coastlines          : acoast
        )


page = plot_page(
        view  : ps_atlantic
        )

dw = plot_superpage(
        custom_width    : 29.7,
        custom_height   : 21.0,
        pages           : page
        )


# Plot the result

plot(dw, fa_diff, pos, neg)
