# HG changeset patch
# User A.M. Thurnherr <athurnherr@yahoo.com>
# Date 1457451784 0
# Node ID 7155adf61d7788628cf8e29b39a2a06de443098c
# Parent  b6ca27a1d19c0fa261e7e29be932df374d04b1cf
V1.4

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: Tue Jan  5 13:22:10 2016
+                    dlm: Tue Mar  8 15:41:38 2016
                     (c) 2012 A.M. Thurnherr
-                    uE-Info: 35 46 NIL 0 0 72 3 2 8 NIL ofnI
+                    uE-Info: 55 28 NIL 0 0 72 3 2 8 NIL ofnI
 ======================================================================
 
 May 15, 2012:
@@ -28,8 +28,30 @@
   - merged with Oct 2 version on Studio desktop, which ignores
     initial garbage in PD0 files
 
-Jan  5, 2015: V1.4
-	- added [ADCP_tools_lib.pl] with compile-time version control
-	- [RDI_Coords.pl] added &velEarthToBeam()
-	- updated [listBins] to use versioned libs and calc w12 & w34
-	  from earth-coordinate data correctly
+Jan  5, 2016: V1.4
+  - added [ADCP_tools_lib.pl] with compile-time version control
+  - [RDI_Coords.pl] added &velEarthToBeam()
+  - updated [listBins] to use versioned libs and calc w12 & w34
+    from earth-coordinate data correctly
+
+Jan 6, 2016:
+  - minor change to [listBins]
+
+Jan 9, 2016:
+  - added &velEarthToBeam(), &velBeamToEarth() to [RDI_Coords.pl]
+  - improvements to [RDI_PD0_IO.pl]
+  - adapted [listHdr] to producer-id in PD0 files & other minor changes
+  - renamed function to solve name conflic in [RDI_Utils.pl]
+
+Jan 9 - Feb 26:
+  - added swap_beams() to [editPD0]
+
+Feb 29, 2016:
+  - improvements to [RDI_PD0_IO.pl]
+  - finished debugging [RDI_Coords.pl]
+
+Mar 8, 2016: V1.4
+  - verified/updated version in [ADCP_tools_lib.pl] [.hg/hgrc]
+  - published V1.4 on server
+
+	  
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: Sun Jan 31 12:42:43 2016
+#                    dlm: Mon Feb 29 18:03:31 2016
 #                    (c) 2003 A.M. Thurnherr
-#                    uE-Info: 185 0 NIL 0 0 72 0 2 4 NIL ofnI
+#                    uE-Info: 44 85 NIL 0 0 72 0 2 4 NIL ofnI
 #======================================================================
 
 # RDI Workhorse Coordinate Transformations
@@ -41,6 +41,8 @@
 #				  - removed unused code from &velBeamToBPInstrument
 #	Jan  5, 2016: - added &velEarthToInstrument(@), &velInstrumentToBeam(@)
 #	Jan  9, 2016: - added &velEarthToBeam(), &velBeamToEarth()
+#	Feb 29, 2016: - debugged & verified velEarthToInstrument(), velInstrumentToBeam()
+#				  - added velBeamToEarth()
 
 use strict;
 use POSIX;
@@ -157,9 +159,17 @@
 	}
 } # STATIC SCOPE
 
+
+sub velBeamToEarth(@)
+{
+	my($dtaR,$e,@v) = @_;
+	return velInstrumentToEarth($dtaR,$e,velBeamToInstrument($dtaR,@v));
+}
+
 #----------------------------------------------------------------------
 # velEarthToInstrument() transforms earth to instrument coordinates
 #	- based on manually inverted rotation matrix M (Sec 5.6 in coord-trans manual)
+#	- code was verified for both down- and uplookers
 #	- missing heading data (IMP) causes undef beam velocities
 #----------------------------------------------------------------------
 
