--- a/RDI_Coords.pl
+++ b/RDI_Coords.pl
@@ -1,9 +1,9 @@
#======================================================================
# R D I _ C O O R D S . P L
# doc: Sun Jan 19 17:57:53 2003
-# dlm: Wed Mar 28 12:30:12 2018
+# dlm: Mon Jun 29 10:59:01 2020
# (c) 2003 A.M. Thurnherr
-# uE-Info: 109 22 NIL 0 0 72 10 2 4 NIL ofnI
+# uE-Info: 61 83 NIL 0 0 72 10 2 4 NIL ofnI
#======================================================================
# RDI Workhorse Coordinate Transformations
@@ -57,6 +57,8 @@
# Nov 26, 2017: - BUG: velBeamtoBPEarth() did not respect missing values
# Nov 27, 2017: - BUG: numbersp() from [antslib.pl] was used
# Mar 28, 2018: - added &loadInstrumentTransformation()
+# Jun 5, 2020: - added sscorr_w & sscorr_w_mooring
+# Jun 29, 2020: - added comments for sscorr_w, which conflicts with LADCP_w_ocean
use strict;
use POSIX;
@@ -673,5 +675,61 @@
}
#----------------------------------------------------------------------
+# Sound Speed Correction for Vertical Velocity
+# - Usage: sscorr_w(<observed w>,<beam_angle>,<ADCP sVel setup>,salin,temp,press,<vertical temp. gradient>
+#----------------------------------------------------------------------
+
+sub sscorr_w($$$$$$$$)
+{
+ my($w,$beamAngle,$ssADCP,$salin,$temp,$press,$dz,$dtdz) = @_;
+
+ return 'nan' unless numberp($w);
+ my($tanSqBeamAngle) = tan(rad($beamAngle))**2;
+
+ my($ssXD) = sVel($salin,$temp,$press);
+ my($ssBin) = sVel($salin,$temp+$dz*$dtdz,$press-$dz);
+ my($Kn) = sqrt(1 + (1 - $ssBin/$ssXD)**2 * $tanSqBeamAngle);
+ return $w * $ssBin/$ssADCP / $Kn;
+}
+
+#----------------------------------------------------------------------
+# Sound Speed Correction for Vertical Velocity
+# - Usage: sscorr_w_mooring(\$ADCP,<ens-idx>,<bin-idx>,<press>,<salin>,<vertical temp. gradient>)
+# - Notes:
+# - RDI Coord. Trans. manual sec. 4.1
+# - manual error: the ^2 applies to the []
+# - difference between pressure and depth over instrument range ignored
+# - Assumptions:
+# - libEOS83.pl loaded
+# - $ADCP{ENSEMBLE}[$ens-idx]->VELOCITY}[2] contains measured w
+# - sound speed variation dominated by temperature
+# - vertical temperature gradient is constant in time (violated
+# at least on superinertial time scales)
+#----------------------------------------------------------------------
+
+sub sscorr_w_mooring($$$$$)
+{
+ my($ADCP,$ei,$bi,$press,$salin,$dtdz) = @_;
+
+ my($w) = $ADCP->{ENSEMBLE}[$ei]->{VELOCITY}[2];
+ return 'nan' unless numberp($w);
+
+ my($tanSqBeamAngle) = tan(rad($ADCP->{BEAM_ANGLE}))**2;
+
+ $Global::P{ITS} = 90;
+ my($ssXD) = sVel($salin,$ADCP->{ENSEMBLE}[$ei]->{TEMPERATURE},$press);
+ my($ssADCP) = $ADCP->{ENSEMBLE}[$ei]->{SPEED_OF_SOUND};
+
+ my($dz) = $ADCP->{ENSEMBLE}[$ei]->{BLANKING_DISTANCE} + $bi * $ADCP->{ENSEMBLE}[$ei]->{BIN_LENGTH};
+ $dz *= -1 if ($ADCP->{ENSEMBLE}[$ei]->{XDUCER_FACING_DOWN}); # z increases upward
+
+ my($ssBin) = sVel($salin,
+ $ADCP->{ENSEMBLE}[$ei]->{TEMPERATURE} + $dz*$dtdz,
+ $press-$dz); # ignore press/depth difference across range
+
+ my($Kn) = sqrt(1 + (1 - $ssBin/$ssXD)**2 * $tanSqBeamAngle); # RDI manual
+ return $w * $ssBin/$ssADCP / $Kn;
+}
+
1;
--- a/listBins
+++ b/listBins
@@ -2,9 +2,9 @@
#======================================================================
# L I S T B I N S
# doc: Fri Aug 25 15:57:05 2006
-# dlm: Thu Feb 13 10:36:24 2020
+# dlm: Fri Jun 5 13:45:23 2020
# (c) 2006 A.M. Thurnherr
-# uE-Info: 133 33 NIL 0 0 72 0 2 4 NIL ofnI
+# uE-Info: 367 103 NIL 0 0 72 0 2 4 NIL ofnI
#======================================================================
# Split data file into per-bin time series.
@@ -74,6 +74,7 @@
# Jun 13, 2018: - adpated to RTI files (disabled BIT error check)
# - BUG: dn did not have sufficient digits
# Feb 13, 2020: - added -z
+# May 11, 2020: - removed -z, added -t -m
# General Notes:
# - everything (e.g. beams) is numbered from 1
@@ -118,19 +119,18 @@
require "$ANTS/ants.pl";
require "$ANTS/libconv.pl";
-die("Usage: $0 [-z) progress dots] [-r)ange <first_ens,last_ens>] [-l)ast <bin>] [-R)enumber ensembles from 1] " .
+die("Usage: $0 [-r)ange <first_ens,last_ens>] [-l)ast <bin>] [-R)enumber ensembles from 1] " .
"[-o)utput <redirection[>bin%d.raw]>] " .
"[output -a)ll ens (not just those with good vels)] " .
"[-M)agnetic <declination>] " .
"[-S)oundspeed correction <salin|*,temp|*,depth|*> " .
"[Instrument -T)ransformation Matrix <file>] " .
+ '[disable bin -m)apping] [use TRDI beam-to-earth -t)ransformation] ' .
"[-P)itch/Roll <bias/bias>] [-B)eamvel <bias/bias/bias/bias>] " .
"[require -4)-beam solutions] [-d)iscard <beam#>] " .
"[-p)ct-good <min>] " .
"<RDI file>\n")
- unless (&getopts("4aB:d:l:M:o:p:r:P:RS:T:z") && @ARGV == 1);
-
-$global::readDataProgress = 10000 if defined($opt_z);
+ unless (&getopts("4aB:d:l:mM:o:p:r:P:RS:tT:") && @ARGV == 1);
($P{pitch_bias},$P{roll_bias}) = split('[,/]',$opt_P);
($P{velbias_b1},$P{velbias_b2},$P{velbias_b3},$P{velbias_b4}) = split('[,/]',$opt_B);
@@ -140,7 +140,9 @@
$opt_p = 0 unless defined($opt_p);
-$RDI_Coords::minValidVels = 4 if ($opt_4); # no 3-beam solutions
+$RDI_Coords::minValidVels = 4 if ($opt_4); # no 3-beam solutions
+$RDI_Coords::binMapping = 'none' if ($opt_m); # 'linterp' is default
+$RDI_Coords::beamTransformation = 'RDI' if ($opt_t); # 'LHR90' is default
print(STDERR "WARNING: magnetic declination not set!\n")
unless defined($opt_M);
@@ -362,7 +364,7 @@
@{$dta{ENSEMBLE}[$e]->{BEAM_VELOCITY}[$b]} = # save beam velocities
@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]};
- @{$dta{ENSEMBLE}[$e]->{BEAMPAIR_VELOCITY}[$b]} = # calculate w12, w34
+ @{$dta{ENSEMBLE}[$e]->{BEAMPAIR_VELOCITY}[$b]} = # calculate v12, w12, v34, w34
velBeamToBPEarth(\%dta,$e,@{$dta{ENSEMBLE}[$e]->{BEAM_VELOCITY}[$b]});
@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]} = # calculate earth velocities
--- a/listEns
+++ b/listEns
@@ -2,9 +2,9 @@
#======================================================================
# L I S T E N S
# doc: Sat Jan 18 18:41:49 2003
-# dlm: Thu Feb 13 10:36:50 2020
+# dlm: Sun Feb 28 15:46:21 2021
# (c) 2003 A.M. Thurnherr
-# uE-Info: 96 54 NIL 0 0 72 2 2 4 NIL ofnI
+# uE-Info: 351 1 NIL 0 0 72 2 2 4 NIL ofnI
#======================================================================
$synopsis = 'list ensemble summaries (default), dump ensembles (-E), time-average ensembles (-T)';
@@ -61,6 +61,7 @@
# Apr 10, 2018: - added -l)ast bin
# May 31, 2018: - BUG: -A was disabled by default
# Feb 13, 2020: - added support for $readDataProgress
+# Feb 19, 2021: - BUG: -T did not handle new years correctly
# Notes:
# - -E/-B outputs data in earth coordinates, unless -b is set also
@@ -311,19 +312,18 @@
} elsif (defined($opt_T)) { # time-average ensembles
my(@tmp) = split(',',$opt_T); # decode -T
- my($Tstart,$deltaT,$Tend,$month,$day,$year);
+ my($Tstart,$deltaT,$Tend,$month,$day); # NB: $yearbase needs to be global!
if (@tmp == 3) {
($Tstart,$deltaT,$Tend) = @tmp;
} elsif (@tmp == 2) {
($Tstart,$deltaT) = @tmp;
- ($month,$day,$year) = split('/',$dta{ENSEMBLE}[$#{$dta{ENSEMBLE}}]->{DATE});
- $Tend = str2dec_time($dta{ENSEMBLE}[$#{$dta{ENSEMBLE}}]->{DATE},$dta{ENSEMBLE}[$#{$dta{ENSEMBLE}}]->{TIME},$year);
+ ($month,$day,$yearbase) = split('/',$dta{ENSEMBLE}[0]->{DATE});
+ $Tend = str2dec_time($dta{ENSEMBLE}[0]->{DATE},$dta{ENSEMBLE}[$#{$dta{ENSEMBLE}}]->{TIME},$yearbase);
} else {
- ($month,$day,$year) = split('/',$dta{ENSEMBLE}[0]->{DATE});
- $Tstart = str2dec_time($dta{ENSEMBLE}[0]->{DATE},$dta{ENSEMBLE}[0]->{TIME},$year);
+ ($month,$day,$yearbase) = split('/',$dta{ENSEMBLE}[0]->{DATE});
+ $Tstart = str2dec_time($dta{ENSEMBLE}[0]->{DATE},$dta{ENSEMBLE}[0]->{TIME},$yearbase);
($deltaT) = @tmp;
- ($month,$day,$year) = split('/',$dta{ENSEMBLE}[$#{$dta{ENSEMBLE}}]->{DATE});
- $Tend = str2dec_time($dta{ENSEMBLE}[$#{$dta{ENSEMBLE}}]->{DATE},$dta{ENSEMBLE}[$#{$dta{ENSEMBLE}}]->{TIME},$year);
+ $Tend = str2dec_time($dta{ENSEMBLE}[$#{$dta{ENSEMBLE}}]->{DATE},$dta{ENSEMBLE}[$#{$dta{ENSEMBLE}}]->{TIME},$yearbase);
}
$Tstart = &{&antsCompileConstExpr($')} if ($Tstart =~ m{^=});
$deltaT = &{&antsCompileConstExpr($')} if ($deltaT =~ m{^=});
@@ -346,11 +346,12 @@
my($e) = @_;
my($b,$i);
- my($month,$day,$year) = ($e >= 0) ? split('/',$dta{ENSEMBLE}[$e]->{DATE}) : (undef,undef,undef);
- my($dn) = ($e >= 0) ? str2dec_time($dta{ENSEMBLE}[$e]->{DATE},$dta{ENSEMBLE}[$e]->{TIME},$year) : undef;
+ my($dn) = ($e >= 0) ? str2dec_time($dta{ENSEMBLE}[$e]->{DATE},$dta{ENSEMBLE}[$e]->{TIME},$yearbase) : undef;
+# print(STDERR "ens#$e at $dn (cbin = $cbin)\n");
if ($e<0 || $dn>=$Tstart+$cbin*$deltaT) { # dump full bin
my($file) = sprintf("$basename.T$fnfmt",$cbin); # file name: <basename>.T0001
+# print(STDERR "dumping average to $file...\n");
do { # skip empy bins
$cbin++;
@@ -530,6 +531,8 @@
}
} # define output format
+die("$yearbase") unless ($yearbase > 1900);
+
#----------------------------------------------------------------------
# Loop Over Ensembles
#----------------------------------------------------------------------