# HG changeset patch
# User Andreas Thurnherr <ant@ldeo.columbia.edu>
# Date 1540166508 14400
# Node ID 5767cbe470a0479734f4bf29c0e6a7234f2d988e
# Parent  b7654ea68af61ded4512415f81e7cc97ef38e805# Parent  b63fa355644c01a01552ae8f3632504366a2f2be
merged whoosher & EN620

diff --git a/ADCP_tools_lib.pl b/ADCP_tools_lib.pl
--- a/ADCP_tools_lib.pl
+++ b/ADCP_tools_lib.pl
@@ -1,9 +1,9 @@
 #======================================================================
 #                    A D C P _ T O O L S _ L I B . P L 
 #                    doc: Tue Jan  5 10:45:47 2016
-#                    dlm: Thu Dec  7 10:43:28 2017
+#                    dlm: Tue Feb  6 21:37:45 2018
 #                    (c) 2016 A.M. Thurnherr
-#                    uE-Info: 17 25 NIL 0 0 72 0 2 4 NIL ofnI
+#                    uE-Info: 16 57 NIL 0 0 72 0 2 4 NIL ofnI
 #======================================================================
 
 # HISTORY:
@@ -13,8 +13,9 @@
 #	Mar 12, 2017: - updated to V1.9 for LADCP_w 1.3
 #	Nov 28, 2017: - updated to V2.0 for LADCP_w 1.4
 #	Dec  7, 2017: - updated to V2.1 for improvements to listHdr
+#	Feb  6, 2018: - updated to V2.2 for changes to PD0_IO
 
-$ADCP_tools_version = 2.1;
+$ADCP_tools_version = 2.2;
 
 die(sprintf("$0: obsolete ADCP_tools V%.1f; V%.1f required\n",
     $ADCP_tools_version,$ADCP_tools_minVersion))
diff --git a/HISTORY b/HISTORY
--- a/HISTORY
+++ b/HISTORY
@@ -1,9 +1,9 @@
 ======================================================================
                     H I S T O R Y 
                     doc: Tue May 15 18:04:39 2012
-                    dlm: Sat Dec 23 16:13:38 2017
+                    dlm: Wed Mar 28 12:21:58 2018
                     (c) 2012 A.M. Thurnherr
-                    uE-Info: 214 15 NIL 0 0 72 3 2 4 NIL ofnI
+                    uE-Info: 237 60 NIL 0 0 72 3 2 4 NIL ofnI
 ======================================================================
 
 --------------------------------------
@@ -212,3 +212,26 @@
 	- updated all tools to use MinVersion 2.1
 	- updated [patchPD0] to use ANTSlib V7.0
 	- PUBLISHED
+
+#---------------------------------------------------------------
+# V2.2
+#	- allow interior garbage in ADCP files
+#	- allow use of individual Instrument Transformation Matrices
+#---------------------------------------------------------------
+
+Feb  6, 2018:
+	- updated to V2.2 [ADCP_tools_lib.pl]
+	- support for partial files in RDI_PD0_IO, listBins
+
+Feb  7, 2018:
+	- added support for garbage inside PD0 files
+	- improvement to listEns
+
+Mar 14-20, 2018:
+	- added consistency check to mk_prof()
+	- fixed bugs in RDI_PD0_IO
+	- updated [HISTORY]
+
+Mar 28, 2018:
+	- added &loadInstrumentTransformation() to [RDI_Coords.pl]
+	- added support for &loadInstrumentTransformation() to [listBins]
diff --git a/RDI_Coords.pl b/RDI_Coords.pl
--- 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: Mon Nov 27 07:13:25 2017
+#                    dlm: Wed Mar 28 12:30:12 2018
 #                    (c) 2003 A.M. Thurnherr
-#                    uE-Info: 58 62 NIL 0 0 72 10 2 4 NIL ofnI
+#                    uE-Info: 109 22 NIL 0 0 72 10 2 4 NIL ofnI
 #======================================================================
 
 # RDI Workhorse Coordinate Transformations
@@ -56,6 +56,7 @@
 #	Oct 12, 2017: - documentation
 #	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()
 
 use strict;
 use POSIX;
@@ -74,7 +75,15 @@
 $RDI_Coords::beamTransformation = 'LHR90';	# set to 'RDI' to use 1st order transformations from RDI manual
 
 #----------------------------------------------------------------------
-# beam to earth transformation 
+# beam to earth transformation
+#	- loadInstrumentTransformation(filename) loads a file that contains the
+#	  output from the PS3 command, which includes the instrument transformation
+#	  matrix as follows:
+#	Instrument Transformation Matrix (Down):    Q14:
+#	  1.4689   -1.4682    0.0030   -0.0035       24067  -24055      49     -58
+#	 -0.0036    0.0029   -1.4664    1.4673         -59      48  -24025   24041
+#	  0.2658    0.2661    0.2661    0.2657        4355    4359    4359    4354
+#	  1.0373    1.0382   -1.0385   -1.0373       16995   17010  -17015  -16995
 #----------------------------------------------------------------------
 
 $RDI_Coords::threeBeam_1 = 0;			# stats from velBeamToInstrument
@@ -88,6 +97,35 @@
 { # STATIC SCOPE
 	my(@B2I);
 
+	sub loadInstrumentTransformation($)
+	{
+		die("loadInstrumentTransformation(): B2I matrix already defined\n")
+			if (@B2I);
+		open(ITF,$_[0]) || die("$_[0]: $!\n");
+		my($row) = 0;
+		while (<ITF>) {
+			if ($row == 0) {
+				next unless m{^Instrument Transformation Matrix \(Down\):};
+				$row = 1;
+			} elsif ($row <= 4) {
+				my(@vals) = split;
+				die("$_[0]: cannot decode row #$row of Instrument Transformation Matrix\n")
+					unless (@vals == 8);
+				for (my($i)=0; $i<4; $i++) {
+					die("$_[0]: cannot decode row #$row of Instrument Transformation Matrix\n")
+						unless numberp($vals[$i]);
+					$B2I[$row-1][$i] = $vals[$i];
+				}
+				$row++;
+			} else {
+				last;
+			}
+		}
+		die("$_[0]: cannot decode Instrument Transformation Matrix (row = $row)\n")
+			unless ($row == 5);
+		close(ITF);
+	}
+
 	sub velBeamToInstrument(@)
 	{
 		my($ADCP,$ens,$v1,$v2,$v3,$v4) = @_;
@@ -134,22 +172,21 @@
 	}
 } # STATIC SCOPE
 
-#----------------------------------------------------------------------
+#--------------------------------------------------------------------------------------------------------------
 # velInstrumentToEarth(\%ADCP,ens,v1,v2,v3,v4) => (u,v,w,e)
 #	- $RDI_Coords::beamTransformation = 'LHR90'
 #		- from Lohrmann, Hackett & Roet (J. Tech., 1990)
 #		- eq A1 maps to RDI matrix M (sec 5.6) with
 #			alpha = roll
 #			beta = gimball_pitch
-#			psi = calculation_pitch
-#			psi = asin{sin(beta) cos(alpha) / sqrt[1- sin(alpha)^2 sin(beta)^2]}
+#			psi (pitch used for calculation) =  asin{sin(beta) cos(alpha) / sqrt[1- sin(alpha)^2 sin(beta)^2]}
 #		- (I only checked for 0 heading, but this is sufficient)
 #	- $RDI_Coords::beamTransformation = 'RDI'
 #		- default prior to LADCP_w V1.3
 #		- from RDI manual
 #		- 99% accurate for p/r<8deg
 #			=> 1cm/s error for 1m/s winch speed!