@@ -170,26 +180,22 @@
 	{
 		my($dta,$ens,$u,$v,$w,$ev) = @_;
 
-		unless (@E2I) {
+		unless (@E2I &&
+				$pitch == $dta->{ENSEMBLE}[$ens]->{PITCH} &&
+				$roll  == $dta->{ENSEMBLE}[$ens]->{ROLL}) {
 			$hdg = $dta->{ENSEMBLE}[$ens]->{HEADING} - $dta->{HEADING_BIAS} 
 				if defined($dta->{ENSEMBLE}[$ens]->{HEADING});
 			$pitch = $dta->{ENSEMBLE}[$ens]->{PITCH};
 			$roll  = $dta->{ENSEMBLE}[$ens]->{ROLL};
 			my($rad_gimbal_pitch) = atan(tan(rad($pitch)) * cos(rad($roll)));
+			my($useRoll) = ($dta->{ENSEMBLE}[$ens]->{XDUCER_FACING_UP}) ? $roll+180 : $roll;
 			my($sh,$ch) = (sin(rad($hdg)),cos(rad($hdg)))
 				if defined($hdg);				
 			my($sp,$cp) = (sin($rad_gimbal_pitch),cos($rad_gimbal_pitch));
-			my($sr,$cr) = (sin(rad($roll)),	cos(rad($roll)));
-			@E2I = $dta->{ENSEMBLE}[$ens]->{XDUCER_FACING_UP}
-				 ? (
-					[$ch*-$cr+$sh*$sp*-$sr,	 $ch*$sp*-$sr-$sh*-$cr,	$cp*-$sr],
-				    [$sh*$cp,				 $ch*$cp,				$sp		],
-				    [$ch*-$sr-$sh*$sp*-$cr,	-$sh*-$sr-$ch*$sp*-$cr,	$cp*-$cr]
-				 ) : (
-					[$ch*$cr+$sh*$sp*$sr,	 $ch*$sp*$sr-$sh*$cr,	$cp*$sr	],
-				    [$sh*$cp,				 $ch*$cp,				$sp		],
-				    [$ch*$sr-$sh*$sp*$cr,	-$sh*$sr-$ch*$sp*$cr,	$cp*$cr	]
-				 );
+			my($sr,$cr) = (sin(rad($useRoll)),	  cos(rad($useRoll)));
+			@E2I = ([$ch*$cr+$sh*$sp*$sr,	 $ch*$sp*$sr-$sh*$cr,	-$cp*$sr],		# M^-1 = R^-1 * P^-1 * R^-1
+				    [$sh*$cp,				 $ch*$cp,				$sp	],
+				    [$ch*$sr-$sh*$sp*$cr,	-$sh*$sr-$ch*$sp*$cr,	$cp*$cr]);
 		}
 
 		return defined($dta->{ENSEMBLE}[$ens]->{HEADING})
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: Sat Jan  9 17:57:01 2016
+#                    dlm: Mon Feb 29 12:30:04 2016
 #                    (c) 2003 A.M. Thurnherr
-#                    uE-Info: 73 64 NIL 0 0 72 10 2 4 NIL ofnI
+#                    uE-Info: 1232 63 NIL 0 0 72 10 2 4 NIL ofnI
 #======================================================================
 
 # Read RDI BroadBand Binary Data Files (*.[0-9][0-9][0-9])
@@ -71,6 +71,11 @@
 #				  - added PRODUCER
 #				  - BUG: writeData() did not work correctly for ECOGIG OC26 moored data (spaces in filename?)
 #				  - added support for patching coordinate system
+#	Feb 16, 2016: - added transducer orientation to WBPens()
+#				  - BUG: most WBPens() error messages used wrong file name
+#	Feb 23, 2016: - changed WBRhdr() to use 2nd ensemble (with correct data-source id)
+#	Feb 26, 2016: - added basic BT data to WBPens(); not BT_RL_* and BT_SIGNAL_STRENGTH
+#	Feb 29, 2016: - LEAP DAY: actually got BT data to work
 
 # FIRMWARE VERSIONS:
 #	It appears that different firmware versions generate different file
@@ -86,6 +91,11 @@
 
 # PD0 FILE FORMAT EXTENSIONS:
 #
+#	- file creator encoded in DATA_SOURCE_ID
+#
+#	- first ensemble uses default RDI DATA_SOURCE_ID because the LDEO_IX
+#	  software assumes this
+#
 #	- DATA_SOURCE_ID = 0x7F						original TRDI PD0 file
 #
 #	- DATA_SOURCE_ID = 0xA0 | PATCHED_MASK		produced by IMP+LADP 
@@ -376,7 +386,14 @@
 	# HEADER
 	#--------------------
 
-	$start_ens = skip_initial_trash();
+	skip_initial_trash();
+	sysread(WBRF,$buf,6) == 6 || die("$WBRcfn: $!");
+	($hid,$did,$dta->{ENSEMBLE_BYTES},$dummy,$dta->{NUMBER_OF_DATA_TYPES})
+		= unpack('CCvCC',$buf);
+	$hid == 0x7f || die(sprintf($FmtErr,$WBRcfn,"Header",$hid,0));
+	$did == 0x7f || die(sprintf($FmtErr,$WBRcfn,"Header",$did,0));
+
+	$start_ens = sysseek(WBRF,$dta->{ENSEMBLE_BYTES}-6+2,1) || die("$WBRcfn: $!");
 	sysread(WBRF,$buf,6) == 6 || die("$WBRcfn: $!");
 	($hid,$did,$dta->{ENSEMBLE_BYTES},$dummy,$dta->{NUMBER_OF_DATA_TYPES})
 		= unpack('CCvCC',$buf);
@@ -910,14 +927,12 @@
 
 		next if ($nxt == $ndt);													# no BT found => next ens
 
-		sysseek(WBRF,14,1) || die("$WBRcfn: $!");								# BT config
-
+		sysseek(WBRF,14,1) || die("$WBRcfn: $!");								# BT range, velocity, corr, %-good, ...
 		sysread(WBRF,$buf,28) == 28 || die("$WBRcfn: $!");
 		@dta = unpack('v4v4C4C4C4',$buf);
-	    
 		for ($beam=0; $beam<4; $beam++) {
-			${$E}[$ens]->{BT_RANGE}[$beam] = $dta[$beam] / 100
-					if ($dta[$beam]);
+			${$E}[$ens]->{BT_RANGE}[$beam] = $dta[$beam] / 100					# lower 2 bytes only! 
+					if ($dta[$beam]);											# (see below for high bytes)
 		}
 		for ($beam=0; $beam<4; $beam++) {
 			${$E}[$ens]->{BT_VELOCITY}[$beam] =
@@ -928,18 +943,16 @@
 			${$E}[$ens]->{BT_CORRELATION}[$beam] = $dta[8+$beam]
 				if ($dta[8+$beam]);
 		}
-		for ($beam=0; $beam<4; $beam++) {
+		for ($beam=0; $beam<4; $beam++) {										# BT filter parameter
 			${$E}[$ens]->{BT_EVAL_AMPLITUDE}[$beam] = $dta[12+$beam];
 		}
 		for ($beam=0; $beam<4; $beam++) {
 			${$E}[$ens]->{BT_PERCENT_GOOD}[$beam] = $dta[16+$beam];
 		}
 
-		sysseek(WBRF,6,1) || die("$WBRcfn: $!");		# BT config
-
+		sysseek(WBRF,6,1) || die("$WBRcfn: $!");								# BT ref level stuff
 		sysread(WBRF,$buf,20) == 20 || die("$WBRcfn: $!");
 		@dta = unpack('v4C4C4C4',$buf);
-
 		for ($beam=0; $beam<4; $beam++) {
 			${$E}[$ens]->{BT_RL_VELOCITY}[$beam] =
 				unpack('s',pack('S',$dta[$beam])) / 1000
@@ -956,17 +969,15 @@
 			${$E}[$ens]->{BT_RL_PERCENT_GOOD}[$beam] = $dta[12+$beam];
 		}
 
-		sysseek(WBRF,2,1) || die("$WBRcfn: $!");		# BT config
-
+		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++) {
+		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]);
 		}
