In the past pressure was recorded on beautiful drum barograph instruments (the meteorological equivalent of the drum seismograph).
Nowadays the instruments are digital, and people plot boring time series of pressure (time on the X axis, pressure on the Y axis). A few years ago I wrote a shell and GMT script to turn digital pressure records into old style barograph charts:
Barograph showing the pressure wave from the 15 January 2022 Hunga Tonga - Hunga Haʻapai volcanic eruption, as recorded at Honolulu International Airport. The pressure wave circumnavigated the world at least three times (in both directions), and was observed at many locations, including Honolulu.
The ksh script is below:
#!/usr/bin/env ksh
function barocoords {
#
# Convert date in format t,p,... to x,y,... on the barograph chart.
#
# Input data are provided on stdin and transformed data written to stdout.
# Ignore lines beginning with ">".
#
# Arguments
#
# $1: length of chart (cm) Default: 30
# $2: height of chart (cm) Default: 8
# $3: mimimum pressure (hPa) Default: 950
# $4: maximum pressure (hPa) Default: 1050
# $5: minimum time (s since 1970-01-01T00:00:00) Default: 7 days ago
# $6: maximum time (s since 1970-01-01T00:00:00) Default: Now
# $7: length of barograph arm (cm) Default: 15
awk -v plen="${1:-63}" -v phgt="${2:-15}" -v pmin="${3:-950}" \
-v pmax="${4:-1050}" -v tmin="${5:-$(($(gdate -u +%s)-604800))}" \
-v tmax="${6:-$(gdate -u +%s)}" -v arm="${7:-25}" '
/^[^>]/{
t = $1;
p = $2;
pmid = (pmin + pmax)/2;
y = phgt*(p - pmid)/(pmax - pmin);
d = arm - sqrt(arm^2 - y^2);
x = plen*(t - tmin)/(tmax - tmin) + d;
printf("%f\t%f\t", x, y);
for (ii=3; ii<=NF; ii++) printf("%s\t", $ii);
printf("\n");
}
/^>/{
printf("%s\n", $0);
}'
}
function decode_us {
#
# Decode 1-minute data from the US NWS.
#
gawk -F, '
/^[[:upper:]]{4}/{
year = substr($3,1,4);
mon = substr($3,6,2);
day = substr($3,9,2);
hour = substr($3,12,2);
min = substr($3,15,2);
utime = mktime(year" "mon" "day" "hour" "min" 00",1);
prs = $4*33.86388666; # Convert inches of Hg to hPa.
if (prs > 0){
printf("%d\t%f\n", utime, $4*33.86388666);
}
}'
}
function find_missing {
#
# Find missing time steps in (t, P) series, and insert GMT "gaps".
#
awk -v tmin=${1:-0} -v tmax=${2:-999999999999} '
BEGIN{
tprev = 0;
nsteps = 0;
ssteps = 0;
avgstep = 86400;
}
($1 >= tmin) && ($1 <= tmax){
if (($1 - avgstep) > tprev){
printf(">\n");
} else {
ssteps += ($1 - tprev);
nsteps ++;
avgstep = ssteps/nsteps;
}
printf("%s\n", $0);
tprev = $1;
}'
}
function helpmsg {
cat <<- END 1>&2
usage: plot_baro -d infile [ options ]
-d infile set data file name (Required)
-w width set plot width (integer, default = 30 cm)
-h height set plot height (integer, default = 8 cm)
-p pmin set minimum prs (integer, default = 950 hPa, 0=auto)
-P pmax set maximum prs (integer, default = 1050 hPa, 0=auto)
-t tmin set start time (gdate parseable, default = 1 week ago)
-n ndays set num of days (integer, default = 7)
-a arm set barograph length (integer, default = 15 cm)
-c ccol set chart colour (GMT parseable, default = 84/189/140)
-i icol set ink colour (GMT parseable, default = 60/60/60)
-s stn set station name (default = Brisbane)
-o outfile set output file name (default = barograph.*)
-h print help message
END
exit 2
}
#
# Read the command line.
#
#
# Default values.
#
pmin=0 # Automatic
pmax=0 # Automatic
arm=15
plen=30
phgt=8
tmin=$(gdate -u -d "2022-01-12 00:00:00" +%s)
ndays=7
chartcol="84/189/140" # A darkish green.
pencol="60/60/60" # Something like India ink.
stn=""
ifile=""
ofile="barograph"
while getopts ":w#h#p#P#n#a#t:c:i:s:o:d:h" option
do
case "${option}" in
w) plen="${OPTARG}";;
h) phgt="${OPTARG}";;
p) pmin="${OPTARG}";;
P) pmax="${OPTARG}";;
n) ndays="${OPTARG}";;
a) arm="${OPTARG}";;
t) tmin=$(gdate -u -d "${OPTARG}" +%s) || helpmsg;;
c) chartcol="${OPTARG}";;
i) pencol="${OPTARG}";;
s) stn="${OPTARG}";;
d) ifile="${OPTARG}";;
o) ofile="${OPTARG}";;
h) helpmsg;;
?) helpmsg;;
esac
done
#
# The data file is mandatory.
#
if [[ ! -f "${ifile}" ]]
then
helpmsg
fi
#
# Make a temporary directory for intermediate files.
#
tmpdir=$(mktemp -d)
#
# If station name not given, try and extract it from the data file.
#
if [[ "${stn}" = "" ]]
then
stn=$(awk -F, '
BEGIN{
stn="";
}
/^[[:upper:]]{4}/{
stn=$2;
}
END{
gsub("_", " ", stn);
print stn;
}' < "${ifile}")
fi
#
# Set some parameters based on length of chart.
#
# tinc: increment between time grid lines (s)
# linc: increment between pressure labels (s)
#
if [[ ${ndays} -eq 1 ]]
then
tinc=3600
linc=86400
else
tinc=10800
linc=172800
fi
#
# Adjust the start and end times of the chart so
# that it starts on multiples of tinc.
#
tmin=$(((${tmin}/${tinc})*${tinc}))
tmax=$((${tmin}+${ndays}*86400))
#
# Extract the data.
#
# This is set up to extract data from an example climate zone output.
#
#gawk -F, -v tmin="${tmin}" -v tmax="${tmax}" '
# /^hd/ && ($13+0 > 500){
# utime = mktime($8" "$9" "$10" "$11" "$12" 00", 1);
# if ((utime >= tmin)&&(utime <= tmax)){
# printf("%d\t%f\n", utime, $13);
# }
# }' < "${ifile}" > "${tmpdir}/data.txt"
decode_us < "${ifile}" | find_missing "${tmin}" "${tmax}" > "${tmpdir}/data.txt"
cp "${tmpdir}/data.txt" data.txt
#
# Settings for the pressure axis.
#
if [[ ${pmin} -eq 0 ]]
then
pmin=$(gmt info -C -I5 "${tmpdir}/data.txt" | cut -f 3)
fi
if [[ ${pmax} -eq 0 ]]
then
pmax=$(gmt info -C -I5 "${tmpdir}/data.txt" | cut -f 4)
fi
if [[ $((${pmax}-${pmin})) -le 40 ]]
then
if [[ $((${pmax}-${pmin})) -le 20 ]]
then
pinc=2
else
pinc=5
fi
else
pinc=10
fi
#
# Compute the coordinates of the lines on the chart.
#
# Bold pressure lines are drawn every ${pinc} hPa
# Fine pressure lines are drawn every 1 hPa
# Fine time lines are drawn every ${tinc} s
# Bold time lines are drawn every 86400 s
#
prs=${pmin}
while [[ ${prs} -le ${pmax} ]]
do
if [[ $((${prs}%${pinc})) -eq 0 ]]
then
echo "> -Wthicker,${chartcol}"
else
echo "> -Wthinnest,${chartcol}"
fi
echo "${tmin} ${prs}" | barocoords ${plen} ${phgt} ${pmin} ${pmax} \
${tmin} ${tmax} ${arm}
echo "${tmax} ${prs}" | barocoords ${plen} ${phgt} ${pmin} ${pmax} \
${tmin} ${tmax} ${arm}
prs=$((${prs}+1))
done > "${tmpdir}/chart.txt"
t=${tmin}
while [[ ${t} -le ${tmax} ]]
do
if [[ $((${t}%86400)) -eq 0 ]]
then
echo "> -Wthicker,${chartcol}"
else
echo "> -Wthinnest,${chartcol}"
fi
prs=${pmin}
while [[ ${prs} -le ${pmax} ]]
do
echo "${t} ${prs}"
prs=$((${prs}+1))
done | barocoords ${plen} ${phgt} ${pmin} ${pmax} ${tmin} ${tmax} ${arm}
t=$((${t}+${tinc}))
done >> "${tmpdir}/chart.txt"
#
# Compute the positions of the labels.
#
tophourpos=${pmax}
bothourpos=${pmin}
datepos=$((${pmin} + (${pmax}-${pmin})*0.95))
firstdate="Y"
t=${tmin}
while [[ ${t} -le ${tmax} ]]
do
thour=$(((${t}%86400)/3600))
printf "%f\t%f\t10p,Times-Italic,%s\t0\tCB\t%02d\n" \
${t} ${tophourpos} ${chartcol} ${thour}
printf "%f\t%f\t10p,Times-Italic,%s\t0\tCB\t%02d\n" \
${t} ${bothourpos} ${chartcol} ${thour}
if [[ ${thour} -eq 12 ]]
then
if [[ ${firstdate} = "Y" ]]
then
datestr=$(gdate -u -d "1970-01-01 ${t} seconds" +"%a %d %b %Y %Z")
firstdate="N"
else
datestr=$(gdate -u -d "1970-01-01 ${t} seconds" +"%a %d %b")
fi
printf "%f\t%f\t12p,Times-Italic,%s\t0\tCM\t%s\n" \
${t} ${datepos} ${chartcol} "${datestr}"
fi
if [[ $((${t}%linc)) -eq 0 && ${t} -lt $((${tmax}-${tinc})) ]]
then
prs=$((${pmin}+${pinc}))
while [[ ${prs} -le $((${pmax}-${pinc})) ]]
do
pnext=$((${prs}+${pinc}))
if [[ ${pnext} -gt $((${pmax}-${pinc})) ]]
then
printf "%f\t%f\t11p,Helvetica,%s\t0\tCM\t%s\n" \
$((${t}+${tinc})) ${prs} ${chartcol} "${prs} hPa"
else
printf "%f\t%f\t11p,Helvetica,%s\t0\tCM\t%s\n" \
$((${t}+${tinc})) ${prs} ${chartcol} "${prs}"
fi
prs=${pnext}
done
fi
t=$((${t}+${tinc}))
done | barocoords ${plen} ${phgt} ${pmin} ${pmax} \
${tmin} ${tmax} ${arm} > "${tmpdir}/labels.txt"
#
# Add the station name and time stamp to the labels file.
#
cat <<- END >> "${tmpdir}/labels.txt"
0 $((${phgt}/2+1.2)) 12p,Helvetica-Oblique,${chartcol} 0 LB STATION: ${stn}
END
#$((${plen}-2)) $((${phgt}/2+1.2)) 8p,Helvetica-Oblique,${chartcol} 0 BR Plotted by Tim Hume
#
# Now create the chart.
#
gmt begin "${ofile}" pdf,png A+m0.5c
gmt gmtset FONT_LOGO auto,Helvetica,${chartcol}
gmt plot "${tmpdir}/chart.txt" -JX${plen}/${phgt} ${chartreg} -N
gmt text "${tmpdir}/labels.txt" -F+f+a+j -N -DJ0/0.2 -Gwhite
barocoords ${plen} ${phgt} ${pmin} ${pmax} ${tmin} ${tmax} \
${arm} < "${tmpdir}/data.txt" | gmt plot -Wthinnest,${pencol}
#gmt logo -DjTR+w1.5+jBR+o0/0.5 -Sn
gmt end
#
# Remove the temporary directory.
#
rm -rf "${tmpdir}"
And the data file has this format:
station,station_name,valid(UTC),pres1,
PHNL,HONOLULU_INTL/OAHU_(ASOS),2022-01-01 10:00,29.823,
PHNL,HONOLULU_INTL/OAHU_(ASOS),2022-01-01 10:01,29.822,
PHNL,HONOLULU_INTL/OAHU_(ASOS),2022-01-01 10:02,29.822,
PHNL,HONOLULU_INTL/OAHU_(ASOS),2022-01-01 10:03,29.821,
PHNL,HONOLULU_INTL/OAHU_(ASOS),2022-01-01 10:04,29.82,
PHNL,HONOLULU_INTL/OAHU_(ASOS),2022-01-01 10:05,29.82,
PHNL,HONOLULU_INTL/OAHU_(ASOS),2022-01-01 10:06,29.819,
PHNL,HONOLULU_INTL/OAHU_(ASOS),2022-01-01 10:07,29.818,
PHNL,HONOLULU_INTL/OAHU_(ASOS),2022-01-01 10:08,29.818,