-#----------------------------------------------------------------------
+#--------------------------------------------------------------------------------------------------------------
 
 { # STATIC SCOPE
 	my($hdg,$pitch,$roll,@I2E);
diff --git a/RDI_PD0_IO.pl b/RDI_PD0_IO.pl
--- a/RDI_PD0_IO.pl
+++ b/RDI_PD0_IO.pl
@@ -1,9 +1,9 @@
 #======================================================================
 #                    R D I _ P D 0 _ I O . P L 
 #                    doc: Sat Jan 18 14:54:43 2003
-#                    dlm: Fri Feb  2 12:49:54 2018
+#                    dlm: Tue Jun 12 19:10:08 2018
 #                    (c) 2003 A.M. Thurnherr
-#					 uE-Info: 772 0 NIL 0 0 72 2 2 4 NIL ofnI
+#					 uE-Info: 115 72 NIL 0 0 72 2 2 4 NIL ofnI
 #======================================================================
 	
 # Read RDI PD0 binary data files (*.[0-9][0-9][0-9])
@@ -97,6 +97,22 @@
 #	Dec 23, 2017: - BUG: could no longer read Anslope II raw files
 #				  - added support for patching ADCP time data
 #				  - added support for RDI_PD0_IO::OVERRIDE_Y2K_CLOCK
+#	Feb  6, 2018: - added supprort for first_ens, last_ens in readData()
+#	Feb  7, 2018: - BUG: new params were not optional
+#				  - made read routines more permissive by allowing
+#					garbage between ensembles
+#	Mar 15, 2018: - BUG: WBPens() did not work for files with garbage
+#	Mar 16, 2018: - added skipped garbage warning
+#				  - BUG: garbage skipping did not work correctly for files w/o BT
+#	Mar 20, 2018: - BUG: garbage skipping STILL did not work correctly for files w/o BT
+#	Apr  9, 2018: - added last_bin argument to readData()
+#	Apr 10, 2018: - slight improvement (parameter range check)
+#	Apr 23, 2018: - make WBRens() work (again?) with BB150 data froim 1996
+#	Apr 24, 2018: - BUG: undefined lat_bin argument to readData() did not work
+#	Apr 30, 2018: - added support for repeated ensembles
+#				  - added warning on wrong ensemble length
+#	Jun  9, 2018: - removed double \n from warnings
+#	Jun 12, 2018: - BUG: IMPed files did not pass the garbage detection 
 	
 # FIRMWARE VERSIONS:
 #	It appears that different firmware versions generate different file
@@ -119,7 +135,7 @@
 #
 #	- DATA_SOURCE_ID = 0x7F 					original TRDI PD0 file
 #
-#	- DATA_SOURCE_ID = 0xA0 | PATCHED_MASK		produced by IMP+LADP 
+#	- DATA_SOURCE_ID = 0xA0 | PATCHED_MASK		produced by IMP+LADCP, KVH+LADCP
 #		PATCHED_MASK & 0x04:						pitch value has been patched
 #		PATCHED_MASK & 0x02:						roll value has been patched
 #		PATCHED_MASK & 0x01:						heading value has been patched
@@ -357,37 +373,53 @@
 	my($FmtErr) = "%s: illegal %s Id 0x%04x at ensemble %d";
 	
 #----------------------------------------------------------------------
-# skip to first valid ensemble (skip over initial garbage)
+# skip to next valid start of ensemble (skip over garbage)
 #----------------------------------------------------------------------
 	
-	sub skip_initial_trash(@)
-	{
-		my($quiet) = @_;
-		my($buf,$dta);
+sub goto_next_ens(@)
+{
+	my($fh,$return_skipped) = @_;							# if return_skipped not set, return file pos
+	my($buf,$dta);
+
+	my($found) = 0; 										# zero consecutive 0x7f found
+	my($skipped) = 0; my($garbage_start);
+	while ($found < 2) {
+		sysread($fh,$buf,1) == 1 || last; 
+		($dta) = unpack('C',$buf);
+		if ($dta == 0x7f) {
+			$found++;
+		} elsif ($found==1 &&
+					($dta==0xE0	||									# from editPD0
+					 (($dta&0xF0)==0xA0 && ($dta&0x0F)<8))) {		# from IMP+LADCP or KVH+LADCP
+			$found++;
+		} elsif ($found == 0) {
+			$garbage_start = sysseek($fh,0,1)-1 unless defined($garbage_start);
+			$skipped++;
+		} else {											# here, found == 1 but 2nd byte was not found
+			$garbage_start = sysseek($fh,0,1)-$found unless defined($garbage_start);
+			$skipped += $found;
+			$found = 0;
+		}
+	}
+	my($fpos) = ($found < 2) ? undef : sysseek($fh,-2,1);
+	return $skipped if ($return_skipped);
 	
-		my($found) = 0; 									# zero consecutive 0x7f found
-		my($skipped) = 0;
-		while ($found < 2) {
-			sysread(WBRF,$buf,1) == 1 || last;
-			($dta) = unpack('C',$buf);
-			if ($dta == 0x7f) {
-				$found++;
-			} elsif ($found==1 && ($dta==0xE0 || ($dta&0xF0==0xA0 && $dta&0x0F<8))) {
-				$found++;
-			} elsif ($found == 0) {
-				$skipped++;
-			} else {
-				$skipped += $found;
-				$found = 0;
-			}
+	if ($skipped) {
+		if (eof($fh)) {
+#				04/18/18: disabled the following line of code because it is very common at
+#						  least with the older RDI instruments I am looking at in the context
+#						  of the SR1b repeat section analysis
+#			print(STDERR "WARNING (RDI_PD0_IO): PD0 file ends with $skipped garbage bytes\n");
+		} elsif ($garbage_start == 0) {
+			print(STDERR "WARNING (RDI_PD0_IO): PD0 file starts with $skipped garbage bytes\n");
+		} else {
+			print(STDERR "WARNING (RDI_PD0_IO): $skipped garbage bytes in PD0 file beginning at byte $garbage_start\n");
 		}
-		die("$WBRcfn: no valid ensemble header found [$!]\n")
-			if ($found < 2);
-		printf(STDERR "WARNING: %d bytes of initial garbage\n",$skipped)
-			if ($skipped > 0 && !$quiet);
-		return sysseek(WBRF,-2,1);
 	}
 	
+	return $fpos;
+}
+	
 #----------------------------------------------------------------------
 # readHeader(file_name,^dta) WBRhdr(^data)
 #	- read header data
@@ -419,7 +451,10 @@
 	# HEADER
 	#--------------------
 
-	skip_initial_trash();
+	my($skipped) = goto_next_ens(\*WBRF,1);	
+   	printf(STDERR "WARNING: %d bytes of initial garbage\n",$skipped)
+		if ($skipped > 0);
+	
 	sysread(WBRF,$buf,6) == 6 || return undef;
 	($hid,$did,$dta->{ENSEMBLE_BYTES},$dummy,$dta->{NUMBER_OF_DATA_TYPES})
 		= unpack('CCvCC',$buf);
@@ -455,7 +490,7 @@
 			|| die("$WBRcfn: $!");
 	@WBRofs = unpack("v$dta->{NUMBER_OF_DATA_TYPES}",$buf);
 #	for ($i=0; $i<$dta->{NUMBER_OF_DATA_TYPES}; $i++) {
-#		printf(STDERR "\nWBRofs[$i] = %d",$WBRofs[$i]);
+#		printf(STDERR "WBRofs[$i] = %d",$WBRofs[$i]);
 #	}
 
 
@@ -698,33 +733,49 @@
 }
 
 #----------------------------------------------------------------------
-# readData(file_name,^data) WBRens(nbins,fixed_leader_bytes,^data)
-# 	- read all ensembles
+# readData(file_name,^data[,first_ens,last_ens[,last_bin]]) 
+# 	- read ensembles
+#	- read all ensembles unless first_ens and last_ens are given
+#	- read all bins unless last_bin is given
 #----------------------------------------------------------------------
 
 sub readData(@)
 {
-	my($fn,$dta) = @_;
+	my($fn,$dta,$fe,$le,$lb) = @_;
 	$WBRcfn = $fn;
     open(WBRF,$WBRcfn) || die("$WBRcfn: $!\n");
     WBRhdr($dta,1) || die("$WBRcfn: Insufficient Data\n");
-	WBRens($dta->{N_BINS},$dta->{FIXED_LEADER_BYTES},
-		   \@{$dta->{ENSEMBLE}});
+    $lb = $dta->{N_BINS}
+		unless (numberp($lb) && $lb>=1 && $lb<=$dta->{N_BINS});
+	WBRens($lb,$dta->{FIXED_LEADER_BYTES},\@{$dta->{ENSEMBLE}},$fe,$le);
 	print(STDERR "$WBRcfn: $BIT_errors built-in-test errors\n")
 		if ($BIT_errors);
 }
 
-sub WBRens($$$)
+sub WBRens(@)
 {
-	my($nbins,$fixed_leader_bytes,$E) = @_;
+	my($nbins,$fixed_leader_bytes,$E,$fe,$le) = @_;
 	my($start_ens,$B1,$B2,$B3,$B4,$I,$id,$bin,$beam,$buf,$dummy,@dta,$i,$cs,@WBRofs);
-	my($ens,$ensNo,$dayStart,$ens_length,$hid,$did,$ndt);
+	my($ens,$ensNo,$dayStart,$ens_length,$hid,$did,$ndt,$el);
 
     sysseek(WBRF,0,0) || die("$WBRcfn: $!");
-	$start_ens = skip_initial_trash(1);
-	for ($ens=0; 1; $ens++,$start_ens+=$ens_length+2) {
-#		print(STDERR "ens = $ens\n");
-#		print(STDERR "start_ens = $start_ens\n");
+ENSEMBLE:
+	for ($ens=0; 1; $ens++) {
+		$start_ens = goto_next_ens(\*WBRF);
+		last unless defined($start_ens);
+
+		#----------------------------------------
+		# Handle first_ens and last_ens
+		#----------------------------------------
+
+		if (defined($fe) && $ens>0 && ${$E}[$ens-1]->{NUMBER}<$fe) {					# delete previous ensemble
+			pop(@{$E}); $ens--;
+		}
+
+		if (defined($le) && $ens>0 && ${$E}[$ens-1]->{NUMBER}>$le) {					# delete previous ensemble and finish
+			pop(@{$E}); $ens--;
+			last;
+		}
 
 		#----------------------------------------
 		# Get ensemble length and # of data types 
@@ -732,7 +783,7 @@
 
 		sysseek(WBRF,$start_ens,0) || die("$WBRcfn: $!");
 		sysread(WBRF,$buf,6) == 6 || last;
-		($hid,$did,$ens_length,$dummy,$ndt) = unpack('CCvCC',$buf);
+		($hid,$did,$el,$dummy,$ndt) = unpack('CCvCC',$buf);
 		$hid == 0x7f || die(sprintf($FmtErr,$WBRcfn,"Header",$hid,defined($ensNo)?$ensNo+1:0));
 		${$E}[$ens]->{DATA_SOURCE_ID} = $did;
 		if ($did == 0x7f) {
@@ -745,7 +796,16 @@
 			${$E}[$ens]->{PRODUCER} = 'unknown';
 	    }
 
-##		printf(STDERR "\n$WBRcfn: WARNING: unexpected number of data types (%d, ens=$ens)\n",$ndt),last
+		if (defined($ens_length) && ($el != $ens_length)) {
+			print(STDERR "WARNING (RDI_PD0_IO): ensemble ${$E}[$#{$E}]->{NUMBER} skipped (unexpected length)\n");
+			pop(@{$E});
+			$ens--;
+			next;
+		}
+		
+		$ens_length = $el;
+
+##		printf(STDERR "$WBRcfn: WARNING: unexpected number of data types (%d, ens=$ens)\n",$ndt),last
 ##				unless ($ndt == 6 || $ndt == 7);
 		sysread(WBRF,$buf,2*$ndt) == 2*$ndt || die("$WBRcfn: $!");
 		@WBRofs = unpack("v$ndt",$buf);
@@ -763,15 +823,16 @@
 		sysseek(WBRF,$start_ens,0) || die("$WBRcfn: $!");
 		unless ((sysread(WBRF,$buf,$ens_length) == $ens_length) &&				# incomplete ensemble
 				(sysread(WBRF,$cs,2) == 2)) {
+#			print(STDERR "INCOMPLETE ENSEMBLE\n");
 			pop(@{$E});
 			last;
 		}
 
 		unless (unpack('%16C*',$buf) == unpack('v',$cs)) {						# bad checksum
-			pop(@{$E});
-			last;
-#			next;																# using this might make the code work
-		}																		# for files with isolated bad ensembles
+#			print(STDERR "BAD CHECKSUM\n");
+			pop(@{$E}); $ens--;
+			next;
+		}		
 
 		#------------------------------
 		# Variable Leader
@@ -803,6 +864,16 @@
 
 		sysread(WBRF,$buf,1) == 1 || die("$WBRcfn: $!");
 		$ensNo += unpack('C',$buf) << 16;
+
+		for (my($i)=$ens; $i>0; $i--) {											# check for duplicate ens; e.g. 2018 S4P 24UL
+			if (${$E}[$i]->{NUMBER} == $ensNo) {									
+				print(STDERR "WARNING (RDI_PD0_IO): duplicate ensemble $ensNo skipped\n");
+				pop(@{$E});
+				$ens--;
+				next ENSEMBLE;
+			}
+		}
+			
 		${$E}[$ens]->{NUMBER} = $ensNo;
 		
 		sysread(WBRF,$buf,30) == 30 || die("$WBRcfn: $!");
@@ -888,14 +959,14 @@
 			= &_dayNo(${$E}[$ens]->{YEAR},${$E}[$ens]->{MONTH},${$E}[$ens]->{DAY},
 					  ${$E}[$ens]->{HOUR},${$E}[$ens]->{MINUTE},${$E}[$ens]->{SECONDS});
 
-		# when analyzing an STA file from an OS75 SADCP (Poseidion),
+		# when analyzing an STA file from an OS75 SADCP (Poseidon),
 		# I noticed that there is no time information. This causes
 		# timegm to bomb. 
 		if (${$E}[$ens]->{MONTH} == 0) {					# no time info
 			${$E}[$ens]->{UNIX_TIME} = 0;
 			${$E}[$ens]->{SECNO} = 0;
         } else {
-#			print(STDERR "\n[$ens]->${$E}[$ens]->{MINUTE}:${$E}[$ens]->{HOUR},${$E}[$ens]->{DAY},${$E}[$ens]->{MONTH},${$E}[$ens]->{YEAR}-<\n");
+#			print(STDERR "[$ens]->${$E}[$ens]->{MINUTE}:${$E}[$ens]->{HOUR},${$E}[$ens]->{DAY},${$E}[$ens]->{MONTH},${$E}[$ens]->{YEAR}-<\n");
 			${$E}[$ens]->{UNIX_TIME}
 				= timegm(0,${$E}[$ens]->{MINUTE},
 						   ${$E}[$ens]->{HOUR},
@@ -988,12 +1059,9 @@
 			die(sprintf($FmtErr,$WBRcfn,"Percent-Good Data",$id,$ensNo));
 
 		for ($i=0,$bin=0; $bin<$nbins; $bin++) {
-#			printf(STDERR "%-GOOD($bin): ");
 			for ($beam=0; $beam<4; $beam++,$i++) {
-#				printf(STDERR "$dta[$i] ");
 				${$E}[$ens]->{PERCENT_GOOD}[$bin][$beam] = $dta[$i];
 			}
-#			printf(STDERR "\n");
 		}
 
 		#-----------------------------------------
@@ -1009,7 +1077,11 @@
 			last if ($id == 0x0600);
 		}
 
-		next if ($nxt == $ndt);													# no BT found => next ens
+		if ($nxt == $ndt) {														# no BT found => next ens
+#			sysseek(WBRF,4,1) || die("$WBRcfn: $!");							# skip over remainder of ensemble
+			sysseek(WBRF,$start_ens+$ens_length+2,0) || die("$WBRcfn: $!");
+			next;
+		}
 
 		sysseek(WBRF,14,1) || die("$WBRcfn: $!");								# BT range, velocity, corr, %-good, ...
 		sysread(WBRF,$buf,28) == 28 || die("$WBRcfn: $!");
@@ -1054,17 +1126,21 @@
 		}
 
 		sysseek(WBRF,2,1) || die("$WBRcfn: $!");								# BT signal strength & BT range high bytes
-		sysread(WBRF,$buf,9) == 9 || die("$WBRcfn: $!");
-		@dta = unpack('C4CC4',$buf);
-		for ($beam=0; $beam<4; $beam++) {
-			${$E}[$ens]->{BT_SIGNAL_STRENGTH}[$beam] = $dta[$beam];
-		}
-		${$E}[$ens]->{HIGH_GAIN} if    ($dta[4]);
-		${$E}[$ens]->{LOW_GAIN} unless ($dta[4]);
-		for ($beam=0; $beam<4; $beam++) {										# high bytes (1 byte per beam)
-			${$E}[$ens]->{BT_RANGE}[$beam] += $dta[5+$beam] * 655.36
-				if ($dta[5+$beam]);
-		}
+#		sysread(WBRF,$buf,9) == 9 || die("$WBRcfn: $!");
+		if (sysread(WBRF,$buf,9) == 9) {										# SR1b JR16 BB150 data files require this
+			@dta = unpack('C4CC4',$buf);
+			for ($beam=0; $beam<4; $beam++) {
+				${$E}[$ens]->{BT_SIGNAL_STRENGTH}[$beam] = $dta[$beam];
+			}
+			${$E}[$ens]->{HIGH_GAIN} if    ($dta[4]);
+			${$E}[$ens]->{LOW_GAIN} unless ($dta[4]);
+			for ($beam=0; $beam<4; $beam++) {									# high bytes (1 byte per beam)
+				${$E}[$ens]->{BT_RANGE}[$beam] += $dta[5+$beam] * 655.36
+					if ($dta[5+$beam]);
+			}
+#			sysseek(WBRF,8,1) || die("$WBRcfn: $!");							# remainder of ensemble
+        }
+        sysseek(WBRF,$start_ens+$ens_length+2,0) || die("$WBRcfn: $!");
 	} # ens loop
 }
 
@@ -1116,9 +1192,12 @@
 {
 	my($nbins,$fixed_leader_bytes,$dta) = @_;
 	my($start_ens,$B1,$B2,$B3,$B4,$I,$id,$bin,$beam,$buf,$dummy,@dta,$i,$cs,@WBPofs);
-	my($ens,$dayStart,$ens_length,$hid,$ndt);
+	my($ens,$dayStart,$ens_length,$hid,$ndt,$el);
 
-	for ($ens=$start_ens=0; $ens<=$#{$dta->{ENSEMBLE}}; $ens++,$start_ens+=$ens_length+2) {
+#	for ($ens=$start_ens=0; $ens<=$#{$dta->{ENSEMBLE}}; $ens++,$start_ens+=$ens_length+2) {
+	for ($ens=0; $ens<=$#{$dta->{ENSEMBLE}}; $ens++) {
+		$start_ens = goto_next_ens(\*WBPF);
+		die("ens = $ens\n") unless defined($start_ens);
 
 		#------------------------------
 		# Patch Header (Data Source Id)
@@ -1134,8 +1213,11 @@
 		$nw == 1 || die("$WBPcfn: $nw bytes written ($!)");
 
 		sysread(WBPF,$buf,4) == 4 || die("$WBPcfn: unexpected EOF");
-		($ens_length,$dummy,$ndt) = unpack('vCC',$buf);
-		printf(STDERR "\n$WBPcfn: WARNING: unexpected number of data types (%d, ens=$ens)\n",$ndt),last
+		($el,$dummy,$ndt) = unpack('vCC',$buf);
+		$ens--,next if (defined($ens_length) && ($el != $ens_length));
+		$ens_length = $el;
+		
+		printf(STDERR "$WBPcfn: WARNING: unexpected number of data types (%d, ens=$ens)\n",$ndt),last
 				unless ($ndt == 6 || $ndt == 7);
 
 		sysread(WBPF,$buf,2*$ndt) == 2*$ndt || die("$WBPcfn: $!");
@@ -1172,6 +1254,20 @@
 		syswrite(WBPF,$buf,1) == 1 || die("$WBPcfn: $!");
 
 		#----------------------------------------------------------------------
+		# Variable Leader #0
+		#	- read ensNo for debugging purposes
+		#----------------------------------------------------------------------
+
+		sysseek(WBPF,$start_ens+$WBPofs[1]+2,0) || die("$WBPcfn: $!");
+		sysread(WBPF,$buf,2) == 2 || die("$WBPcfn: $!");
+		my($ensNo) = unpack("v",$buf);											# only lower two bytes!!!
+		sysseek(WBPF,$start_ens+$WBPofs[1]+13,0) || die("$WBPcfn: $!");			# jump to high byte
+		sysread(WBPF,$buf,1) == 1 || die("$WBPcfn: $!");
+		$ensNo += unpack('C',$buf) << 16;
+		die("ensNo = $ensNo (should be $dta->{ENSEMBLE}[$ens]->{NUMBER})\n")
+			unless ($ensNo == $dta->{ENSEMBLE}[$ens]->{NUMBER});
+
+		#----------------------------------------------------------------------
 		# Variable Leader #1
 		#	- if $RDI_PD0_IO::OVERRIDE_Y2K_CLOCK is set, the data from the pre-Y2K
 		#	  clock are used to override the ADCP clock values; this allows
@@ -1181,7 +1277,7 @@
 
 		if ($RDI_PD0_IO::OVERRIDE_Y2K_CLOCK) {
 			sysseek(WBPF,$start_ens+$WBPofs[1]+4,0) || die("$WBPcfn: $!");		# jump to RTC_YEAR
-			sysread(WBPF,$buf,7) == 7 || die("$WBRcfn: $!");					# read pre-Y2K clock
+			sysread(WBPF,$buf,7) == 7 || die("$WBPcfn: $!");					# read pre-Y2K clock
 			($dta->{ENSEMBLE}[$ens]->{YEAR},
 			 $dta->{ENSEMBLE}[$ens]->{MONTH},
 			 $dta->{ENSEMBLE}[$ens]->{DAY},
@@ -1195,6 +1291,7 @@
 		
 		#----------------------------------------------------------------------
 		# Variable Leader #2
+		#   - read ensemble number for debugging purposes
 		#	- patch everything from SPEED_OF_SOUND to TEMPERATURE
 		# 	- at one stage, IMP allowed for missing values in pitch/roll and heading;
 		#	  on 12/23/2017 the corresponding code was disabled (replaced by assertion)
diff --git a/RDI_Utils.pl b/RDI_Utils.pl
--- a/RDI_Utils.pl
+++ b/RDI_Utils.pl
@@ -1,9 +1,9 @@
 #======================================================================
 #                    R D I _ U T I L S . P L 
 #                    doc: Wed Feb 12 10:21:32 2003
-#                    dlm: Mon Nov 27 10:32:50 2017
+#                    dlm: Sat Jun  9 12:11:01 2018
 #                    (c) 2003 A.M. Thurnherr
-#                    uE-Info: 59 75 NIL 0 0 72 2 2 4 NIL ofnI
+#                    uE-Info: 61 58 NIL 0 0 72 2 2 4 NIL ofnI
 #======================================================================
 
 # miscellaneous RDI-specific utilities
@@ -57,6 +57,8 @@
 #	Aug  7, 2017: - added abmiguity velocity
 #	Aug  8, 2017: - changed transducer frequency to kHz
 #	Nov 27, 2017: - BUG: profile-restart heuristic did not work with P6#001
+#	Mar 18, 2018: - added -ve dt consistency check
+#	Jun  9, 2018: - added support for ENV{NO_GAP_WARNINGS}
 
 use strict;
 
@@ -466,6 +468,10 @@
 	
 		my($dt) = $dta->{ENSEMBLE}[$e]->{UNIX_TIME} -		# time step since
 				  $dta->{ENSEMBLE}[$lastgood]->{UNIX_TIME}; # ... last good ens
+
+		die(sprintf("PANIC: negative dt = $dt between ensembles %d and %d\n",
+			$dta->{ENSEMBLE}[$lastgood]->{NUMBER},$dta->{ENSEMBLE}[$e]->{NUMBER}))
+				if ($dt < 0);
 	
 		if ($dt > $max_gap) {
 			# 2nd heuristic test added Nov 2017 for P06 profile #001
@@ -475,7 +481,8 @@
 #					printf(STDERR "\t[#ens = %d, end-of-gap = $e]\n",scalar(@{$dta->{ENSEMBLE}}));
 					last;
 			}
-			printf(STDERR "WARNING: %.1f-s gap in first half of profile is too long; profile restarted at ensemble $e\n",$dt);
+			printf(STDERR "WARNING: %.1f-s gap beginning at ens#%d in first half of profile is too long; profile restarted at ensemble %d\n",
+				$dt,$dta->{ENSEMBLE}[$lastgood+1]->{NUMBER},$dta->{ENSEMBLE}[$e]->{NUMBER});
 			$firstgood = $lastgood = $e;
 			$dta->{ENSEMBLE}[$e]->{ELAPSED_TIME} = 0;
 			$z = $zErr = $maxz = 0;
@@ -495,7 +502,9 @@
 			$rms_heave_accel_ssq += (($dta->{ENSEMBLE}[$e]->{W}-$dta->{ENSEMBLE}[$lastgood]->{W})/$dt)**2;
 			$rms_heave_accel_nsamp++;
 		} elsif ($dt > 15) {
-	       	printf(STDERR "WARNING: long-ish w gap (dt=%.1fs)\n",$dt);
+	       	printf(STDERR "WARNING: long-ish w gap at ens#%d-%d (dt=%.1fs)\n",
+				$dta->{ENSEMBLE}[$lastgood+1]->{NUMBER},$dta->{ENSEMBLE}[$e-1]->{NUMBER},$dt)
+					unless defined($ENV{NO_GAP_WARNINGS});
 		}
 	
 		$dta->{ENSEMBLE}[$e]->{DEPTH} = $z;
diff --git a/editPD0 b/editPD0
--- a/editPD0
+++ b/editPD0
@@ -2,9 +2,9 @@
 #======================================================================
 #                    E D I T P D 0 
 #                    doc: Mon Nov 25 20:24:31 2013
-#                    dlm: Sat Dec 23 16:10:02 2017
+#                    dlm: Wed Mar 14 21:15:51 2018
 #                    (c) 2013 A.M. Thurnherr
-#                    uE-Info: 83 28 NIL 0 0 72 2 2 4 NIL ofnI
+#                    uE-Info: 417 0 NIL 0 0 72 2 2 4 NIL ofnI
 #======================================================================
 
 # edit RDI PD0 file, e.g. to replace pitch/roll/heading with external values
diff --git a/listBT b/listBT
--- a/listBT
+++ b/listBT
@@ -2,9 +2,9 @@
 #======================================================================
 #                    L I S T B T 
 #                    doc: Sat Jan 18 18:41:49 2003
-#                    dlm: Mon Nov 25 18:30:11 2013
+#                    dlm: Mon Apr  2 18:31:18 2018
 #                    (c) 2003 A.M. Thurnherr
-#                    uE-Info: 42 44 NIL 0 0 72 11 2 4 NIL ofnI
+#                    uE-Info: 43 61 NIL 0 0 72 11 2 4 NIL ofnI
 #======================================================================
 
 # Extract Bottom-Track Data
@@ -40,6 +40,7 @@
 #	Nov  1, 2008: - BUG: sig(u) was reported instead of sig(v)
 #	Jul 30, 2009: - NaN => nan
 #	Nov 25, 2013: - checkEnsemble() expunged
+#	Apr  2, 2018: - BUG: velBeamToInstrument() used old usage
 
 # NOTES:
 #	- the RDI BT data contains ranges that are greater than the
@@ -164,7 +165,7 @@
 					 $dta{ENSEMBLE}[$ens]->{PERCENT_GOOD}[$b][2] < 100 ||
 					 $dta{ENSEMBLE}[$ens]->{PERCENT_GOOD}[$b][3] < 100);
 			@v = velInstrumentToEarth(\%dta,$ens,
-					velBeamToInstrument(\%dta,
+					velBeamToInstrument(\%dta,$ens,
 	                    @{$dta{ENSEMBLE}[$ens]->{VELOCITY}[$b]}));
 	    } else {
 			next if ($dta{ENSEMBLE}[$ens]->{PERCENT_GOOD}[$b][0] > 0 ||
@@ -298,7 +299,7 @@
 		last if ($dta{ENSEMBLE}[$e]->{BT_MEAN_RANGE}-$dz <= $slc);
 		my(@v) = $beamCoords
 			   ? velInstrumentToEarth(\%dta,$e,
-					velBeamToInstrument(\%dta,
+					velBeamToInstrument(\%dta,$e,
 						@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]}))
 			   : velApplyHdgBias(\%dta,$e,@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]});
 		next unless defined($v[0]);
@@ -338,7 +339,7 @@
 		last if ($dta{ENSEMBLE}[$e]->{BT_MEAN_RANGE}-$dz <= $slc);
 		my(@v) = $beamCoords
 			   ? velInstrumentToEarth(\%dta,$e,
-					velBeamToInstrument(\%dta,
+					velBeamToInstrument(\%dta,$e,
 						@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]}))
 			   : velApplyHdgBias(\%dta,$e,@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]});
 		next unless defined($v[0]);
@@ -484,7 +485,7 @@
 
 	@{$dta{ENSEMBLE}[$e]->{BT_VELOCITY}} = $beamCoords	# xform BT vel
 	    ? velInstrumentToEarth(\%dta,$e,
-			velBeamToInstrument(\%dta,@{$dta{ENSEMBLE}[$e]->{BT_VELOCITY}}))
+			velBeamToInstrument(\%dta,$e,@{$dta{ENSEMBLE}[$e]->{BT_VELOCITY}}))
 		: velApplyHdgBias(\%dta,$e,@{$dta{ENSEMBLE}[$e]->{BT_VELOCITY}});
 
  	$dta{ENSEMBLE}[$e]->{BT_MEAN_RANGE} =				# mean vals
