listBins
changeset 28 7c7da52363c2
parent 24 1a761865f839
child 31 b6ca27a1d19c
--- 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: Tue Jun 16 09:28:17 2015
+#                    dlm: Wed Jan  6 16:08:54 2016
 #                    (c) 2006 A.M. Thurnherr
-#                    uE-Info: 53 74 NIL 0 0 72 2 2 4 NIL ofnI
+#                    uE-Info: 153 28 NIL 0 0 72 10 2 4 NIL ofnI
 #======================================================================
 
 # Split data file into per-bin time series.
@@ -51,6 +51,10 @@
 #	Nov 24, 2014: - enabled -w always
 #	Mar 22, 2015: - replaced -f by -o (allowing for pipes)
 #	Jun 16, 2015: - BUG: velocity bias code did not respect bad velocities
+#	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)
 
 # General Notes:
 #	- everything (e.g. beams) is numbered from 1
@@ -79,10 +83,10 @@
 #	  data, where one of the beams performed much worse than the others
 
 require "getopts.pl";
-$0 =~ m{(.*/)[^/]+};
-require "$1RDI_BB_Read.pl";
-require "$1RDI_Coords.pl";
-require "$1RDI_Utils.pl";
+
+$ADCP_tools_minVersion = 1.4;
+($ADCP_TOOLS) = ($0 =~ m{(.*/)[^/]+});
+require "$ADCP_TOOLS/ADCP_tools_lib.pl";
 
 die("Usage: $0 [-r)ange <first_ens,last_ens>] [-R)enumber ensembles] " .
 			  "[-o)utput <redirection[>bin%d.raw]>] " .
@@ -92,9 +96,8 @@
 			  "[-P)itch/Roll <bias/bias>] [-B)eamvel <bias/bias/bias/bias>] " .
 		 	  "[require -4)-beam solutions] [-d)iscard <beam#>] " .
 		 	  "[-p)ct-good <min>] " .
-		 	  "[output -b)eam coordinates] " .
 			  "<RDI file>\n")
-	unless (&Getopts("4abB:d:M:o:p:r:P:RS:") && @ARGV == 1);
+	unless (&Getopts("4aB:d:M:o:p:r:P:RS:") && @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);
@@ -102,15 +105,12 @@
 die("$0: -4 and -d are mutually exclusive\n")
 	if ($opt_4 && defined($opt_d));
 
-die("$0: -p and -b are mutually exclusive\n")
-	if ($opt_b && defined($opt_p));
-
 $opt_p = 0 unless defined($opt_p);
 
 $RDI_Coords::minValidVels = 4 if ($opt_4);			# no 3-beam solutions
 
 print(STDERR "WARNING: magnetic declination not set!\n")
-	unless defined($opt_M) || defined($opt_b);
+	unless defined($opt_M);
 
 $opt_o = '>bin%d.raw' unless defined($opt_o);
 $ifn = $ARGV[0];
@@ -146,13 +146,14 @@
 	foreach my $k (keys(%P)) {
 		print(P "$k\{$P{$k}\} ");
 	}
-	if ($opt_b) {
-		printf(STDERR "%.0f%% ",100*$good_vels[$b]/($le-$fe+1));
-	} else {
+	if ($beamCoords) {
 		my($pct3b) = ($good_vels[$b] > 0) ? 100*$three_beam[$b]/$good_vels[$b] : nan;
-		printf(STDERR "%.0f%%/%.0f%% ",100*$good_vels[$b]/($le-$fe+1),$pct3b);
-		printf(P "pct_3_beam{%.0f} ",$pct3b);
+		printf(STDERR "%02d:%.0f%%/%.0f%% ",$b+1,100*$good_vels[$b]/($le-$fe+1),$pct3b);
+	} else {
+		printf(STDERR "%02d:%.0f%% ",$b+1,100*$good_vels[$b]/($le-$fe+1));
 	}
+
+	printf(P "pct_3_beam{%.0f} ",$pct3b);
 	printf(P "pct_good_vels{%.0f} ",100*$good_vels[$b]/($le-$fe+1));
 	printf(P "bin{%d}",$b+1);
 	printf(P " soundspeed_correction{%s}",defined($opt_S) ? $opt_S : 'NONE!');
@@ -167,11 +168,11 @@
 			"{sig_heading} {sig_pitch} {sig_roll} " .
 			"{xmit_current} {xmit_voltage} " .
 			"{temp} " .
-			($opt_b ? "{v1} {v2} {v3} {v4} " : "{u} {v} {w} {err_vel} ") .
+			"{bv1} {bv2} {bv3} {bv4} {u} {v} {w} {err_vel} " .
 			"{v12} {w12} {v34} {w34} " .
 			"{corr1} {corr2} {corr3} {corr4} " .
 			"{amp1} {amp2} {amp3} {amp4} " .
-			($opt_b ? "{pcg1} {pcg2} {pcg3} {pcg4}" : "{3_beam} {min_pcg}")
+			"{pcg1} {pcg2} {pcg3} {pcg4} {3_beam} {min_pcg}"
 	);
 	print(P " {dz}") if ($variable_ssCorr);
 	print(P "\n");
