/usr/bin/bat2gts is in gerris 20131206+dfsg-17.
This file is owned by root:root, with mode 0o755.
The actual contents of the file can be viewed below.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 | #!/bin/sh
# reference longitude (degrees)
rlong=174
# reference latitude (degrees)
rlat=-41
# width of the domain (kilometers)
width=500
# coastline reference depth (m)
coast=0.1
# verbose output
verbose=""
# angle of rotation
angle=0
# relative error
relative=0.05
# reference depth
H=5000
usage()
{
cat <<EOF
Usage: bat2gts [OPTIONS] < BATHY
Takes a bathymetry file (three columns: longitude, latitude (degree),
depth (meters)) and generates a GTS depth file for the given
domain.
Options:
[--long=V] set reference longitude to V (default 174 deg)
[--lat=V] set reference latitude to V (default -41 deg)
[--width=V] set domain width (default is 500 km)
[--depth=V] set reference depth to V
(default is 5000 meters)
[--coast=V] set coastline reference depth to V
(default is 0.1 meters)
[--rel=T] set relative error allowed on bathymetry (default is 0.05)
[--angle=V] rotation of V degrees (default is 0)
[--verbose] display info about the process
[--help] display this message and exits
EOF
exit $1
}
if test $# -lt 1; then
usage 1 1>&2
fi
command="bat2gts $*"
while test $# -gt 0; do
case "$1" in
-*=*) optarg=`echo "$1" | sed 's/[-_a-zA-Z0-9]*=//'` ;;
*) optarg= ;;
esac
case $1 in
--long=*)
rlong=$optarg
;;
--lat=*)
rlat=$optarg
;;
--width=*)
width=$optarg
;;
--depth=*)
H=$optarg
;;
--coast=*)
coast=$optarg
;;
--angle=*)
angle=$optarg;
;;
--rel=*)
relative=$optarg
;;
--verbose)
verbose="-v"
;;
--help)
usage 0 1>&2
;;
*)
usage 0 1>&2
;;
esac
shift
done
proj=-Ja$rlong/$rlat/1:`echo $width | awk '{print $1*1e5}'`c
area=-R`echo $rlong $rlat $width | awk '{
Pi = 3.14159265359
dlon = 0.71*$3/(6378.*cos ($2*Pi/180))*180./Pi
dlat = 0.71*$3/6378.*180./Pi
print $1 - dlon "/" $1 + dlon "/" $2 - dlat "/" $2 + dlat
}'`
cat <<EOF
# Created by bat2gts
# wdir: $PWD
# command line: $command
# reference depth: $H m
EOF
mapproject -Dc -C $proj $area |\
awk -v H=$H -v coast=$coast '{
if ($1 < 0.71 && $1 > -0.71 && $2 < 0.71 && $2 > -0.71)
print $1 " " $2 " " (coast - $3)/H;
}
' | sort | uniq | happrox -f $verbose -r `awk -v H=$H 'BEGIN{print 1./H}'` -c $relative |\
transform --rz $angle | transform --revert --tz 0.5
cat <<EOF > lolat2xy
gmtset D_FORMAT = %.12lf
mapproject -Dc -C $proj $area | awk '
BEGIN {
angle = $angle * 3.14159265359/180.;
cosa = cos (angle);
sina = sin (angle);
}
{
if (NF >= 2) {
printf ("%f %f", cosa*\$1 - sina*\$2, sina*\$1 + cosa*\$2);
for (i = 3; i <= NF; i++)
printf (" %s", \$i);
printf ("\n");
}
else
print \$0;
}'
EOF
chmod +x lolat2xy
cat <<EOF > xy2lolat
gmtset D_FORMAT = %.12lf
awk '
BEGIN {
angle = - $angle * 3.14159265359/180.;
cosa = cos (angle);
sina = sin (angle);
}
{
if (NF >= 2) {
printf ("%f %f", cosa*\$1 - sina*\$2, sina*\$1 + cosa*\$2);
for (i = 3; i <= NF; i++)
printf (" %s", \$i);
printf ("\n");
}
else
print \$0;
}' | mapproject -Dc -I -C $proj $area
EOF
chmod +x xy2lolat
|