@@ -555,7 +556,7 @@
 		my($btm_e) = int($btm_e[0]/4+$btm_e[1]/4+$btm_e[2]/4+$btm_e[3]/4+0.5);
 		@{$dta{ENSEMBLE}[$btm_e]->{BTFWT_VELOCITY}} = $beamCoords # BT from WT
 			? velInstrumentToEarth(\%dta,$e,
-				velBeamToInstrument(\%dta,@{$dta{ENSEMBLE}[$e]->{BT_VELOCITY}}))
+				velBeamToInstrument(\%dta,$e,@{$dta{ENSEMBLE}[$e]->{BT_VELOCITY}}))
 	        : velApplyHdgBias(\%dta,$e,@{$dta{ENSEMBLE}[$e]->{BT_VELOCITY}});
 	}
 
diff --git a/listBins b/listBins
--- a/listBins
+++ b/listBins
@@ -66,6 +66,10 @@
 #				  - BUG: 3-beam %ages were incorrect: 1) they were based on goodvels instead of the
 #						 entire ensemble range; 2) pings-per-ensemble were not considered
 #				  - BUG: output layout was all messed up for non-valid velocities
+#	Feb  6, 2018: - added support for PD0_IO first_ens, last_ens
+#	Apr 10, 2018: - added day number to output
+#				  - added -l)ast bin
+#				  - activate output files
 #	Aug 29, 2018: - added error message on -r decoding failures
 
 # General Notes:
@@ -102,20 +106,26 @@
 
 use Getopt::Std;
 
-$ADCP_tools_minVersion = 2.1;
+$ADCP_tools_minVersion = 2.2;
 ($ADCP_TOOLS) = ($0 =~ m{(.*/)[^/]+});
 require "$ADCP_TOOLS/ADCP_tools_lib.pl";
 
-die("Usage: $0 [-r)ange <first_ens,last_ens>] [-R)enumber ensembles] " .
+$antsMinLibVersion = 7.0;
+($ANTS) = (`which ANTSlib` =~ m{^(.*)/[^/]*$});
+require "$ANTS/ants.pl";
+require "$ANTS/libconv.pl";
+
+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>] " .
 			  "[-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:M:o:p:r:P:RS:") && @ARGV == 1);
+	unless (&getopts("4aB:d:l:M:o:p:r:P:RS:T:") && @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);
@@ -146,6 +156,9 @@
 	print(STDERR "WARNING: no soundspeed correction applied!\n");
 }
 
+loadInstrumentTransformation($opt_T)				# load instrument-transformation matrix
+	if (defined($opt_T));
+
 #----------------------------------------------------------------------
 
 sub min(@)								# return minimum
@@ -163,6 +176,8 @@
 	my($out) = sprintf($opt_o,$b+1);
 
 	open(P,"$out") || die("$out: $!\n");