@@ -196,6 +197,7 @@
 		print(P "$dta{ENSEMBLE}[$e]->{ADC_XMIT_VOLTAGE} ");
 		print(P "$dta{ENSEMBLE}[$e]->{TEMPERATURE} ");
 		if ($dta{ENSEMBLE}[$e]->{GOOD_VEL}[$b]) {
+			printf(P "%g %g %g %g ",@{$dta{ENSEMBLE}[$e]->{BEAM_VELOCITY}[$b]});
 			printf(P "%g ",$dta{ENSEMBLE}[$e]->{VELOCITY}[$b][0] * $ssCorr);
 			printf(P "%g ",$dta{ENSEMBLE}[$e]->{VELOCITY}[$b][1] * $ssCorr);
 			printf(P "%g ",$dta{ENSEMBLE}[$e]->{VELOCITY}[$b][2] * $ssCorr);
@@ -218,11 +220,14 @@
 		}
 		print(P "@{$dta{ENSEMBLE}[$e]->{CORRELATION}[$b]} ");
 		print(P "@{$dta{ENSEMBLE}[$e]->{ECHO_AMPLITUDE}[$b]} ");
-		if ($opt_b) {
+		if ($beamCoords) {
 			print(P "@{$dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b]} ");
-		} else {
 			printf(P "%d ",$dta{ENSEMBLE}[$e]->{THREE_BEAM}[$b]);
 			printf(P "%s ",min(@{$dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b]}));
+		} else {
+			print(P "nan nan nan nan ");
+			print(P "$dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][0] ");
+			print(P "nan ");
 		}
 		printf(P "%g ",$dz[$b]*$ssCorr) if ($variable_ssCorr);
 		print(P "\n");
@@ -244,10 +249,14 @@
 if ($dta{BEAM_COORDINATES}) {						# coords
 	$beamCoords = 1;
 } else {
-	die("$0: -b requires input in beam coordinates\n")
-		if ($opt_b);
-	die("$ifn: only beam and earth coordinates implemented so far\n")
+	die("$ifn: only beam and earth coordinates supported\n")
 		if (!$dta{EARTH_COORDINATES});
+	die("$ifn: -p requires beam-coordinate data\n")
+		if ($opt_p > 0);
+	die("$ifn: -d requires beam-coordinate data\n")
+		if defined($opt_d);
+	die("$ifn: -B requires beam-coordinate data\n")
+		if defined($opt_B);
 }
 
 for (my($b)=0; $b<$dta{N_BINS}; $b++) {				# calc dz
@@ -280,34 +289,63 @@
         if ($dta{ENSEMBLE}[$e]->{LOW_GAIN});
 
 	for (my($b)=0; $b<$dta{N_BINS}; $b++) {
-		$dta{ENSEMBLE}[$e]->{VELOCITY}[$b][0] -= $P{velbias_b1} if defined($dta{ENSEMBLE}[$e]->{VELOCITY}[$b][0]);
-		$dta{ENSEMBLE}[$e]->{VELOCITY}[$b][1] -= $P{velbias_b2} if defined($dta{ENSEMBLE}[$e]->{VELOCITY}[$b][1]);
-		$dta{ENSEMBLE}[$e]->{VELOCITY}[$b][2] -= $P{velbias_b3} if defined($dta{ENSEMBLE}[$e]->{VELOCITY}[$b][2]);
-		$dta{ENSEMBLE}[$e]->{VELOCITY}[$b][3] -= $P{velbias_b4} if defined($dta{ENSEMBLE}[$e]->{VELOCITY}[$b][3]);
-		if (defined($opt_d)) {
-			undef($dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][$opt_d-1]);
-			undef($dta{ENSEMBLE}[$e]->{VELOCITY}[$b][$opt_d-1]);
+		if ($beamCoords) {
+			for (my($i)=0; $i<4; $i++) {										# percent-good editing (-p)
+				if ($dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][$i] < $opt_p) {
+					undef($dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][$i]);
+					undef($dta{ENSEMBLE}[$e]->{VELOCITY}[$b][$i]);
+				}
+	        }
+
+			$dta{ENSEMBLE}[$e]->{VELOCITY}[$b][0] -= $P{velbias_b1}			# beam-velocity biases (-B)
+				if defined($dta{ENSEMBLE}[$e]->{VELOCITY}[$b][0]);
+			$dta{ENSEMBLE}[$e]->{VELOCITY}[$b][1] -= $P{velbias_b2}
+				if defined($dta{ENSEMBLE}[$e]->{VELOCITY}[$b][1]);
+			$dta{ENSEMBLE}[$e]->{VELOCITY}[$b][2] -= $P{velbias_b3}
+				if defined($dta{ENSEMBLE}[$e]->{VELOCITY}[$b][2]);
+			$dta{ENSEMBLE}[$e]->{VELOCITY}[$b][3] -= $P{velbias_b4}
+				if defined($dta{ENSEMBLE}[$e]->{VELOCITY}[$b][3]);
+
+			if (defined($opt_d)) {											# discard data from given beam (-d)
+				undef($dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][$opt_d-1]);
+				undef($dta{ENSEMBLE}[$e]->{VELOCITY}[$b][$opt_d-1]);
+			}
+
+			@{$dta{ENSEMBLE}[$e]->{BEAM_VELOCITY}[$b]} =					# save beam velocities
+				@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]};
+
+			@{$dta{ENSEMBLE}[$e]->{BEAMPAIR_VELOCITY}[$b]} =				# calculate w12, w34
+				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]})
+			  	);
+			$dta{ENSEMBLE}[$e]->{THREE_BEAM}[$b] = $RDI_Coords::threeBeamFlag;
+			$three_beam[$b] += $RDI_Coords::threeBeamFlag;
+
+			unless (defined($dta{ENSEMBLE}[$e]->{VELOCITY}[$b][0])) {
+				undef(@{$dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b]});			# not sure when this can happen
+				next;
+			}
+		} else { 															# Earth coordinates
+			@{$dta{ENSEMBLE}[$e]->{BEAM_VELOCITY}[$b]} =					# calculate beam velocities
+				velInstrumentToBeam(\%dta,
+					&velEarthToInstrument(\%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]});
+
+			@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]} = 						# correct for heading bias
+				velApplyHdgBias(\%dta,$e,@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]});
+
+			unless (defined($dta{ENSEMBLE}[$e]->{VELOCITY}[$b][0])) {
+				undef(@{$dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b]});			# not sure when/if this can happen
+				next;
+			}
+
 		}
