#* Draw a TS diagram, possibly choosing isopycnals automatically
#* @param $rho_inc (optional), density step between isopycnals
cmd draw TS diagram (;$) { # [$rho_inc]
my ($rho, $rho_min, $rho_max, $rho_inc, $rho_inc_range, $S, $dS);
set x name "Salinity, PSU";
set y name "Temperature, \(\circ\)C";
$rho_min = eos_dens($_xleft, $_ytop, 0); # min density on diagram
$rho_max = eos_dens($_xright, $_ybottom, 0); # max density on diagram
if ($#_ == 0) {
# Use specified density increment.
$rho_inc = abs($_[0]);
} else {
# Didn't give an increment, so pick one, using 1-2-5 scaling
# to give of the order of 5 isopycnals in range of data.
$rho_inc = ($rho_max - $rho_min) / 5;
$rho_inc_range = exp10(floor(log10($rho_inc))); # a power of 10
if ($rho_inc / $rho_inc_range < 2) {
$rho_inc = $rho_inc_range;
} else {
if ($rho_inc / $rho_inc_range < 5) {
$rho_inc = $rho_inc_range * 2;
} else {
$rho_inc = $rho_inc_range * 5;
}
}
}
$rho_min = $rho_inc * floor($rho_min / $rho_inc);
$rho_max = $rho_inc * ceil($rho_max / $rho_inc);
for ($rho = $rho_min; $rho <= $rho_max; $rho += $rho_inc) {
draw isopycnal $rho;
}
if (eos_freezing_temp($_xleft, 0) > $_ybottom) {
draw line connecting $_xleft eos_freezing_temp($_xleft, 0)
$_xright eos_freezing_temp($_xright, 0);
}
draw symbols "bullet";
# Draw freezing-point line
$dS = ($_xright - $_xleft) / 50;
for ($S = $_xleft; $S <= $_xright; $S += $dS) {
draw line connecting ($S-$dS) eos_freezing_temp($S-$dS, 0)
$S eos_freezing_temp($S, 0);
}
}
$file = "/users/dek/bmay/thesis/eubex/data/stn308";
open(IN, "$file") or die "cannot access $file\n";
read columns IN x=3 y=2;
draw TS diagram;
#draw TS diagram 1;