+	print(P "#!/usr/bin/perl -S list\n");
+	chmod(0777&~umask,*P);																	# activate output
 	print(P "#ANTS#PARAMS# ");
 	foreach my $k (keys(%P)) {
 		print(P "$k\{$P{$k}\} ");
@@ -182,7 +197,7 @@
 	print( P "\n");
 
 	print(P "#ANTS#FIELDS# " .
-			"{ensemble} {date} {time} {elapsed} " .
+			"{ensemble} {date} {time} {elapsed} {dn} " .
 			"{heading} {pitch} {roll} " .
 			"{sig_heading} {sig_pitch} {sig_roll} " .
 			"{xmit_current} {xmit_voltage} " .
@@ -205,7 +220,8 @@
 		print(P "$dta{ENSEMBLE}[$e]->{NUMBER} ");
 		print(P "$dta{ENSEMBLE}[$e]->{DATE} ");
 		print(P "$dta{ENSEMBLE}[$e]->{TIME} ");
-		printf(P "%d ",$dta{ENSEMBLE}[$e]->{UNIX_TIME}-$t0);		# elapsed time
+		printf(P "%d ",$dta{ENSEMBLE}[$e]->{UNIX_TIME}-$t0);												# elapsed time
+		printf(P "%g ",str2dec_time($dta{ENSEMBLE}[$e]->{DATE},$dta{ENSEMBLE}[$e]->{TIME}));				# decimal day 
 		print(P defined($dta{ENSEMBLE}[$e]->{HEADING}) ? "$dta{ENSEMBLE}[$e]->{HEADING} " : 'nan ');
 		print(P defined($dta{ENSEMBLE}[$e]->{PITCH}) ? "$dta{ENSEMBLE}[$e]->{PITCH} " : 'nan ');
 		print(P defined($dta{ENSEMBLE}[$e]->{ROLL}) ? "$dta{ENSEMBLE}[$e]->{ROLL} " : 'nan ');
@@ -269,7 +285,7 @@
 $P{RDI_file} = $ifn;
 $P{mag_decl} = $opt_M if defined($opt_M);
 
-readData($ifn,\%dta);
+readData($ifn,\%dta,$first_ens,$last_ens,$opt_l);
 printf(STDERR "%d complete ensembles...\n",scalar(@{$dta{ENSEMBLE}}));
 $dta{HEADING_BIAS} = -$opt_M;						# magnetic declination
 
diff --git a/listEns b/listEns
--- a/listEns
+++ b/listEns
@@ -2,13 +2,12 @@
 #======================================================================
 #                    L I S T E N S 
 #                    doc: Sat Jan 18 18:41:49 2003
-#                    dlm: Wed Nov  9 12:27:32 2016
+#                    dlm: Thu May 31 11:10:51 2018
 #                    (c) 2003 A.M. Thurnherr
-#                    uE-Info: 115 0 NIL 0 0 72 2 2 4 NIL ofnI
+#                    uE-Info: 62 51 NIL 0 0 72 2 2 4 NIL ofnI
 #======================================================================
 
-# Print useful info from the ensemble list or dump ensembles to
-# separate files.
+$synopsis = 'list ensemble summaries (default), dump ensembles (-E), time-average ensembles (-T)';
 
 # HISTORY:
 #	Jan 18, 2003: - created
@@ -48,40 +47,80 @@
 #   Mar 17, 2016: - adapted to new Getopt library
 #	Apr 19, 2016: - added %date, %time to -E output
 #	Nov  9, 2016: - BUG: no error on missing files
+#	Feb  7, 2018: - removed 3-beam error
+#	Apr  1, 2018: - improved usage message
+#				  - removed -Q option (errcheck only, which is not necessary; can use mkProfile to check for errors)
+#				  - added -T (time averaging)
+#	Apr  2, 2018: - made it work
+#				  - BUG: velBeamToInstrument() was using old usage
+#	Apr  3, 2018: - BUG: typo
+#				  - added -S from [listBins]
+#				  - removed -B and an BT data (current version did not treat beam-coord BT data correctly)
+#	Apr  4, 2018: - added support for first_ens and last_ens in [RDI_PD0_IO.pl]
+#				  - removed support for multiple input files
+#	Apr 10, 2018: - added -l)ast bin
+#	May 31, 2018: - BUG: -A was disabled by default
 
 # Notes:
-#	- -E outputs data in earth coordinates, unless -b is set also
-#	- -E output is always in ANTS format, ignoring -A
-#	- no soundspeed correction
+#	- -E/-B outputs data in earth coordinates, unless -b is set also
+#	- -E/-T output is always in ANTS format
 
 use Getopt::Std;
-$0 =~ m{(.*/)[^/]+};
-require "$1RDI_BB_Read.pl";
-require "$1RDI_Coords.pl";
+
+$ADCP_tools_minVersion = 2.2;
+($ADCP_TOOLS) = ($0 =~ m{(.*/)[^/]+});
+require "$ADCP_TOOLS/ADCP_tools_lib.pl";
+
+$antsMinLibVersion = 7.0;
+($ANTS) = (`which ANTSlib` =~ m{^(.*)/[^/]*$});
+require "$ANTS/ants.pl";
+require "$ANTS/libconv.pl";
+
+($self) = ($0 =~ m{.*/([^/]+)});
+$cmdline = "$self @ARGV";
 
-die("Usage: $0 [-A)nts] [-Q)uiet (errcheck only)] " .
-			  "[-f)ields <[name=]FIELD[,...]>] " .
-			  "[require -4)-beam solutions] [-d)iscard <beam#>] " .
-			  "[write -E)nsemples <.suff> [use -B)T] [-M)agnetic <declination>] [min -p)ercent-good <#>] [keep -b)eam coords]] " .
-			  "[-r)ange <first_ens,last_ens>] [in-w)ater ensembles only] " .
-			  "<RDI file...>\n")
-	unless (&getopts("4ABbd:E:f:M:p:Qr:w") && $#ARGV >= 0);
+die("\n$self [-4ABbdEfiMprSTw] <PD0 file> -- $synopsis\n\nCommand-Line Options:\n\t" .
+			  "Output Ensemble Summary (default mode):\n\t\t" .
+				  "[-A)NTS format]\n\t\t" .
+			  "Dump Ensembles (-E|-T; ANTS Format):\n\t\t" .
+				  "[dump individual -E)nsemples <.suff>]\n\t\t" .
+				  "[-T)ime-average ensembles [time-series start (decimal day),]<averaging interval (days)>[,time-series end (decimal day)]]\n\t\t\t" .
+				  	"[-i)gnore bins with <fraction> of max samples (-T only)]\n\t\t\t" .
+				  	"[output -B)asename <bn>]\n\t\t" .
+				  "[-M)agnetic <declination>]\n\t\t" .
+				  "[-S)oundspeed correction <salin|*,temp|*,depth|*>\n\t\t" .
+				  "[require min -p)ercent-good <#>]\n\t\t" .
+				  "[keep -b)eam coords (do not transform to earth coordinates)]\n\t" .
+			  "Common Options:\n\t\t" .
+				  "[add -f)ields <[name=]FIELD[,...]>]\n\t\t" .
+				  "[require -4)-beam solutions] [-d)iscard <beam#>]\n\t\t" .
+				  "[-r)ange <first_ens,last_ens>] [-l)ast <bin>]\n\t\t" .
+				  "[in-w)ater ensembles only]\n")
+	unless (&getopts("4AB:bd:E:f:i:l:M:p:r:S:T:w") && @ARGV == 1);
 
-print(STDERR "WARNING: no soundspeed correction applied!\n");
+die("$ARGV[0]: no such file\n")
+	unless (-f $ARGV[0]);
+
+$dump_ens = defined($opt_E) + defined($opt_T);
+die("$self: cannot combine -E with -T\n") if ($dump_ens > 1);
+
+if (defined($opt_S)) {
+	($SS_salin,$SS_temp,$SS_depth) = split(',',$opt_S);
+} else {
+	print(STDERR "WARNING: no soundspeed correction applied!\n")
+		if ($dump_ens);
+}
 
 print(STDERR "WARNING: magnetic declination not set!\n")
-	if defined($opt_E) && !defined($opt_M);
+	if ($dump_ens && !defined($opt_M));
 
-die("$0: illegal option combination\n")
-	if ($opt_Q && $opt_A) || ((defined($opt_M) || defined($opt_p) || defined($opt_b) || $opt_B) && !defined($opt_E));
+die("$self: illegal option combination\n")
+	if ((defined($opt_M) || defined($opt_p) || defined($opt_b)) && !defined($dump_ens));
 
-die("$0: -4 and -d are mutually exclusive\n")
+die("$self: -4 and -d are mutually exclusive\n")
 	if ($opt_4 && defined($opt_d));
 
-($first_ens,$last_ens) = split(',',$opt_r)
-	if defined($opt_r);
-
-undef($opt_A) if defined($opt_E);
+#undef($opt_A) if defined($dump_ens);
 
 $opt_p = 0 unless defined($opt_p);
 
@@ -107,239 +146,406 @@
 # MAIN
 #----------------------------------------------------------------------
 
-while ($ARGV[0] ne '') {
-	die("$ARGV[0]: No such file or directory\n")
-		unless (-f $ARGV[0]);
+printf(STDERR "Reading $ARGV[0]...");
+if (defined($opt_r)) {								# read selected range
+	my($fe,$le) = split(',',$opt_r);
+	readData(@ARGV,\%dta,$fe,$le,$opt_l);
+} else {											# read entire file (possibly selected bins)
+	readData(@ARGV,\%dta,undef,undef,$opt_l);
+}
+printf(STDERR "\n\t%d complete ensembles\n",scalar(@{$dta{ENSEMBLE}}));
+
+$dta{HEADING_BIAS} = -$opt_M;						# magnetic declination
+
+if ($dta{BEAM_COORDINATES}) {						# coords used
+	$beamCoords = 1;
+} elsif (!$dta{EARTH_COORDINATES}) {
+	die("$ARGV[0]: beam or earth coordinates required (implementation restriction)\n");
+}
+die("$ARGV[0]: -b requires beam-coordinate data\n")
+	if ($opt_b && !$beamCoords);
+
+($basename) = defined($opt_B)						# set basename of output files
+			? $opt_B
+			: ($ARGV[0] =~ m{([^\./]+)\.[^\.]+});
+
+#----------------------------------------------------------------------
+# define &dumpEns() routine for different output formats:
+#	-A			ANTS format ensemble summary
+#	-E			create one file per ensemble
+#	-T			time-average multiple ensembles
+#	default:	ASCII ensemble summary (human readable)
+#----------------------------------------------------------------------
 
-	readData(@ARGV,\%dta);
-	if ($opt_A && !$opt_E) {
-		print("#ANTS#PARAMS# RDI_file{$ARGV[0]}\n");
-	} elsif (!$opt_Q) {
-		print(STDERR "$ARGV[0]: ");
-	}
-	printf(STDERR "%d complete ensembles...\n",scalar(@{$dta{ENSEMBLE}}))
-		unless ($opt_Q);
-	$dta{HEADING_BIAS} = -$opt_M;						# magnetic declination
-	shift;
+if ($opt_A) {										# select output fmt: ANTS
+	print("#ANTS#PARAMS# PD0_file{$ARGV[0]}\n");
+	printf("#ANTS#PARAMS# N_ensembles{%d}\n",scalar(@{$dta{ENSEMBLE}}));
+	print('#ANTS#FIELDS# {ens} {date} {time} {unix-time} {xducer_up} {temp} {hdg} {pitch} {roll} {XMIT_VOLTAGE} {XMIT_CURRENT}');
+	print(' {ESW}') if ($dta{FIXED_LEADER_BYTES} >= 53);
+	print("$addLayout\n");
 
-	if ($dta{BEAM_COORDINATES}) {						# coords used
-		$beamCoords = 1;
-	} elsif (!$dta{EARTH_COORDINATES}) {
-		die("$ARGV[0]: only beam and earth coordinates implemented so far\n");
+	$dumpEns = sub ($)
+	{
+		my($e) = @_;
+	    
+		printf('%d %s %s %lf %d %g',
+			$dta{ENSEMBLE}[$e]->{NUMBER},
+			$dta{ENSEMBLE}[$e]->{DATE},
+			$dta{ENSEMBLE}[$e]->{TIME},
+			$dta{ENSEMBLE}[$e]->{UNIX_TIME},
+			$dta{ENSEMBLE}[$e]->{XDUCER_FACING_UP} ? 1 : 0,
+			$dta{ENSEMBLE}[$e]->{TEMPERATURE},
+		);
+		if (defined($dta{ENSEMBLE}[$e]->{HEADING})) { printf(' %g',$dta{ENSEMBLE}[$e]->{HEADING}); }
+		else { printf(' nan'); }
+		if (defined($dta{ENSEMBLE}[$e]->{PITCH})) { printf(' %g',$dta{ENSEMBLE}[$e]->{PITCH}); }
+		else { printf(' nan'); }
+		if (defined($dta{ENSEMBLE}[$e]->{ROLL})) { printf(' %g',$dta{ENSEMBLE}[$e]->{ROLL}); }
+		else { printf(' nan'); }
+		printf(' %g %g',
+			$dta{ENSEMBLE}[$e]->{ADC_XMIT_VOLTAGE},
+			$dta{ENSEMBLE}[$e]->{ADC_XMIT_CURRENT},
+		);
+		printf(' %08X',$dta{ENSEMBLE}[$e]->{ERROR_STATUS_WORD})
+			if ($dta{FIXED_LEADER_BYTES} >= 53);
+		foreach my $f (@addFields) {
+			my($fn,$fi) = ($f =~ m{([^[]*)(\[.*)});
+			$fn = $f unless defined($fn);
+			my($v) = eval("\$dta{ENSEMBLE}[$e]->{$fn}$fi");
+			print(defined($v) ? " $v" : " nan");
+		}
+		print("\n");
 	}
 
-	die("$ARGV[0]: -b only makes sense for beam-coordinate data\n")
-		if ($opt_b && !$beamCoords);
+} elsif ($opt_E) {										# dump each ensemble in separate file 
 
-	if ($opt_A) {										# select output fmt: ANTS
-		unless ($opt_Q) {
-			printf("#ANTS#PARAMS# N_ensembles{%d}\n",scalar(@{$dta{ENSEMBLE}}));
-			print('#ANTS#FIELDS# {ens} {time} {xducer_up} {temp} {hdg} {pitch} {roll} {XMIT_VOLTAGE} {XMIT_CURRENT}');
-			print(' {ESW}') if ($dta{FIXED_LEADER_BYTES} >= 53);
-	        print("$addLayout\n");
-	    }
+	$dumpEns = sub ($)
+	{
+		my($e) = @_;
+		my($b,$i);
+		my($file) = "$dta{ENSEMBLE}[$e]->{NUMBER}$opt_E";
+    
+		my($ssCorr) = defined($opt_S) ? ssCorr($dta{ENSEMBLE}[$e],$SS_salin,$SS_temp,$SS_depth) : 1;
 
-		$dumpEns = sub ($)
-		{
-			my($e) = @_;
-			
-			printf('%d %lf %d %g',
-				$dta{ENSEMBLE}[$e]->{NUMBER},
-				$dta{ENSEMBLE}[$e]->{UNIX_TIME},
-				$dta{ENSEMBLE}[$e]->{XDUCER_FACING_UP} ? 1 : 0,
-				$dta{ENSEMBLE}[$e]->{TEMPERATURE},
+		open(P,">$file") || die("$file: $!\n");
+		print(P "#!/usr/bin/perl -S list\n");
+		print(P "#ANTS#PARAMS# " .
+				"date{$dta{ENSEMBLE}[$e]->{DATE}} " .
+				"time{$dta{ENSEMBLE}[$e]->{TIME}} " .
+				"soundspeed_correction{%s} " .
+				"magnetic_declination{%g} " .
+				"\n",
+					(defined($opt_S) ? $opt_S : "NONE!"),
+					$opt_M
+		);
+		print(P "#ANTS#FIELDS# " .
+				"{bin} {dz} {u} {v} {w} {e} {w12} {w34} {corr1} {corr2} {corr3} {corr4} " .
+				"{amp1} {amp2} {amp3} {amp4} "
+		);
+		if ($beamCoords) {
+			print(P "{pcg1} {pcg2} {pcg3} {pcg4}");
+		} else {
+			print(P "{pc3beam} {pcBadErrVel} {pc1or2beam} {pc4beam}");
+		}
+		print(P "$addLayout\n");
+		    
+		my($ssCorr) = defined($opt_S) ? ssCorr($dta{ENSEMBLE}[$e],$SS_salin,$SS_temp,$SS_depth) : 1;
+		for (my($b)=0; $b<$dta{N_BINS}; $b++) {
+			my(@v,$w12,$w34);
+			my($dz) = $ssCorr * ($dta{DISTANCE_TO_BIN1_CENTER} + $b*$dta{BIN_LENGTH});
+    
+			if ($beamCoords) {
+				undef($dta{ENSEMBLE}[$e]->{VELOCITY}[$b][0])
+					if (($dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][0] < $opt_p) || ($opt_d == 1));
+				undef($dta{ENSEMBLE}[$e]->{VELOCITY}[$b][1])
+					if (($dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][1] < $opt_p) || ($opt_d == 2));
+				undef($dta{ENSEMBLE}[$e]->{VELOCITY}[$b][2])
+					if (($dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][2] < $opt_p) || ($opt_d == 3));
+				undef($dta{ENSEMBLE}[$e]->{VELOCITY}[$b][3])
+					if (($dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][3] < $opt_p) || ($opt_d == 4));
+				($dummy,$w12,$dummy,$w34) =
+					velBeamToBPEarth(\%dta,$e,@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]});
+				@v = $opt_b ? @{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]} :
+							  velBeamToEarth(\%dta,$e,@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]});
+			} else {
+				@v = velApplyHdgBias(\%dta,$e,@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]});
+			}
+			$v[0] *= $ssCorr if defined($v[0]);
+			$v[1] *= $ssCorr if defined($v[1]);
+			$v[2] *= $ssCorr if defined($v[2]);
+			$v[3] *= $ssCorr if defined($v[3]);
+			$w12 *= $ssCorr if defined($w12);
+			$w34 *= $ssCorr if defined($w34);
+    
+			$v[0] = nan unless defined($v[0]);		# u
+			$v[1] = nan unless defined($v[1]);		# v
+			$v[2] = nan unless defined($v[2]);		# w
+			$v[3] = nan unless defined($v[3]);		# err_vel
+			$w12 = nan unless defined($w12);		# w from beams 1&2
+			$w34 = nan unless defined($w34);		# w from beams 3&4
+
+			my(@out) = (
+				$b+1,$dz,$v[0],$v[1],$v[2],$v[3],$w12,$w34,
+				@{$dta{ENSEMBLE}[$e]->{CORRELATION}[$b]},
+				@{$dta{ENSEMBLE}[$e]->{ECHO_AMPLITUDE}[$b]},
+				@{$dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b]}
 			);
-			if (defined($dta{ENSEMBLE}[$e]->{HEADING})) { printf(' %g',$dta{ENSEMBLE}[$e]->{HEADING}); }
-			else { printf(' nan'); }
-			if (defined($dta{ENSEMBLE}[$e]->{PITCH})) { printf(' %g',$dta{ENSEMBLE}[$e]->{PITCH}); }
-			else { printf(' nan'); }
-			if (defined($dta{ENSEMBLE}[$e]->{ROLL})) { printf(' %g',$dta{ENSEMBLE}[$e]->{ROLL}); }
-			else { printf(' nan'); }
-			printf(' %g %g',
-				$dta{ENSEMBLE}[$e]->{ADC_XMIT_VOLTAGE},
-				$dta{ENSEMBLE}[$e]->{ADC_XMIT_CURRENT},
-			);
-			printf(' %08X',$dta{ENSEMBLE}[$e]->{ERROR_STATUS_WORD})
-				if ($dta{FIXED_LEADER_BYTES} >= 53);
 			foreach my $f (@addFields) {
 				my($fn,$fi) = ($f =~ m{([^[]*)(\[.*)});
 				$fn = $f unless defined($fn);
-				my($v) = eval("\$dta{ENSEMBLE}[$e]->{$fn}$fi");
-				print(defined($v) ? " $v" : " nan");
+				push(@out,eval("\$dta{ENSEMBLE}[$e]->{$fn}$fi"));
+			}
+			for ($i=0; $i<19+@addFields; $i++) {
+				$out[$i] = nan unless defined($out[$i]);
 			}
-			print("\n");
+			print(P "@out\n");
 		}
+		chmod(0777&~umask,*P);
+		close(P);
+	}
 
-	} elsif ($opt_E) {										# one file per ens
+} elsif (defined($opt_T)) { 									# time-average ensembles
+
+	my(@tmp) = split(',',$opt_T);								# decode -T 
+	my($Tstart,$deltaT,$Tend,$month,$day,$year);
+	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);
+	} else {
+		($month,$day,$year) = split('/',$dta{ENSEMBLE}[0]->{DATE});
+		$Tstart = str2dec_time($dta{ENSEMBLE}[0]->{DATE},$dta{ENSEMBLE}[0]->{TIME},$year);
+		($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);
+	}
+	$Tstart = &{&antsCompileConstExpr($')} if ($Tstart =~ m{^=});
+	$deltaT = &{&antsCompileConstExpr($')} if ($deltaT =~ m{^=});
+	$Tend	= &{&antsCompileConstExpr($')} if ($Tend   =~ m{^=});
+	$deltaT = 9e99 unless ($deltaT > 0);
+																# format string for zero-padded filenames
+	my($fnfmt) = sprintf('%%0%dd',length(sprintf('%d',($Tend-$Tstart)/$deltaT+1)));
+
+	$cbin		= 1;											# current tim-bin number; used inside dumpEns
+	$max_nens	= 0;											# max number of samples
+
+	sub dde($$) 												# divide allowing for zero; used inside dumpEns
+	{
+		my($sum,$n) = @_;
+		return $n ? $sum/$n : nan;
+	}
 
-		$dumpEns = sub ($)
-		{
-			my($e) = @_;
-			my($b,$i);
-			my($file) = "$dta{ENSEMBLE}[$e]->{NUMBER}$opt_E";
-		
-			open(P,">$file") || die("$file: $!\n");
-			print(P "#!/usr/bin/perl -S list\n");
-			print(P "#ANTS#PARAMS# " .
-					"date{$dta{ENSEMBLE}[$e]->{DATE}} " .
-					"time{$dta{ENSEMBLE}[$e]->{TIME}} " .
-					"BT_u{$dta{ENSEMBLE}[$e]->{BT_VELOCITY}[0]} " .
-					"BT_v{$dta{ENSEMBLE}[$e]->{BT_VELOCITY}[1]} " .
-					"BT_w{$dta{ENSEMBLE}[$e]->{BT_VELOCITY}[2]} " .
-					"BT_e{$dta{ENSEMBLE}[$e]->{BT_VELOCITY}[3]} " .
-					"BT_cor1{$dta{ENSEMBLE}[$e]->{BT_CORRELATION}[0]} " .
-					"BT_cor2{$dta{ENSEMBLE}[$e]->{BT_CORRELATION}[1]} " .
-					"BT_cor3{$dta{ENSEMBLE}[$e]->{BT_CORRELATION}[2]} " .
-					"BT_cor4{$dta{ENSEMBLE}[$e]->{BT_CORRELATION}[3]} " .
-					"BT_amp1{$dta{ENSEMBLE}[$e]->{BT_AMPLITUDE}[0]} " .
-					"BT_amp2{$dta{ENSEMBLE}[$e]->{BT_AMPLITUDE}[1]} " .
-					"BT_amp3{$dta{ENSEMBLE}[$e]->{BT_AMPLITUDE}[2]} " .
-					"BT_amp4{$dta{ENSEMBLE}[$e]->{BT_AMPLITUDE}[3]} " .
-					"soundspeed_correction{NONE!} " .
-					"\n"
-			);
-			print(P "#ANTS#FIELDS# " .
-					"{bin} {dz} {u} {v} {w} {e} {w12} {w34} {cor1} {cor2} {cor3} {cor4} " .
-					"{amp1} {amp2} {amp3} {amp4} "
-			);
-			if ($beamCoords) {
-				print(P "{pcg1} {pcg2} {pcg3} {pcg4}");
-			} else {
-				print(P "{pc3beam} {pcBadErrVel} {pc1or2beam} {pc4beam}");
-			}
-			print(P "$addLayout\n");
+	$dumpEns = sub ($)											# time average and output when next bin is started
+	{
+		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;
+
+		if ($e<0 || $dn>=$Tstart+$cbin*$deltaT) {				# dump full bin
+			my($file) = sprintf("$basename.T$fnfmt",$cbin); 	# file name: <basename>.T0001
+
+			do {												# skip empy bins
+				$cbin++;
+			} while ($dn>=$Tstart+$cbin*$deltaT);
+
+			$max_nens = $a1_n[0] if ($a1_n[0] > $max_nens); 	# update max number of samples in time bin
+			if ($a1_n[0] >= $opt_i*$max_nens) { 				# write file only if sufficient samples (-i)
 			    
-			for (my($b)=0; $b<$dta{N_BINS}; $b++) {
-				my(@v,$w12,$w34);
-				my($dz) = $dta{DISTANCE_TO_BIN1_CENTER} + $b*$dta{BIN_LENGTH};
-		
+				open(P,">$file") || die("$file: $!\n"); 		# open file and output metadata
+				print(P "#!/usr/bin/perl -S list\n");
+				print(P "#ANTS# $cmdline\n");
+				printf(P "#ANTS#PARAMS# " .
+						"PD0_file{$ARGV[0]} " .
+						"dn{%s} " .
+						"N_ensembles{%d} ensemble_range{%d,%d} " .
+						"delta-T{%g} " .
+						"soundspeed_correction{%s} " .
+						"magnetic_declination{%g} " .
+						"\n",
+							&dde($dn_s,$dn_n),
+							$a1_n[0],$feib,$leib,$deltaT,
+							(defined($opt_S) ? $opt_S : "NONE!"),
+							$opt_M
+
+				);
+			    
+				print(P "#ANTS#FIELDS# " .													# Layout
+						"{bin} {dz} {u} {v} {w} {e} {w12} {w34} {corr1} {corr2} {corr3} {corr4} " .
+						"{amp1} {amp2} {amp3} {amp4} "
+				);
 				if ($beamCoords) {
-					undef($dta{ENSEMBLE}[$e]->{VELOCITY}[$b][0])
-						if (($dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][0] < $opt_p) || ($opt_d == 1));
-					undef($dta{ENSEMBLE}[$e]->{VELOCITY}[$b][1])
-						if (($dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][1] < $opt_p) || ($opt_d == 2));
-					undef($dta{ENSEMBLE}[$e]->{VELOCITY}[$b][2])
-						if (($dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][2] < $opt_p) || ($opt_d == 3));
-					undef($dta{ENSEMBLE}[$e]->{VELOCITY}[$b][3])
-						if (($dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][3] < $opt_p) || ($opt_d == 4));
-					($dummy,$w12,$dummy,$w34) =
-						velBeamToBPEarth(\%dta,$e,@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]});
-					@v = $opt_b ? @{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]} :
-							velInstrumentToEarth(\%dta,$e,
-								velBeamToInstrument(\%dta,
-									@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]}));
+					print(P "{pcg1} {pcg2} {pcg3} {pcg4} ");
 				} else {
-					@v = velApplyHdgBias(\%dta,$e,@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]});
+					print(P "{pc3beam} {pcBadErrVel} {pc1or2beam} {pc4beam} ");
 				}