-		@{$dta{ENSEMBLE}[$e]->{BEAMPAIR_VELOCITY}[$b]} =
-			velBeamToBPEarth(\%dta,$e,@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]});
-		for (my($i)=0; $i<4; $i++) {
-			if ($dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][$i] < $opt_p) {
-				undef($dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b][$i]);
-				undef($dta{ENSEMBLE}[$e]->{VELOCITY}[$b][$i]);
-			}
-        }
-		@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]} = $beamCoords
-			? velInstrumentToEarth(\%dta,$e,
-				  velBeamToInstrument(\%dta,@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]})
-			  )
-			: velApplyHdgBias(\%dta,$e,@{$dta{ENSEMBLE}[$e]->{VELOCITY}[$b]})
-				unless ($opt_b);
-		unless (defined($dta{ENSEMBLE}[$e]->{VELOCITY}[$b][0])) {
-			undef(@{$dta{ENSEMBLE}[$e]->{PERCENT_GOOD}[$b]});
-			next;
-		}
-		$dta{ENSEMBLE}[$e]->{THREE_BEAM}[$b] = $RDI_Coords::threeBeamFlag;
-		$three_beam[$b] += $RDI_Coords::threeBeamFlag;
+
 		$dta{ENSEMBLE}[$e]->{GOOD_VEL}[$b] = 1;
 		$good_vels[$b]++; 
 		$lastGoodBin = $b if ($b > $lastGoodBin);
@@ -329,15 +367,15 @@
 print( STDERR "Start      : $dta{ENSEMBLE}[$fe]->{DATE} $dta{ENSEMBLE}[$fe]->{TIME}\n");
 print( STDERR "End        : $dta{ENSEMBLE}[$le]->{DATE} $dta{ENSEMBLE}[$le]->{TIME}\n");
 printf(STDERR "Bins       : %d-%d\n",$firstBin+1,$lastBin+1);
-if ($opt_b) {
-	print(STDERR "Good       : ");
-} else {
+if ($beamCoords) {
 	printf(STDERR "3-Beam     : %d %d %d %d\n",$RDI_Coords::threeBeam_1,
 											   $RDI_Coords::threeBeam_2,
 											   $RDI_Coords::threeBeam_3,
 											   $RDI_Coords::threeBeam_4);
 	
 	print(STDERR "Good/3-beam: ");
+} else {
+	print(STDERR "Good       : ");
 }
 for ($b=$firstBin; $b<=$lastBin; $b++) {				# generate output
 	dumpBin($b,$fe,$le);