@@ -1051,7 +1062,18 @@
 		# Fixed Leader
 		#--------------------
 	
-		sysseek(WBPF,$start_ens+$WBPofs[0]+25,0) || die("$WBPcfn: $!");		# jump to EX/Coord-Transform
+		sysseek(WBPF,$start_ens+$WBPofs[0]+4,0)								# system config (transducer orientation)
+			|| die("$WBPcfn: $!");
+		sysread(WBPF,$buf,1) == 1 || die("$WBPcfn: $!");
+		$B1 = unpack('C',$buf);
+		$B1 |= 0x80 if ($dta->{ENSEMBLE}[$ens]->{XDUCER_FACING_UP});
+		$B1 &= 0x7F if ($dta->{ENSEMBLE}[$ens]->{XDUCER_FACING_DOWN});
+		$buf = pack('C',$B1);
+		sysseek(WBPF,$start_ens+$WBPofs[0]+4,0)
+			|| die("$WBPcfn: $!");
+		syswrite(WBPF,$buf,1) == 1 || die("$WBPcfn: $!");
+
+		sysseek(WBPF,$start_ens+$WBPofs[0]+25,0) || die("$WBPcfn: $!");		# EX / coord-transformation
 		sysread(WBPF,$buf,1) == 1 || die("$WBPcfn: $!");
 		my($EX) = unpack('C',$buf);
 		if ($dta->{BEAM_COORDINATES}) {
@@ -1073,9 +1095,9 @@
 		
 		$dta->{ENSEMBLE}[$ens]->{XDUCER_DEPTH} = round($dta->{ENSEMBLE}[$ens]->{XDUCER_DEPTH}*10);
 
-		#-----------------------------
-		# IMP allows for missing value
-		#-----------------------------
+		#---------------------------------
+		# NB: IMP allows for missing value
+		#---------------------------------
 
 		$dta->{ENSEMBLE}[$ens]->{HEADING} = defined($dta->{ENSEMBLE}[$ens]->{HEADING})
 							   ? round($dta->{ENSEMBLE}[$ens]->{HEADING}*100)
@@ -1104,7 +1126,7 @@
 		# Velocity Data
 		#--------------------
 
-		sysseek(WBPF,$start_ens+$WBPofs[2]+2,0) || die("$WBRcfn: $!");	# skip velocity data id (assume it is correct)
+		sysseek(WBPF,$start_ens+$WBPofs[2]+2,0) || die("$WBPcfn: $!");	# skip velocity data id (assume it is correct)
 		for ($bin=0; $bin<$nbins; $bin++) {
 			for ($beam=0; $beam<4; $beam++) {
 				$dta->{ENSEMBLE}[$ens]->{VELOCITY}[$bin][$beam] = defined($dta->{ENSEMBLE}[$ens]->{VELOCITY}[$bin][$beam])
@@ -1120,7 +1142,7 @@
 		# Correlation Data
 		#--------------------
 
-		sysseek(WBPF,$start_ens+$WBPofs[3]+2,0) || die("$WBRcfn: $!");
+		sysseek(WBPF,$start_ens+$WBPofs[3]+2,0) || die("$WBPcfn: $!");
 		for ($bin=0; $bin<$nbins; $bin++) {
 			for ($beam=0; $beam<4; $beam++) {
 				$buf = pack('C',$dta->{ENSEMBLE}[$ens]->{CORRELATION}[$bin][$beam]);
@@ -1133,7 +1155,7 @@
 		# Echo Intensity Data
 		#--------------------
 
-		sysseek(WBPF,$start_ens+$WBPofs[4]+2,0) || die("$WBRcfn: $!");
+		sysseek(WBPF,$start_ens+$WBPofs[4]+2,0) || die("$WBPcfn: $!");
 
 		for ($bin=0; $bin<$nbins; $bin++) {
 			for ($beam=0; $beam<4; $beam++) {
@@ -1147,7 +1169,7 @@
 		# Percent Good Data
 		#--------------------
 
-		sysseek(WBPF,$start_ens+$WBPofs[5]+2,0) || die("$WBRcfn: $!");
+		sysseek(WBPF,$start_ens+$WBPofs[5]+2,0) || die("$WBPcfn: $!");
 
 		for ($i=0,$bin=0; $bin<$nbins; $bin++) {
 			for ($beam=0; $beam<4; $beam++,$i++) {
@@ -1164,17 +1186,55 @@
 
 		my($nxt);
 		for ($nxt=6; $nxt<$ndt; $nxt++) {										# scan until BT found
-			sysseek(WBPF,$start_ens+$WBPofs[$nxt],0) || die("$WBRcfn: $!");
-			sysread(WBPF,$buf,2) == 2 || die("$WBRcfn: $!");
+			sysseek(WBPF,$start_ens+$WBPofs[$nxt],0) || die("$WBPcfn: $!");
+			sysread(WBPF,$buf,2) == 2 || die("$WBPcfn: $!");
 			$id = unpack('v',$buf);
 			last if ($id == 0x0600);
 		}
 
 		unless ($nxt == $ndt) {													# BT found
-			sysseek(WBPF,14,1) || die("$WBRcfn: $!");							# skip BT config
-			# NOT YET IMPLEMENTED
-		}
+			sysseek(WBPF,14,1) || die("$WBPcfn: $!");							# skip BT config
+			for ($beam=0; $beam<4; $beam++) {									# BT range low bytes (2 per beam)
+				$buf = pack('v',round($dta->{ENSEMBLE}[$ens]->{BT_RANGE}[$beam] * 100) & 0xFFFF);
+				my($nw) = syswrite(WBPF,$buf,2);
+				$nw == 2 || die("$WBPcfn: $nw bytes written ($!)");
+			}
+			
+			for ($beam=0; $beam<4; $beam++) {									# BT velocities
+				$buf = pack('v',unpack('S',pack('s',
+						defined($dta->{ENSEMBLE}[$ens]->{BT_VELOCITY}[$beam])
+							? round($dta->{ENSEMBLE}[$ens]->{BT_VELOCITY}[$beam]*1000)
+							: 0x8000)));
+				my($nw) = syswrite(WBPF,$buf,2);
+				$nw == 2 || die("$WBPcfn: $nw bytes written ($!)");
+			}
+			
+			for ($beam=0; $beam<4; $beam++) {									# BT correlation
+				$buf = pack('C',$dta->{ENSEMBLE}[$ens]->{BT_CORRELATION}[$beam]);
+				my($nw) = syswrite(WBPF,$buf,1);
+				$nw == 1 || die("$WBPcfn: $nw bytes written ($!)");
+			}
+			
+			for ($beam=0; $beam<4; $beam++) {									# BT evaluation amp of matching filter
+                $buf = pack('C',$dta->{ENSEMBLE}[$ens]->{BT_EVAL_AMPLITUDE}[$beam]);
+                my($nw) = syswrite(WBPF,$buf,1);
+                $nw == 1 || die("$WBPcfn: $nw bytes written ($!)");
+            }
+            
+            for ($beam=0; $beam<4; $beam++) {									# BT percent good
+                $buf = pack('C',$dta->{ENSEMBLE}[$ens]->{BT_PERCENT_GOOD}[$beam]);
+                my($nw) = syswrite(WBPF,$buf,1);
+                $nw == 1 || die("$WBPcfn: $nw bytes written ($!)");
+            }
 
+			sysseek(WBPF,33,1) || die("$WBPcfn: $!");							# BT range high bytes (1 per beam)
+			for ($beam=0; $beam<4; $beam++) {
+				$buf = pack('C',(round($dta->{ENSEMBLE}[$ens]->{BT_RANGE}[$beam]*100) & 0xFF0000) >> 16);
+				my($nw) = syswrite(WBPF,$buf,1);
+				$nw == 1 || die("$WBPcfn: $nw bytes written ($!)");
+			}
+
+        } # if BT found
 
 		#----------------
 		# Update Checksum
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 Jan  9 20:07:59 2016
+#                    dlm: Fri Feb 26 17:24:55 2016
 #                    (c) 2013 A.M. Thurnherr
-#                    uE-Info: 128 0 NIL 0 0 72 2 2 4 NIL ofnI
+#                    uE-Info: 26 98 NIL 0 0 72 2 2 4 NIL ofnI
 #======================================================================
 
 # edit RDI PD0 file, e.g. to replace pitch/roll/heading with external values
@@ -23,11 +23,15 @@
 #		swap_beams(<b1>,<b2>)	swap data from beams b1 and b2
 #								input in beam coords required
 #								beam rotation is equivalent to 3 consecutive beam swaps
+#								basic BT data are swapped as well (not RL and not SIGNAL_STRENGTH)
 #
 #		earth2beam()			transform beam to earth coordinates
 #								does not handle bin-remapping
 #								input in beam coords required
 #
+#		ensure_UL()				overwrite transducer-orientation flag in data file
+#		ensure_DL()
+#
 #	- -x notes:
 #		- multiple perl expressions can be combined with ,
 #
@@ -47,6 +51,9 @@
 #				  - change output data-source ID from 0x7F to 0xE0
 #				  - updated getopts to current perl version
 #				  - adapted to [ADCP_tools_lib.pl]
+#	Feb 15, 2016: - added ensure_UL() ensure_DL()
+#	Feb 23, 2016: - added -k
+#	Feb 26, 2016: - added basic BT data to swap_beams()
 
 use Getopt::Std;
 
@@ -57,8 +64,9 @@
 $USAGE = "$0 @ARGV";
 die("Usage: $0 " .
     '-e) <edit-file> | -x) <expr> ' .
+    '-k)eep original data-source id' .
     "<input file> <output file>\n")
-        unless (&getopts('e:x:') && @ARGV == 2);
+        unless (&getopts('ke:x:') && @ARGV == 2);
 
 die("$0: -e <edit-file> or -x <expr> required\n")
     unless (defined($opt_x) || -r $opt_e);
@@ -75,21 +83,66 @@
 # Data Editing Library
 #--------------------------------------------------
 
+#
+# override pitch/roll/heading
+#
 sub p($) { $dta{ENSEMBLE}[$e]->{PITCH} = $_[0]; return 1; }
 sub r($) { $dta{ENSEMBLE}[$e]->{ROLL} = $_[0]; return 1; }
 sub h($) { $dta{ENSEMBLE}[$e]->{HEADING} = $_[0]; return 1;}
 
+#
+# override transducer orientation
+#
+sub ensure_DL()
+{
+	$dta{ENSEMBLE}[$e]->{XDUCER_FACING_DOWN} = 1;
+	$dta{ENSEMBLE}[$e]->{XDUCER_FACING_UP} = undef;
+	return 1;
+}
+sub ensure_UL()
+{
+	$dta{ENSEMBLE}[$e]->{XDUCER_FACING_DOWN} = undef;
+	$dta{ENSEMBLE}[$e]->{XDUCER_FACING_UP} = 1;
+	return 1;
+}
+
+#
+# swap data from two mis-connected beams
+#
 sub swap_beams($$)
 {
 	my($b1,$b2) = @_;
+	my($tmp);
 
 #	print(STDERR "\n entering swap_beams($b1,$b2) for ens = $e...");
 
 	die("$ARGV[0]: beam-coordinate data required\n")
 		unless ($dta{BEAM_COORDINATES});
 
+	if ($dta{BT_PRESENT}) {
+		$tmp = $dta{ENSEMBLE}[$e]->{BT_RANGE}[$b1-1];
+		$dta{ENSEMBLE}[$e]->{BT_RANGE}[$b1-1] = $dta{ENSEMBLE}[$e]->{BT_RANGE}[$b2-1];
+		$dta{ENSEMBLE}[$e]->{BT_RANGE}[$b2-1] = $tmp;
+		                
+		$tmp = $dta{ENSEMBLE}[$e]->{BT_VELOCITY}[$b1-1];
+		$dta{ENSEMBLE}[$e]->{BT_VELOCITY}[$b1-1] = $dta{ENSEMBLE}[$e]->{BT_VELOCITY}[$b2-1];
+		$dta{ENSEMBLE}[$e]->{BT_VELOCITY}[$b2-1] = $tmp;
+		                
+		$tmp = $dta{ENSEMBLE}[$e]->{BT_CORRELATION}[$b1-1];
+		$dta{ENSEMBLE}[$e]->{BT_CORRELATION}[$b1-1] = $dta{ENSEMBLE}[$e]->{BT_CORRELATION}[$b2-1];
+		$dta{ENSEMBLE}[$e]->{BT_CORRELATION}[$b2-1] = $tmp;
+		                
+		$tmp = $dta{ENSEMBLE}[$e]->{BT_EVAL_AMPLITUDE}[$b1-1];
+		$dta{ENSEMBLE}[$e]->{BT_EVAL_AMPLITUDE}[$b1-1] = $dta{ENSEMBLE}[$e]->{BT_EVAL_AMPLITUDE}[$b2-1];
+		$dta{ENSEMBLE}[$e]->{BT_EVAL_AMPLITUDE}[$b2-1] = $tmp;
+		                
+		$tmp = $dta{ENSEMBLE}[$e]->{BT_PERCENT_GOOD}[$b1-1];
+		$dta{ENSEMBLE}[$e]->{BT_PERCENT_GOOD}[$b1-1] = $dta{ENSEMBLE}[$e]->{BT_PERCENT_GOOD}[$b2-1];
+		$dta{ENSEMBLE}[$e]->{BT_PERCENT_GOOD}[$b2-1] = $tmp;
+    }		                
+
 	for (my($bin)=0; $bin<$dta{N_BINS}; $bin++) {
-		my($tmp) = $dta{ENSEMBLE}[$e]->{VELOCITY}[$bin][$b1-1];
+		$tmp = $dta{ENSEMBLE}[$e]->{VELOCITY}[$bin][$b1-1];
 		$dta{ENSEMBLE}[$e]->{VELOCITY}[$bin][$b1-1] = $dta{ENSEMBLE}[$e]->{VELOCITY}[$bin][$b2-1];
 		$dta{ENSEMBLE}[$e]->{VELOCITY}[$bin][$b2-1] = $tmp;
 
@@ -109,6 +162,9 @@
 }
 
 
+#
+# transform earth to beam coordinates
+#
 { my($checked);
 
 	sub earth2beam()
@@ -160,7 +216,8 @@
 }
 
 for (local($e)=my($eei)=0; $e<@{$dta{ENSEMBLE}}; $e++) {						# local() needed for p(), r(), h()
-	$dta{ENSEMBLE}[$e]->{DATA_SOURCE_ID} = 0xE0;								# mark all ensembles
+	$dta{ENSEMBLE}[$e]->{DATA_SOURCE_ID} = 0xE0									# mark all ensembles except first
+		unless ($opt_k || $e==0);
 	if ($EE[$eei] eq '*' || $EE[$eei] == $dta{ENSEMBLE}[$e]->{NUMBER}) {		# match => edit
 		eval($EX[$eei]) || die("$@ while executing <$EX[$eei]>\n");
 	} elsif ($EE[$eei] > $dta{ENSEMBLE}[$e]->{NUMBER}) {						# next edit later in file => skip
diff --git a/listBins b/listBins
--- 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: Sun Jan 31 12:43:28 2016
+#                    dlm: Mon Feb 29 18:10:25 2016
 #                    (c) 2006 A.M. Thurnherr
-#                    uE-Info: 340 0 NIL 0 0 72 10 2 4 NIL ofnI
+#                    uE-Info: 59 63 NIL 0 0 72 10 2 4 NIL ofnI
 #======================================================================
 
 # Split data file into per-bin time series.
@@ -54,7 +54,9 @@
 #	Jan  5, 2016: - adapted to [ANTS_tools_lib.pl]
 #				  - adapted to calculation of w12, w34 from earth-coordinate data
 #				  - several other changes to the code that should not affect the results
-#	Jan  6, 2015: - -b removed (always output beamvels)
+#	Jan  6, 2016: - -b removed (always output beamvels)
+#	Jan 31, 2016: - started debugging the obviously wrong Earth2Beam() transformation
+#	Feb 29, 2016: - continued debugging; removed debugging code
 
 # General Notes:
 #	- everything (e.g. beams) is numbered from 1
@@ -318,9 +320,7 @@
 				velBeamToBPEarth(\%dta,$e,@{$dta{ENSEMBLE}[$e]->{BEAM_VELOCITY}[$b]});
 
 			@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]} = 						# calculate earth velocities
-				velInstrumentToEarth(\%dta,$e,
-					velBeamToInstrument(\%dta,@{$dta{ENSEMBLE}[$e]->{BEAM_VELOCITY}[$b]})
-			  	);
+				velBeamToEarth(\%dta,@{$dta{ENSEMBLE}[$e]->{BEAM_VELOCITY}[$b]});
 			$dta{ENSEMBLE}[$e]->{THREE_BEAM}[$b] = $RDI_Coords::threeBeamFlag;
 			$three_beam[$b] += $RDI_Coords::threeBeamFlag;
 
@@ -329,18 +329,8 @@
 				next;
 			}
 		} else { 															# Earth coordinates
-
-			my(@iVels) = &velEarthToInstrument(\%dta,$e,@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]});
-			my(@eVels) = &velInstrumentToEarth(\%dta,$e,@iVels);
-#			printf(STDERR "err($e,$b) = %g,%g/%g,%g/%g,%g/%g,%g\n",
-#												$dta{ENSEMBLE}[$e]->{VELOCITY}[$b][0],$iVels[0],
-#												$dta{ENSEMBLE}[$e]->{VELOCITY}[$b][1],$iVels[1],
-#												$dta{ENSEMBLE}[$e]->{VELOCITY}[$b][2],$iVels[2],
-#												$dta{ENSEMBLE}[$e]->{VELOCITY}[$b][3],$iVels[3]);
-
 			@{$dta{ENSEMBLE}[$e]->{BEAM_VELOCITY}[$b]} =					# calculate beam velocities
-				velInstrumentToBeam(\%dta,
-					&velEarthToInstrument(\%dta,$e,@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]}));
+				velEarthToBeam(\%dta,$e,@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]});
 				                                            
 			@{$dta{ENSEMBLE}[$e]->{BEAMPAIR_VELOCITY}[$b]} =				# calculate w12, w34
 				velBeamToBPEarth(\%dta,$e,@{$dta{ENSEMBLE}[$e]->{BEAM_VELOCITY}[$b]});