-		
-				$v[0] = nan unless defined($v[0]);		# u
-				$v[1] = nan unless defined($v[1]);		# v
-				$v[2] = nan unless defined($v[2]);		# w
-				$v[3] = nan unless defined($v[3]);		# err_vel
-				$w12 = nan unless defined($w12);		# w from beams 1&2
-				$w34 = nan unless defined($w34);		# w from beams 3&4
+				print(P "{uvw.nsamp} {e.nsamp} ");
+				print(P "$addLayout\n");
+																							# ssCorr for dz based on first ensemble in bin
+				my($ssCorr) = defined($opt_S)
+							? ssCorr($dta{ENSEMBLE}[$feib-$dta{ENSEMBLE}[0]->{NUMBER}],$SS_salin,$SS_temp,$SS_depth)
+							: 1;
+				for (my($b)=0; $b<$dta{N_BINS}; $b++) { 									# output data
+					my($dz) = $ssCorr * ($dta{DISTANCE_TO_BIN1_CENTER} + $b*$dta{BIN_LENGTH});
+					printf(P "%d %g  %g %g %g %g  %g %g  %d %d %d %d  %d %d %d %d  %d %d %d %d	%d %d",
+						$b+1,$dz,
+						&dde($v1[$b],$v1_n[$b]),&dde($v2[$b],$v2_n[$b]),&dde($v3[$b],$v3_n[$b]),&dde($v4[$b],$v4_n[$b]),
+						&dde($w12[$b],$w12_n[$b]),&dde($w34[$b],$w34_n[$b]),
+						&dde($c1[$b],$c1_n[$b]),&dde($c2[$b],$c2_n[$b]),&dde($c3[$b],$c3_n[$b]),&dde($c4[$b],$c4_n[$b]),
+						&dde($a1[$b],$a1_n[$b]),&dde($a2[$b],$a2_n[$b]),&dde($a3[$b],$a3_n[$b]),&dde($a4[$b],$a4_n[$b]),
+						&dde($p1[$b],$p1_n[$b]),&dde($p2[$b],$p2_n[$b]),&dde($p3[$b],$p3_n[$b]),&dde($p4[$b],$p4_n[$b]),
+						$v1_n[$b],$v4_n[$b]
+					);
 
-				my(@out);
-				if ($opt_B) {
-					$v[0] = nan unless defined($dta{ENSEMBLE}[$e]->{BT_VELOCITY}[0]);
-					$v[1] = nan unless defined($dta{ENSEMBLE}[$e]->{BT_VELOCITY}[1]);
-					@out = (
-						$b,$dz,$v[0]-$dta{ENSEMBLE}[$e]->{BT_VELOCITY}[0],
-						$v[1]-$dta{ENSEMBLE}[$e]->{BT_VELOCITY}[1],$v[2],$v[3],$w12,$w34,
-						@{$dta{ENSEMBLE}[$e]->{CORRELATION}[$b]},
-						@{$dta{ENSEMBLE}[$e]->{ECHO_AMPLITUDE}[$b]},
-						@{$dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b]}
-					);
-				} else {
-					@out = (
-						$b,$dz,$v[0],$v[1],$v[2],$v[3],$w12,$w34,
-						@{$dta{ENSEMBLE}[$e]->{CORRELATION}[$b]},
-						@{$dta{ENSEMBLE}[$e]->{ECHO_AMPLITUDE}[$b]},
-						@{$dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b]}
-					);
+					for (my($i)=0; $i<@af; $i++) {
+						printf(P "%g ",&dde($af[$i][$b],$af_n[$i][$b]));
+					}
+					print(P "\n");
 				}
