Tephigram from ODB Example

#Metview Macro

#  **************************** LICENSE START ***********************************
# 
#  Copyright 2020 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 ************************************
# 

# station id
statid="40417"

# read db
db = read("temp.odb")

# run queries to get profile data
f_t = odb_filter(
	odb_query	:	"select obsvalue, vertco_reference_1 " & 
	                "  where varno=2 and statid='" & statid & "'",
	odb_data	:	db
	)
	
f_q = odb_filter(
	odb_query	:	"select obsvalue, vertco_reference_1 " & 
	                "  where varno=7 and statid='" & statid & "'",
	odb_data	:	db
	)

f_u = odb_filter(
	odb_query	:	"select obsvalue, vertco_reference_1 " & 
	                "  where varno=3 and statid='" & statid & "'",
	odb_data	:	db
	)

f_v = odb_filter(
	odb_query	:	"select obsvalue, vertco_reference_1 " & 
	                "  where varno=4 and statid='" & statid & "'",
	odb_data	:	db
	)

# define query for metadata 
f_meta=odb_filter(odb_query : "select lat@hdr as lat, lon@hdr as lon, " &
     "andate, antime where statid='" & statid & "'" ,
	odb_data	:	db
	)

# load temperature and humidity profiles into vectors
t = values(f_t,"obsvalue@body")
p_t = values(f_t,"vertco_reference_1@body")
q = values(f_q,"obsvalue@body")
p_q = values(f_q,"vertco_reference_1@body")

# compute dewpoint from specific humidity (q)
td = dewpoint_from_specific_humidity(q, p_q)

# scale data for visualisation
t = t - 273.16
td = td - 273.16
p_t = p_t/100
p_q = p_q/100

# define t curve
t_vis = xy_curve(t, p_t, "red", "solid", 8)

# define td curve
td_vis = xy_curve(td, p_q, "red", "dash", 6)

# load wind profiles into vectors
u = values(f_u,"obsvalue@body")
p = values(f_u,"vertco_reference_1@body")
v = values(f_v,"obsvalue@body")

# define wind x position in tephigram
wind_x_pos = u
for i=1 to count(wind_x_pos) do
    wind_x_pos[i] = 1010
end for
   
# define wind visualisation    
wp = mwind(
	wind_field_type  : "flags",
	wind_flag_colour : "evergreen"
    )

wind_vis = input_visualiser(
		input_plot_type          : "xy_vectors",
		input_x_values           : wind_x_pos,
		input_y_values           : p/100,
		input_x_component_values : u,
		input_y_component_values : v
		)


# define the title
dd=values(f_meta,"andate")[1]
tt=values(f_meta,"antime")[1]

title = mtext(text_line_1 : "TEMP sounding date=" & dd & " time=" & tt & " station=" & statid,       
         text_font_size: 0.4
        )

# define the view
view = cartesianview(
	map_projection	:	"tephigram",
	y_automatic_reverse	:	"on"
	)
 
# define the output plot file
setoutput(pdf_output(output_name : 'tephigram_odb')) 

# generate plot
plot(view, t_vis, td_vis, wind_vis, wp, title)

import metview as mv

#  **************************** LICENSE START ***********************************
#
#  Copyright 2020 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 ************************************
#

# station id
statid = "40417"

# read db
db = mv.read("temp.odb")

# run queries to get profile data
f_t = mv.odb_filter(
    odb_query="""select obsvalue, vertco_reference_1 
	                   where varno=2 and statid='{}'""".format(
        statid
    ),
    odb_data=db
)

f_q = mv.odb_filter(
    odb_query="""select obsvalue, vertco_reference_1  
	                  where varno=7 and statid='{}'""".format(
        statid
    ),
    odb_data=db
)

f_u = mv.odb_filter(
    odb_query="""select obsvalue, vertco_reference_1 
	                  where varno=3 and statid='{}'""".format(
        statid
    ),
    odb_data=db
)

f_v = mv.odb_filter(
    odb_query="""select obsvalue, vertco_reference_1 
	                  where varno=4 and statid='{}'""".format(
        statid
    ),
    odb_data=db
)

# define query for metadata
f_meta = mv.odb_filter(
    odb_query="""select lat@hdr as lat, lon@hdr as lon, 
     andate, antime where statid='{}'""".format(
        statid
    ),
    odb_data=db
)

# load temperature and humidity profiles into vectors
t = mv.values(f_t, "obsvalue@body")
p_t = mv.values(f_t, "vertco_reference_1@body")
q = mv.values(f_q, "obsvalue@body")
p_q = mv.values(f_q, "vertco_reference_1@body")

# compute dewpoint from specific humidity (q)
td = mv.dewpoint_from_specific_humidity(q, p_q)

# scale data for visualisation
t = t - 273.16
td = td - 273.16
p_t = p_t / 100
p_q = p_q / 100

# define t curve
t_vis = mv.xy_curve(t, p_t, "red", "solid", 8)

# define td curve
td_vis = mv.xy_curve(td, p_q, "red", "dash", 6)

# load wind profiles into vectors
u = mv.values(f_u, "obsvalue@body")
p = mv.values(f_u, "vertco_reference_1@body")
v = mv.values(f_v, "obsvalue@body")

# define wind x position in tephigram
wind_x_pos = [1010 for i in u]

# define wind visualisation
wp = mv.mwind(wind_field_type="flags", wind_flag_colour="evergreen")

wind_vis = mv.input_visualiser(
    input_plot_type="xy_vectors",
    input_x_values=wind_x_pos,
    input_y_values=p / 100,
    input_x_component_values=u,
    input_y_component_values=v
)


# define the title
dd = mv.values(f_meta, "andate")[1]
tt = mv.values(f_meta, "antime")[1]

title = mv.mtext(
    text_line_1="TEMP sounding date={} time={} station={}".format(dd, tt, statid),
    text_font_size=0.4,
)

# define the view
view = mv.cartesianview(map_projection="tephigram", y_automatic_reverse="on")

# define the output plot file
mv.setoutput(mv.pdf_output(output_name="tephigram_odb"))

# generate plot
mv.plot(view, t_vis, td_vis, wind_vis, wp, title)