-				foreach my $f (@addFields) {
-					my($fn,$fi) = ($f =~ m{([^[]*)(\[.*)});
-					$fn = $f unless defined($fn);
-					push(@out,eval("\$dta{ENSEMBLE}[$e]->{$fn}$fi"));
+				chmod(0777&~umask,*P);														# activate output
+				close(P);
+			} # if -i check okay
+
+			for (my($b)=0; $b<$dta{N_BINS}; $b++) { 										# reset stats	 
+				$v1[$b] = $v1_n[$b] = $v2[$b] = $v2_n[$b] = $v3[$b] = $v3_n[$b] = $v4[$b] = $v4_n[$b] = 0;
+				$w12[$b] = $w12_n[$b] = $w34[$b] = $w34_n[$b] = 0;
+				$c1[$b] = $c1_n[$b] = $c2[$b] = $c2_n[$b] = $c3[$b] = $c3_n[$b] = $c4[$b] = $c4_n[$b] = 0;							    
+				$a1[$b] = $a1_n[$b] = $a2[$b] = $a2_n[$b] = $a3[$b] = $a3_n[$b] = $a4[$b] = $a4_n[$b] = 0;
+				$p1[$b] = $p1_n[$b] = $p2[$b] = $p2_n[$b] = $p3[$b] = $p3_n[$b] = $p4[$b] = $p4_n[$b] = 0;
+				for (my($i)=0; $i<@af; $i++) {
+					$af[$i][$b] = $af_n[$i][$b] = 0;
 				}
-				for ($i=0; $i<19+@addFields; $i++) {
-					$out[$i] = nan unless defined($out[$i]);
-				}
-				print(P "@out\n");
+				$dn_s = $dn_n = 0;															# day number
 			}
-			close(P);
-		}
+				    
+			undef($feib);																	# make sure first ensemble in bin will be updated below
+		} # if time bin is full
 
-	} else {											# neither ANTS nor ens files
-		unless ($opt_Q) {
-			if ($dta{FIXED_LEADER_BYTES} >= 53) {
-				printf("    # Date       Time        XD  Temp Headng Pitch  Roll #vv DSID ESW$addLayout\n");
-				printf("-----------------------------------------------------------------------\n");
-			} else {
-				printf("    # Date       Time        XD  Temp Headng Pitch  Roll #vv DSID$addLayout\n");
-				printf("-------------------------------------------------------------------\n");
-			}
-		}
-
-		$dumpEns = sub ($)
-		{
-			my($e) = @_;
+		$feib = $dta{ENSEMBLE}[$e]->{NUMBER} unless defined($feib); 						# update first and last ensembles in current bin
+		$leib = $dta{ENSEMBLE}[$e]->{NUMBER};
 
-			printf('%5d %s %s %s %5.1f',
-				$dta{ENSEMBLE}[$e]->{NUMBER},
-				$dta{ENSEMBLE}[$e]->{DATE},
-				$dta{ENSEMBLE}[$e]->{TIME},
-				$dta{ENSEMBLE}[$e]->{XDUCER_FACING_UP} ? "UP" : "DN",
-				$dta{ENSEMBLE}[$e]->{TEMPERATURE},
-			);
-			if (defined($dta{ENSEMBLE}[$e]->{HEADING})) { printf(' %6.1f',$dta{ENSEMBLE}[$e]->{HEADING}); }
-			else { printf('    nan'); }
-			if (defined($dta{ENSEMBLE}[$e]->{PITCH})) { printf(' %5.1f',$dta{ENSEMBLE}[$e]->{PITCH}); }
-			else { printf('   nan'); }
-			if (defined($dta{ENSEMBLE}[$e]->{ROLL})) { printf(' %5.1f',$dta{ENSEMBLE}[$e]->{ROLL}); }
-			else { printf('   nan'); }
-			printf(' %3d 0x%02X',
-				$dta{ENSEMBLE}[$e]->{BIN1VELS},
-				$dta{ENSEMBLE}[$e]->{DATA_SOURCE_ID},
-			);
-			printf(' 0x%08X',$dta{ENSEMBLE}[$e]->{ERROR_STATUS_WORD})
-				if ($dta{FIXED_LEADER_BYTES} >= 53);
+		$dn_s += $dn; $dn_n++;																# day number
+		    
+		my($ssCorr) = defined($opt_S) ? ssCorr($dta{ENSEMBLE}[$e],$SS_salin,$SS_temp,$SS_depth) : 1;
+		for (my($b)=0; $b<$dta{N_BINS}; $b++) {
+			my(@v,$this_w12,$this_w34);
+			if ($beamCoords) {																# convert to earth coordinates
+				undef($dta{ENSEMBLE}[$e]->{VELOCITY}[$b][0])
+					if (($dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][0] < $opt_p) || ($opt_d == 1));
+				undef($dta{ENSEMBLE}[$e]->{VELOCITY}[$b][1])
+					if (($dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][1] < $opt_p) || ($opt_d == 2));
+				undef($dta{ENSEMBLE}[$e]->{VELOCITY}[$b][2])
+					if (($dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][2] < $opt_p) || ($opt_d == 3));
+				undef($dta{ENSEMBLE}[$e]->{VELOCITY}[$b][3])
+					if (($dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][3] < $opt_p) || ($opt_d == 4));
+				($dummy,$this_w12,$dummy,$this_w34) =
+					velBeamToBPEarth(\%dta,$e,@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]});
+				@v = $opt_b ? @{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]} :
+							  velBeamToEarth(\%dta,$e,@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]});
+			} else {
+				@v = velApplyHdgBias(\%dta,$e,@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]});
+			}
+			$v[0] *= $ssCorr if defined($v[0]); $v[1] *= $ssCorr if defined($v[1]); 		# apply sound-speed correction
+			$v[2] *= $ssCorr if defined($v[2]); $v[3] *= $ssCorr if defined($v[3]);
+			$w12  *= $ssCorr if defined($w12);	$w34  *= $ssCorr if defined($w34);
+    
+			$v1[$b]+=$v[0],$v1_n[$b]++ if defined($v[0]);									# update sums and nsamps
+			$v2[$b]+=$v[1],$v2_n[$b]++ if defined($v[1]);
+			$v3[$b]+=$v[2],$v3_n[$b]++ if defined($v[2]);
+			$v4[$b]+=$v[3],$v4_n[$b]++ if defined($v[3]);
+			$w12[$b]+=$this_w12,$w12_n[$b]++ if defined($this_w12);
+			$w34[$b]+=$this_w34,$w34_n[$b]++ if defined($this_w34);
+			$c1[$b]+=$dta{ENSEMBLE}[$e]->{CORRELATION}[$b][0],$c1_n[$b]++;
+			$c2[$b]+=$dta{ENSEMBLE}[$e]->{CORRELATION}[$b][1],$c2_n[$b]++;
+			$c3[$b]+=$dta{ENSEMBLE}[$e]->{CORRELATION}[$b][2],$c3_n[$b]++;
+			$c4[$b]+=$dta{ENSEMBLE}[$e]->{CORRELATION}[$b][3],$c4_n[$b]++;
+			$a1[$b]+=$dta{ENSEMBLE}[$e]->{ECHO_AMPLITUDE}[$b][0],$a1_n[$b]++;
+			$a2[$b]+=$dta{ENSEMBLE}[$e]->{ECHO_AMPLITUDE}[$b][1],$a2_n[$b]++;
+			$a3[$b]+=$dta{ENSEMBLE}[$e]->{ECHO_AMPLITUDE}[$b][2],$a3_n[$b]++;
+			$a4[$b]+=$dta{ENSEMBLE}[$e]->{ECHO_AMPLITUDE}[$b][3],$a4_n[$b]++;
+			$p1[$b]+=$dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][0],$p1_n[$b]++;
+			$p2[$b]+=$dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][1],$p2_n[$b]++;
+			$p3[$b]+=$dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][2],$p3_n[$b]++;
+			$p4[$b]+=$dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][3],$p4_n[$b]++;
+
+			my($fi) = 0;
 			foreach my $f (@addFields) {
 				my($fn,$fi) = ($f =~ m{([^[]*)(\[.*)});
 				$fn = $f unless defined($fn);
-				my($v) = eval("\$dta{ENSEMBLE}[$e]->{$fn}$fi");
-				print(defined($v) ? " $v" : " nan");
+				my($val) = eval("\$dta{ENSEMBLE}[$e]->{$fn}$fi");
+				$af[$fi][$b]+=$val,$af_n[$fi][$b]++ if defined($val);
+				$fi++;
 			}
-			print(" BUILT-IN-TEST ERROR")
-				if defined($dta{ENSEMBLE}[$e]->{BUILT_IN_TEST_ERROR});
-			print("\n");
 		}
-		
+	}
+
+} else {											# neither ANTS nor ens files (DEFAULT OUTPUT)
+	if ($dta{FIXED_LEADER_BYTES} >= 53) {
+		printf("Ens # Date		 Time		 XD  Temp Headng Pitch	Roll #vv DSID ESW$addLayout\n");
+		printf("----------------------------------------------------------------------------\n");
+	} else {
+		printf("Ens # Date		 Time		 XD  Temp Headng Pitch	Roll #vv DSID$addLayout\n");
+		printf("------------------------------------------------------------------------\n");
 	}
 
-	for ($e=0; $e<=$#{$dta{ENSEMBLE}}; $e++) {
-		next if (defined($first_ens) &&
-				 $dta{ENSEMBLE}[$e]->{NUMBER} < $first_ens);
-		last if (defined($last_ens) &&
-				 $dta{ENSEMBLE}[$e]->{NUMBER} > $last_ens);
-		 $dta{ENSEMBLE}[$e]->{BIN1VELS} =
-				defined($dta{ENSEMBLE}[$e]->{VELOCITY}[1][0]) +
-				defined($dta{ENSEMBLE}[$e]->{VELOCITY}[1][1]) +
-				defined($dta{ENSEMBLE}[$e]->{VELOCITY}[1][2]) +
-	            defined($dta{ENSEMBLE}[$e]->{VELOCITY}[1][3]);
-		next if ($opt_w && $dta{ENSEMBLE}[$e]->{BIN1VELS}<3);
+	$dumpEns = sub ($)
+	{
+		my($e) = @_;
 
-		die("3-beams used in ensemble #$dta{ENSEMBLE}[$e]->{NUMBER}\n")
-			if ($dta{ENSEMBLE}[$e]->{N_BEAMS_USED} < 4);
-		die("BIT error in ensemble $dta{ENSEMBLE}[$e]->{NUMBER}\n")
-			if ($opt_Q || $opt_A || $opt_E) && defined($dta{ENSEMBLE}[$e]->{BUILT_IN_TEST_ERROR});
-		die("Low gain in ensemble #$dta{ENSEMBLE}[$e]->{NUMBER}\n")
-			if ($dta{ENSEMBLE}[$e]->{LOW_GAIN});
+		printf('%5d %s %s %s %5.1f',
+			$dta{ENSEMBLE}[$e]->{NUMBER},
+			$dta{ENSEMBLE}[$e]->{DATE},
+			$dta{ENSEMBLE}[$e]->{TIME},
+			$dta{ENSEMBLE}[$e]->{XDUCER_FACING_UP} ? "UP" : "DN",
+			$dta{ENSEMBLE}[$e]->{TEMPERATURE},
+		);
+		if (defined($dta{ENSEMBLE}[$e]->{HEADING})) { printf(' %6.1f',$dta{ENSEMBLE}[$e]->{HEADING}); }
+		else { printf('    nan'); }
+		if (defined($dta{ENSEMBLE}[$e]->{PITCH})) { printf(' %5.1f',$dta{ENSEMBLE}[$e]->{PITCH}); }
+		else { printf('   nan'); }
+		if (defined($dta{ENSEMBLE}[$e]->{ROLL})) { printf(' %5.1f',$dta{ENSEMBLE}[$e]->{ROLL}); }
+		else { printf('   nan'); }
+		printf(' %3d 0x%02X',
+			$dta{ENSEMBLE}[$e]->{BIN1VELS},
+			$dta{ENSEMBLE}[$e]->{DATA_SOURCE_ID},
+		);
+		printf(' 0x%08X',$dta{ENSEMBLE}[$e]->{ERROR_STATUS_WORD})
+			if ($dta{FIXED_LEADER_BYTES} >= 53);
+		foreach my $f (@addFields) {
+			my($fn,$fi) = ($f =~ m{([^[]*)(\[.*)});
+			$fn = $f unless defined($fn);
+			my($v) = eval("\$dta{ENSEMBLE}[$e]->{$fn}$fi");
+			print(defined($v) ? " $v" : " nan");
+		}
+		print(" BUILT-IN-TEST ERROR")
+			if defined($dta{ENSEMBLE}[$e]->{BUILT_IN_TEST_ERROR});
+		print("\n");
+	}
+} # define output format
 
-		&$dumpEns($e)
-			unless ($opt_Q);
-	}
+#----------------------------------------------------------------------
+# Loop Over Ensembles
+#----------------------------------------------------------------------
+
+for ($e=0; $e<=$#{$dta{ENSEMBLE}}; $e++) {
+	 $dta{ENSEMBLE}[$e]->{BIN1VELS} =
+			defined($dta{ENSEMBLE}[$e]->{VELOCITY}[1][0]) +
+			defined($dta{ENSEMBLE}[$e]->{VELOCITY}[1][1]) +
+			defined($dta{ENSEMBLE}[$e]->{VELOCITY}[1][2]) +
+			defined($dta{ENSEMBLE}[$e]->{VELOCITY}[1][3]);
+	next if ($opt_w && $dta{ENSEMBLE}[$e]->{BIN1VELS}<3);
+
+	die("BIT error in ensemble $dta{ENSEMBLE}[$e]->{NUMBER}\n")
+		if ($opt_A || $dump_ens) && defined($dta{ENSEMBLE}[$e]->{BUILT_IN_TEST_ERROR});
+	die("Low gain in ensemble #$dta{ENSEMBLE}[$e]->{NUMBER}\n")
+		if ($dta{ENSEMBLE}[$e]->{LOW_GAIN});
+
+	&$dumpEns($e);
 }
 
+&$dumpEns(-1) if defined($opt_T);										# dump final bin
+
 exit(0);	
diff --git a/listHdr b/listHdr
--- a/listHdr
+++ b/listHdr
@@ -2,9 +2,9 @@
 #======================================================================
 #                    L I S T H D R 
 #                    doc: Sat Jan 18 18:41:49 2003
-#                    dlm: Thu Dec  7 10:45:31 2017
+#                    dlm: Tue Mar 20 12:18:03 2018
 #                    (c) 2003 A.M. Thurnherr
-#                    uE-Info: 79 60 NIL 0 0 72 10 2 4 NIL ofnI
+#                    uE-Info: 74 9 NIL 0 0 72 10 2 4 NIL ofnI
 #======================================================================
 
 # Print useful info from the RDI BB header
@@ -68,8 +68,10 @@
 
 	if ($opt_s) {														# summary ANTS output
 		my($id) = $ARGV[0];
-		$id =~ s/00[0-9]\.000//;										# leave just deployment name for std RDI files
-		$id =~ s@^.*/([^/]+)@\1@;
+		if ($id =~ /^\w{5}\d{3}\.\d{3}/) {								# leave just deployment name for std RDI files
+			$id =~ s/00[0-9]\.000//;										
+			$id =~ s@^.*/([^/]+)@\1@;
+		}
 		if ($valid) {
 			printf("%s %d %.1f %d %g %d %.1f\n",
 				$id,$hdr{SERIAL_NUMBER},$hdr{BEAM_FREQUENCY},
diff --git a/listVels b/listVels
--- a/listVels
+++ b/listVels
@@ -2,15 +2,16 @@
 #======================================================================
 #                    L I S T V E L S 
 #                    doc: Mon Apr 25 21:12:54 2016
-#                    dlm: Sat Dec 23 16:10:33 2017
+#                    dlm: Thu Oct  4 16:03:35 2018
 #                    (c) 2016 A.M. Thurnherr
-#                    uE-Info: 21 28 NIL 0 0 72 10 2 4 NIL ofnI
+#                    uE-Info: 14 61 NIL 0 0 72 10 2 4 NIL ofnI
 #======================================================================
 
 # list water-track velocity samples as ANTS records (PD02ANTS)
 
 # HISTORY:
 #	Apr 25, 2016: - created from [listBins]
+#	Oct  4, 2016: - removed pointless transducer config check
 
 # General Notes:
 #	- everything (e.g. beams) is numbered from 1
@@ -95,8 +96,6 @@
 	$P{last_ens} = $dta{ENSEMBLE}[$e]->{NUMBER};
 	$le = $e;
 
-	die("3-beams used in ensemble #$dta{ENSEMBLE}[$e]->{NUMBER}\n")
-		if ($dta{ENSEMBLE}[$e]->{N_BEAMS_USED} < 4);
 	die("BIT error in ensemble $dta{ENSEMBLE}[$e]->{NUMBER}\n")
 		if defined($dta{ENSEMBLE}[$e]->{BUILT_IN_TEST_ERROR});
 	die("Low gain in ensemble #$dta{ENSEMBLE}[$e]->{NUMBER}\n")
diff --git a/listW b/listW
--- a/listW
+++ b/listW
@@ -2,9 +2,9 @@
 #======================================================================
 #                    L I S T W 
 #                    doc: Wed Mar 24 06:45:09 2004
-#                    dlm: Thu Mar 17 07:45:02 2016
+#                    dlm: Mon Apr  2 18:32:03 2018
 #                    (c) 2004 A.M. Thurnherr
-#                    uE-Info: 24 49 NIL 0 0 72 0 2 4 NIL ofnI
+#                    uE-Info: 25 61 NIL 0 0 72 0 2 4 NIL ofnI
 #======================================================================
 
 # dump vertical velocities
@@ -22,6 +22,7 @@
 #   Nov 25, 2013: - checkEnsemble() expunged
 #	Mar 17, 2016: - removed warning
 #				  - updated ancient library names
+#	Apr  2, 2018: - BUG: velBeamToInstrument() used old usage
 
 $0 =~ m{(.*)/[^/]+}; 
 require "$1/RDI_PD0_IO.pl";
@@ -102,7 +103,7 @@
 					 $dta{ENSEMBLE}[$ens]->{PERCENT_GOOD}[$i][2] < 100 ||
 					 $dta{ENSEMBLE}[$ens]->{PERCENT_GOOD}[$i][3] < 100);
 			@v = velInstrumentToEarth(\%dta,$ens,
-					velBeamToInstrument(\%dta,
+					velBeamToInstrument(\%dta,$ens,
 						@{$dta{ENSEMBLE}[$ens]->{VELOCITY}[$i]}));
 		} else {
 			next if ($dta{ENSEMBLE}[$ens]->{PERCENT_GOOD}[$i][0] > 0 ||
diff --git a/loadANTS.m b/loadANTS.m
old mode 100644
new mode 100755
diff --git a/meanProf b/meanProf
--- a/meanProf
+++ b/meanProf
@@ -12,9 +12,15 @@
 # HISTORY:
 #	Feb 22, 2008: - created from [listBins]
 #	Mar 16, 2016: - adapted to new Getopt library
+#	Apr  2, 2018: - BUG: velBeamToInstrument() used old usage
+#	Apr  9, 2018: - adapted to "new" readData() ensemble limits
+#				  - added -l to set final bin
+#				  - BUG: division by zero in empty bins
+#	Apr 10, 2018: - activate output
 #	Aug 24, 2018: - BUG: code bombed when there are no 4-beam solutions
 
 # Soundspeed Correction:
+#	- based on first ensemble only
 #	- applied as described in the RDI coord-trans manual
 #	- sound-speed variation over range is ignored (valid for small gradients)
 #	=> - same simple correction for all velocity components
@@ -22,12 +28,11 @@
 
 use Getopt::Std;
 
-$0 =~ m{(.*/)[^/]+};
-require "$1RDI_BB_Read.pl";
-require "$1RDI_Coords.pl";
-require "$1RDI_Utils.pl";
+$ADCP_tools_minVersion = 2.2;
+($ADCP_TOOLS) = ($0 =~ m{(.*/)[^/]+});
+require "$ADCP_TOOLS/ADCP_tools_lib.pl";
 
-die("Usage: $0 [-r)ange <first_ens,last_ens>] " .
+die("Usage: $0 [-r)ange <first_ens,last_ens>] [-l)ast <bin>] " .
 			  "[-Q)uiet (stats only)] " .
 			  "[-S)oundspeed correction <salin|*,temp|*,depth|*> " .
 		 	  "[require -4)-beam solutions] [-d)iscard <beam#>] " .
@@ -36,7 +41,7 @@
 			  "[-M)agnetic <declination>] " .
 			  "[-D)epth <depth>] " .
 			  "<RDI file>\n")
-	unless (&getopts("4bd:D:M:p:r:QS:") && @ARGV == 1);
+	unless (&getopts("4bd:D:l:M:p:r:QS:") && @ARGV == 1);
 
 die("$0: -4 and -d are mutually exclusive\n")
 	if ($opt_4 && defined($opt_d));
@@ -53,9 +58,6 @@
 
 $ifn = $ARGV[0];
 
-($first_ens,$last_ens) = split(',',$opt_r)
-	if defined($opt_r);
-
 ($SS_salin,$SS_temp,$SS_depth) = split(',',$opt_S)
 	if defined($opt_S);
 die("$0: Cannot do variable soundspeed correction (implementation restriction)\n")
@@ -69,7 +71,12 @@
 $P{mag_decl} = $opt_M if defined($opt_M);
 
 print(STDERR "reading $ifn: ");
-readData($ifn,\%dta);
+if (defined($opt_r)) {								# read selected range
+	my($fe,$le) = split(',',$opt_r);
+	readData($ifn,\%dta,$fe,$le,$opt_l);
+} else {											# read entire file
+	readData($ifn,\%dta,undef,undef,$opt_l);
+}
 printf(STDERR "%d complete ensembles.\n",scalar(@{$dta{ENSEMBLE}}));
 $dta{HEADING_BIAS} = -$opt_M;						# magnetic declination
 
@@ -88,12 +95,8 @@
 
 $lastGoodBin = 0;
 for ($e=0; $e<=$#{$dta{ENSEMBLE}}; $e++) {			# check/transform velocities
-	next if (defined($first_ens) &&
-			 $dta{ENSEMBLE}[$e]->{NUMBER} < $first_ens);
 	$P{first_ens} = $dta{ENSEMBLE}[$e]->{NUMBER},$fe = $e
 		unless defined($P{first_ens});
-	last if (defined($last_ens) &&
-			 $dta{ENSEMBLE}[$e]->{NUMBER} > $last_ens);
 	$P{last_ens} = $dta{ENSEMBLE}[$e]->{NUMBER};
 	$le = $e;
 
@@ -117,7 +120,7 @@
         }
 		@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]} = $beamCoords
 			? velInstrumentToEarth(\%dta,$e,
-				  velBeamToInstrument(\%dta,@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]})
+				  velBeamToInstrument(\%dta,$e,@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]})
 			  )
 			: velApplyHdgBias(\%dta,$e,@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]})
 				unless ($opt_b);
@@ -191,12 +194,18 @@
 	$mean_corr3[$b] = $sum_corr3[$b] / $nEns; $mean_corr4[$b] = $sum_corr4[$b] / $nEns;
 	$mean_amp1[$b] = $sum_amp1[$b] / $nEns; $mean_amp2[$b] = $sum_amp2[$b] / $nEns;
 	$mean_amp3[$b] = $sum_amp3[$b] / $nEns; $mean_amp4[$b] = $sum_amp4[$b] / $nEns;
-	$mean_pcg1[$b] = $sum_pcg1[$b] / $n_pcg1[$b]; $mean_pcg2[$b] = $sum_pcg2[$b] / $n_pcg2[$b];
-	$mean_pcg3[$b] = $sum_pcg3[$b] / $n_pcg3[$b]; $mean_pcg4[$b] = $sum_pcg4[$b] / $n_pcg4[$b];
-	$mean_u[$b] = $sum_u[$b] / $good_vels[$b]; $mean_v[$b] = $sum_v[$b] / $good_vels[$b];
+	
+	$mean_pcg1[$b] = $sum_pcg1[$b] / $n_pcg1[$b] if ($n_pcg1[$b] > 0);
+    $mean_pcg2[$b] = $sum_pcg2[$b] / $n_pcg2[$b] if ($n_pcg2[$b] > 0);
+	$mean_pcg3[$b] = $sum_pcg3[$b] / $n_pcg3[$b] if ($n_pcg3[$b] > 0); 
+	$mean_pcg4[$b] = $sum_pcg4[$b] / $n_pcg4[$b] if ($n_pcg4[$b] > 0);
+
+	next unless ($good_vels[$b] > 0);
+
+	$mean_u[$b] = $sum_u[$b] / $good_vels[$b]; 
+	$mean_v[$b] = $sum_v[$b] / $good_vels[$b];
 	$mean_w[$b] = $sum_w[$b] / $good_vels[$b];
-	$mean_e[$b] = $sum_e[$b] / ($good_vels[$b] - $three_beam[$b])
-		if ($good_vels[$b] - $three_beam[$b] > 0);
+	$mean_e[$b] = ($good_vels[$b] - $three_beam[$b] > 0) ? $sum_e[$b] / ($good_vels[$b] - $three_beam[$b]) : undef;
 }
 
 for ($e=$fe; $e<=$le; $e++) {
@@ -212,13 +221,13 @@
 		$sumsq_amp4[$b] += ($mean_amp4[$b] - $dta{ENSEMBLE}[$e]->{ECHO_AMPLITUDE}[$b][3])**2;
 
 		$sumsq_pcg1[$b] += ($mean_pcg1[$b] - $dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][0])**2
-			if defined($dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][0]);
+			if defined($mean_pcg1[$b]) && defined($dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][0]);
 		$sumsq_pcg2[$b] += ($mean_pcg2[$b] - $dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][1])**2
-			if defined($dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][1]);
+			if defined($mean_pcg2[$b]) && defined($dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][1]);
 		$sumsq_pcg3[$b] += ($mean_pcg3[$b] - $dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][2])**2
-			if defined($dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][2]);
+			if defined($mean_pcg3[$b]) && defined($dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][2]);
 		$sumsq_pcg4[$b] += ($mean_pcg4[$b] - $dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][3])**2
-			if defined($dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][3]);
+			if defined($mean_pcg4[$b]) && defined($dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][3]);
 
 		next unless ($dta{ENSEMBLE}[$e]->{GOOD_VEL}[$b]);
 
@@ -228,7 +237,8 @@
 
 		next if ($dta{ENSEMBLE}[$e]->{THREE_BEAM}[$b]);
 
-		$sumsq_e[$b] += ($mean_e[$b] - $dta{ENSEMBLE}[$e]->{VELOCITY}[$b][3])**2;
+		$sumsq_e[$b] += ($mean_e[$b] - $dta{ENSEMBLE}[$e]->{VELOCITY}[$b][3])**2
+			if defined($mean_e[$b]);
 	}
 }
 
@@ -253,6 +263,8 @@
 # Calculate Beam Statistics
 #----------------------------------------------------------------------
 
+# not implemented yet
+
 #----------------------------------------------------------------------
 # Produce Output
 #----------------------------------------------------------------------
@@ -262,6 +274,8 @@
 				? ssCorr($dta{ENSEMBLE}[$fe],$SS_salin,$SS_temp,$SS_depth)
 				: 1;
 
+	print("#!/usr/bin/perl -S list\n");
+	chmod(0777&~umask,*STDOUT);
 	print("#ANTS#PARAMS# ");
 	foreach my $k (keys(%P)) {
 		print("$k\{$P{$k}\} ");
diff --git a/mkProfile b/mkProfile
--- a/mkProfile
+++ b/mkProfile
@@ -2,9 +2,9 @@
 #======================================================================
 #                    M K P R O F I L E 
 #                    doc: Sun Jan 19 18:55:26 2003
-#                    dlm: Fri Oct 13 13:39:55 2017
+#                    dlm: Tue May  1 10:12:56 2018
 #                    (c) 2003 A.M. Thurnherr
-#                    uE-Info: 230 37 NIL 0 0 72 0 2 4 NIL ofnI
+#                    uE-Info: 395 0 NIL 0 0 72 0 2 4 NIL ofnI
 #======================================================================
 
 # Make an LADCP Profile by Integrating W (similar to Firing's scan*).
@@ -93,6 +93,8 @@
 #				  - removed warning
 #	Sep 12, 2016: - added %PD0_file
 #	Oct 13, 2017: - added instrument orientation
+#	Apr  2, 2018: - BUG: velBeamToInstrument() used old usage
+#	Apr 24, 2018: - BUG: bin1 was used even with zero blanking
 
 # NOTES:
 #	- the battery values are based on transmission voltages (different
@@ -256,6 +258,7 @@
 # Step 1: Integrate w & determine water depth 
 #======================================================================
 
+$minb = 2 if ($dta{BLANKING_DISTANCE} == 0);
 ($firstgood,$lastgood,$atbottom,$w_gap_time,$zErr,$maxz,$rms_heave_accel) =
 	mk_prof(\%dta,0,$opt_F,$minb,$maxb,$opt_c,$opt_e,$opt_g,$opt_p);
 
@@ -324,7 +327,7 @@
 					 $dta{ENSEMBLE}[$ens]->{PERCENT_GOOD}[$i][2] < 100 ||
 					 $dta{ENSEMBLE}[$ens]->{PERCENT_GOOD}[$i][3] < 100);
 			@v = velInstrumentToEarth(\%dta,$ens,
-					velBeamToInstrument(\%dta,
+					velBeamToInstrument(\%dta,$ens,
 						@{$dta{ENSEMBLE}[$ens]->{VELOCITY}[$i]}));
 		} else {
 			next if ($dta{ENSEMBLE}[$ens]->{PERCENT_GOOD}[$i][0] > 0 ||
diff --git a/patchPD0 b/patchPD0
--- a/patchPD0
+++ b/patchPD0
@@ -2,9 +2,9 @@
 #======================================================================
 #                    P A T C H P D 0 
 #                    doc: Tue Aug 23 20:00:15 2016
-#                    dlm: Sat Dec 23 16:12:01 2017
+#                    dlm: Wed Jun 13 20:35:05 2018
 #                    (c) 2010 A.M. Thurnherr
-#                    uE-Info: 51 24 NIL 0 0 72 0 2 4 NIL ofnI
+#                    uE-Info: 184 88 NIL 0 0 72 2 2 4 NIL ofnI
 #======================================================================
 
 $antsSummary = 'patch TRDI PD0 file with external attitude data';
@@ -17,6 +17,8 @@
 #	Dec  9, 2017: - added $antsSuppressCommonOptions = 1;
 #	Dec 23, 2017: - added support for -c
 #				  - BUG: not backward compatible with old IMP files any more
+#	Jun 13, 2017: - added pitch and roll to -o
+#				  - BUG: ??? does -o handle pitch and roll ANOMALIES correctly?
 
 # PATCH-FILE REQUIREMENTS (ANTS format)
 #	- %LADCP_pitch.mu %LADCP_roll.mu		mean LADCP pitch and roll
@@ -60,8 +62,8 @@
 $antsSuppressCommonOptions = 1;
 &antsUsage('cdhko:pr',2,
 	'[patch -p)itch] [-r)oll] [-h)eading] (none patches all)',
-	'[patch -c)lock with pre-Y2K RTC calues]',
-	'[-o) <heading-offset>] [-k)eep velocities of unpatched ensembles]',
+	'[patch -c)lock with pre-Y2K RTC values]',
+	'[-o) <[pitch,roll,]heading-offset>] [-k)eep velocities of unpatched ensembles]',
 	'[keep original -d)ata-source id]',
 	'<original PD0 file> <patched PD0 file> [external attitude file]');
 
@@ -77,7 +79,7 @@
 # Step 1: Read LADCP Data
 #----------------------------------------------------------------------
 
-readData($LADCP_file,\%LADCP);
+readData($LADCP_file,\%LADCP);												# TRDI PD0 file
 
 #----------------------------------------------------------------------
 # Step 2: Process External Attidue Input to Patch PD0 file
@@ -85,7 +87,7 @@
 
 my($pr_missing,$hdg_missing) = (0,0);
 
-&antsIn();
+&antsIn();																	# load first IMP record
 
 my($ensF) 	= &fnr('LADCP_ens');
 my($pitchF) = &fnr('pitch');
@@ -94,13 +96,36 @@
 my($LADCP_pitch_mean)	= &antsRequireParam('LADCP_pitch.mu');
 my($LADCP_roll_mean)	= &antsRequireParam('LADCP_roll.mu');
 
-my($rho,$crho,$srho);
+
+my($pofs,$rofs) = (0,0);													# apply externally supplied offset(s)
+my($rho,$crho,$srho);												
 if (defined($opt_o)) {
+	my($pofs,$rofs,$hofs) = split(/,/,$opt_o);
+
+	if (defined($pofs)) {													# pitch and roll offsets supplied
+		croak("$0: cannot decode -o $opt_o\n")
+			unless numbersp($pofs,$rofs);
+	} else {																# no pitch and roll, only heading
+		$hofs = $pofs; 
+		$pofs = undef;
+    }
+	croak("$0: cannot decode -o $opt_o\n")									
+		unless numberp($hofs);
+																			# set up heading correction
 	&antsAddParams('IMU_hdg_offset',$P{IMP_hdg_offset})						# backward compatibility
 		if defined($P{IMP_hdg_offset});
-	$rho  = $opt_o - &antsRequireParam('IMU_hdg_offset');
+	$rho  = $hofs - &antsRequireParam('IMU_hdg_offset');					# calculate correction relative to already applied one 
 	$crho = cos(rad($rho));
 	$srho = sin(rad($rho));
+
+	if (defined($pofs)) {													# rotate IMP pitch and roll offsets into new LADCP frame
+		my($IMP_pitch_mean) =  &antsRequireParam('IMP_pitch.mu') * $crho
+									+ &antsRequireParam('IMP_roll.mu') * $srho;
+		my($IMP_roll_mean)  = -&antsRequireParam('IMP_pitch.mu') * $srho
+									+ &antsRequireParam('IMP_roll.mu') * $crho;
+		$LADCP_pitch_mean = $IMP_pitch_mean - $pofs;						# apply externally supplied offsets
+		$LADCP_roll_mean  = $IMP_roll_mean  - $rofs;
+	}
 }
 
 do {
@@ -109,8 +134,8 @@
 		unless ($ants_[0][$ensF] == $LADCP{ENSEMBLE}[$ens]->{NUMBER});
 	$LADCP{ENSEMBLE}[$ens]->{DATA_SOURCE_ID} = 0xA0;
 	
-	if (numbersp($ants_[0][$pitchF],$ants_[0][$rollF])) {
-		if (defined($opt_o)) {
+	if (numbersp($ants_[0][$pitchF],$ants_[0][$rollF])) {					# valid IMP data -> patch LADCP ensemble
+		if (defined($opt_o)) {												# -o set: rotate pitch and roll into correct coordinates
 			my($rot_p) = ($ants_[$r][$pitchF]  * $crho +
 						  $ants_[$r][$rollF]   * $srho);
 			my($rot_r) = (-$ants_[$r][$pitchF] * $srho +
@@ -118,16 +143,16 @@
 			$ants_[$r][$pitchF] = $rot_p;
 			$ants_[$r][$rollF]	= $rot_r;
         } 
-		if ($opt_p) {
+		if ($opt_p) {														# patch pitch
 			$LADCP{ENSEMBLE}[$ens]->{DATA_SOURCE_ID} |= ($opt_p<<2);
 			$LADCP{ENSEMBLE}[$ens]->{PITCH} = RDI_pitch($LADCP_pitch_mean + $ants_[0][$pitchF],
 														$LADCP_roll_mean  + $ants_[0][$rollF]);
 		}
-		if ($opt_r) {
+		if ($opt_r) {														# patch roll
 			$LADCP{ENSEMBLE}[$ens]->{DATA_SOURCE_ID} |= ($opt_r<<1);
 			$LADCP{ENSEMBLE}[$ens]->{ROLL} = $LADCP_roll_mean + $ants_[0][$rollF];
 		}
-    } else {
+    } else {																# no valid IMP pitch and roll => invalidate LADCP data
     	$pr_missing++;
 		unless ($opt_k)  {
 	    	clearEns(\%LADCP,$ens);
@@ -135,15 +160,15 @@
 	    }
     }
 	    
-    if (numberp($ants_[0][$hdgF])) {
+    if (numberp($ants_[0][$hdgF])) {										# valid IMP heading
     	$LADCP{ENSEMBLE}[$ens]->{DATA_SOURCE_ID} |= $opt_h;
-   		if (defined($opt_o)) {
+   		if (defined($opt_o)) {												# apply offset on -o; otherwise, data are correctly rotated
 	    	$ants_[0][$hdgF] -= $rho;
 	    	$ants_[0][$hdgF] += 360 if ($ants_[0][$hdgF] < 0);
 	    }
-		$LADCP{ENSEMBLE}[$ens]->{HEADING} = $ants_[0][$hdgF]
+		$LADCP{ENSEMBLE}[$ens]->{HEADING} = $ants_[0][$hdgF]				# patch heading
 			if $opt_h;
-	} else {
+	} else {																# no valid IMP heading => invalidate LADCP data
 		$hdg_missing++;
 		unless ($opt_k)  {
 	    	clearEns(\%LADCP,$ens);
@@ -152,11 +177,11 @@
 	}
 } while (&antsIn());
 
-$LADCP{ENSEMBLE}[0]->{DATA_SOURCE_ID} = 0x7F;				# ensure correct DSID (1st ens: orig; 2nd ens: this prog)
+$LADCP{ENSEMBLE}[0]->{DATA_SOURCE_ID} = 0x7F;								# ensure correct DSID (1st ens: orig; 2nd ens: this prog)
 $LADCP{ENSEMBLE}[1]->{DATA_SOURCE_ID} = 0xA0
 	unless ($LADCP{ENSEMBLE}[1]->{DATA_SOURCE_ID}&0xF0 == 0xA0);
 
-writeData($outPD0,\%LADCP);
+writeData($outPD0,\%LADCP);													# write new PD0
 
 my($verb) = $opt_k ? 'retained' : 'cleared';
 printf(STDERR "$outPD0: %d pitch/roll & %d heading values $verb\n",$pr_missing,$hdg_missing)
diff --git a/splitPD0 b/splitPD0
--- a/splitPD0
+++ b/splitPD0
@@ -2,9 +2,9 @@
 #======================================================================
 #                    S P L I T P D 0 
 #                    doc: Sat Aug 21 22:20:27 2010
-#                    dlm: Sat Jul 30 18:50:43 2016
+#                    dlm: Mon Apr  2 15:49:09 2018
 #                    (c) 2010 A.M. Thurnherr
-#                    uE-Info: 69 60 NIL 0 0 72 2 2 4 NIL ofnI
+#                    uE-Info: 28 36 NIL 0 0 72 2 2 4 NIL ofnI
 #======================================================================
 
 # split RDI files based on list of ensemble numbers (e.g. from yoyo -t)
@@ -19,11 +19,13 @@
 #	Jul 30, 2016: - modified -o default
 #				  - added code to set DSID of first ensemble of each
 #					output file to 0x7f7f
+#	Apr  2, 1018: - BUG in error messages
+#				  - added header id check
 
 # NOTES:
 #   - it is assumed that the input file begins with ensemble #1
-#   - turning-point ensembles are written to preceding profile,
-#     for compatibility with [yoyo]
+#   - turning-point ensembles are written to next profile,
+#     for compatibility with [yoyo]?
 
 # FILE NAME CONVENTION:
 #   - in order to assign individual yoyo casts numerical station numbers,
@@ -38,7 +40,7 @@
 
 die("Usage: $0 " .
 	"[-o)ut-file <fmt[e.g. 017DL_%02d.000]>] " .
-	"[-f)irst <cast #>] " .
+	"[-f)irst output <cast #>] " .
 	"[require -m)in <ens> to produce output] " .
 	"<RDI file> <ens> <ens[...]>\n")
 		unless (&getopts('f:o:m:') && @ARGV>=3);
@@ -50,7 +52,8 @@
 
 $opt_m = 0 unless defined($opt_m);								# default: produce tiny files as well
 	
-readHeader($ARGV[0],\%hdr); shift;								# get length of ensembles
+$fn = $ARGV[0]; shift;
+readHeader($fn,\%hdr); 											# get length of ensembles
 $ens_len = $hdr{ENSEMBLE_BYTES} + 2;
 
 $first_ens = $ARGV[0]+1; shift;									# initialize loop
@@ -59,19 +62,21 @@
 
 do {															# split data
 	sysseek(WBRF,($first_ens-2)*$ens_len,0) ||					# begin next block of ensembles
-		die("$WBRcfn: $!");
-	$last_ens++ unless defined($ARGV[0]);
+		die("$fn: $!");
+	$last_ens++ unless defined($fn);
 
-	sysread(WBRF,$ids,2) || die("$WBRcfn: file truncated");		# read 1st ensemble & ensure DSID is 0x7f
+	sysread(WBRF,$ids,2) || die("$fn: file truncated");			# read 1st ensemble & ensure DSID is 0x7f
+	die("$fn: illegal header id [0x" . unpack('H4',$ids) . "]")	# require 1st byte to be 0x7f
+		unless (substr(unpack('H4',$ids),0,2) eq '7f');
 	$ids = pack('H4','7f7f');									# reset DSID
 	sysread(WBRF,$febuf,$ens_len-4) == $ens_len-4 ||
-		die("$WBRcfn: file truncated");
-	sysread(WBRF,$csum,2) || die("$WBRcfn: file truncated");
+		die("$fn: file truncated");
+	sysread(WBRF,$csum,2) || die("$fn: file truncated");
 	$csum = pack('v',unpack('%16C*',$ids.$febuf));				# re-calculate checksum
 	
 	$nBytes = ($last_ens-$first_ens) * $ens_len;				# read remaining ensembles in block
 	sysread(WBRF,$buf,$nBytes) == $nBytes ||
-		die("$WBRcfn: file truncated");
+		die("$fn: file truncated (ends before ens#$last_ens)");
 
 	if ($last_ens-$first_ens+1 >= $opt_m) {						# write output only if sufficient # of ensembles
 		$fn = sprintf($opt_o,$opt_f+$cnr++);