# HG changeset patch # User A.M. Thurnherr # Date 1285633511 14400 # Node ID a3b6a908dec54db4a899e32687a16e5b7d07d4bf # Parent 229a0d72d2abfa76ec9a15857d3edbc1a316a5fe after P302 cruise diff --git a/RDI_BB_Read.pl b/RDI_BB_Read.pl --- a/RDI_BB_Read.pl +++ b/RDI_BB_Read.pl @@ -1,9 +1,9 @@ #====================================================================== # R D I _ B B _ R E A D . P L # doc: Sat Jan 18 14:54:43 2003 -# dlm: Wed Jun 4 09:43:15 2008 +# dlm: Sun Aug 15 16:35:54 2010 # (c) 2003 A.M. Thurnherr -# uE-Info: 44 25 NIL 0 0 72 0 2 4 NIL ofnI +# uE-Info: 47 72 NIL 0 0 72 0 2 4 NIL ofnI #====================================================================== # Read RDI BroadBand Binary Data Files (*.[0-9][0-9][0-9]) @@ -42,6 +42,9 @@ # Sep 18, 2007: - modified readHeader() readDta() WBRhdr() WBRens() to # conserve memory (no large arrays as retvals) # Jun 4, 2008: - BUG: BB150 code was not considered on Sep 18, 2007 +# Aug 15, 2010: - downgraded "unexpected number of data types" from error to warning +# - BUG: WBRcfn had not been set correctly +# - modified to allow processing files without time info # FIRMWARE VERSIONS: # It appears that different firmware versions generate different file @@ -294,8 +297,8 @@ = unpack('CCvCC',$buf); $hid == 0x7f || die(sprintf($FmtErr,$WBRcfn,"Header",$hid,0)); $did == 0x7f || die(sprintf($FmtErr,$WBRcfn,"Data Source",$did,0)); - die(sprintf("\n$WBRcfn: WARNING: unexpected number of data types (%d)\n", - $dta->{NUMBER_OF_DATA_TYPES})) + printf(STDERR "\n$WBRcfn: WARNING: unexpected number of data types (%d)\n", + $dta->{NUMBER_OF_DATA_TYPES}) unless ($dta->{NUMBER_OF_DATA_TYPES} == 6 || $dta->{NUMBER_OF_DATA_TYPES} == 7); $dta->{BT_PRESENT} = ($dta->{NUMBER_OF_DATA_TYPES} == 7); @@ -600,19 +603,27 @@ = &dayNo(${$E}[$ens]->{YEAR},${$E}[$ens]->{MONTH},${$E}[$ens]->{DAY}, ${$E}[$ens]->{HOUR},${$E}[$ens]->{MINUTE},${$E}[$ens]->{SECONDS}); - ${$E}[$ens]->{UNIX_TIME} - = timegm(0,${$E}[$ens]->{MINUTE}, - ${$E}[$ens]->{HOUR}, - ${$E}[$ens]->{DAY}, - ${$E}[$ens]->{MONTH}-1, # NB!!! - ${$E}[$ens]->{YEAR}) - + ${$E}[$ens]->{SECONDS}; - - $dayStart = timegm(0,0,0,${$E}[$ens]->{DAY}, - ${$E}[$ens]->{MONTH}-1, # NB!!! - ${$E}[$ens]->{YEAR}) - unless defined($dayStart); - ${$E}[$ens]->{SECNO} = ${$E}[$ens]->{UNIX_TIME} - $dayStart; + # when analyzing an STA file from an OS75 SADCP (Poseidion), + # 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 { + ${$E}[$ens]->{UNIX_TIME} + = timegm(0,${$E}[$ens]->{MINUTE}, + ${$E}[$ens]->{HOUR}, + ${$E}[$ens]->{DAY}, + ${$E}[$ens]->{MONTH}-1, # timegm jan==0!!! + ${$E}[$ens]->{YEAR}) + + ${$E}[$ens]->{SECONDS}; + + $dayStart = timegm(0,0,0,${$E}[$ens]->{DAY}, + ${$E}[$ens]->{MONTH}-1, + ${$E}[$ens]->{YEAR}) + unless defined($dayStart); + ${$E}[$ens]->{SECNO} = ${$E}[$ens]->{UNIX_TIME} - $dayStart; + } seek(WBRF,$start_ens+$WBRofs[0]+4,0) # System Config / Fixed Leader || die("$WBRcfn: $!"); @@ -774,14 +785,16 @@ sub readHeader(@) { - my($WBRcfn,$dta) = @_; + my($fn,$dta) = @_; + $WBRcfn = $fn; open(WBRF,$WBRcfn) || die("$WBRcfn: $!\n"); WBRhdr($dta); } sub readData(@) { - my($WBRcfn,$dta) = @_; + my($fn,$dta) = @_; + $WBRcfn = $fn; open(WBRF,$WBRcfn) || die("$WBRcfn: $!\n"); WBRhdr($dta); WBRens($dta->{N_BINS},$dta->{ENSEMBLE_BYTES}, diff --git a/WorkhorseUtils.pl b/WorkhorseUtils.pl deleted file mode 100644 --- a/WorkhorseUtils.pl +++ /dev/null @@ -1,162 +0,0 @@ -#====================================================================== -# W O R K H O R S E U T I L S . P L -# doc: Wed Feb 12 10:21:32 2003 -# dlm: Sun May 23 16:32:37 2010 -# (c) 2003 A.M. Thurnherr -# uE-Info: 142 42 NIL 0 0 72 2 2 4 NIL ofnI -#====================================================================== - -# miscellaneous Workhorse-specific utilities - -# History: -# Feb 12, 2003: - created -# Feb 14, 2003: - added check for short (bad) data files -# Feb 26, 2004: - added Earth-coordinate support -# - added ensure_BT_RANGE() -# Mar 17, 2004: - set bad BT ranges to undef in ensure_BT_RANGE -# - calc mean/stddev in ensure_BT_RANGE -# Mar 20, 2004: - BUG: find_seabed() could bomb when not enough -# bottom guesses passed the mode_width test -# Mar 30, 2004: - added &soundSpeed() -# Nov 8, 2005: - WATER_DEPTH => Z_BT -# May 23, 2010: - Z* -> DEPTH* - -use strict; - -#====================================================================== -# fake_BT_RANGE(dta ptr) -#====================================================================== - -# During cruise NBP0204 it was found that one of Visbeck's RDIs consistently -# returns zero as the bottom-track range, even thought the BT velocities -# are good. This functions calculates the ranges if they are missing. - -sub cBTR($$$) -{ - my($d,$e,$b) = @_; - my($maxamp) = -9e99; - my($maxi); - - for (my($i)=0; $i<$d->{N_BINS}; $i++) { - next unless ($d->{ENSEMBLE}[$e]->{ECHO_AMPLITUDE}[$i][$b] > $maxamp); - $maxamp = $d->{ENSEMBLE}[$e]->{ECHO_AMPLITUDE}[$i][$b]; - $maxi = $i; - } - $d->{ENSEMBLE}[$e]->{BT_RANGE}[$b] = - $d->{DISTANCE_TO_BIN1_CENTER} + $maxi * $d->{BIN_LENGTH}; -} - -sub ensure_BT_RANGE($) -{ - my($d) = @_; - for (my($e)=0; $e<=$#{$d->{ENSEMBLE}}; $e++) { - my($sum) = my($n) = 0; - if (defined($d->{ENSEMBLE}[$e]->{BT_VELOCITY}[0])) { - for (my($b)=0; $b<4; $b++) { - cBTR($d,$e,$b) - unless defined($d->{ENSEMBLE}[$e]->{BT_RANGE}[$b]); - $sum += $d->{ENSEMBLE}[$e]->{BT_RANGE}[$b]; $n++; - } - } else { - for (my($b)=0; $b<4; $b++) { - $d->{ENSEMBLE}[$e]->{BT_RANGE}[$b] = undef; - } - } - $d->{ENSEMBLE}[$e]->{BT_MEAN_RANGE} = $sum/$n if ($n == 4); - } -} - -#====================================================================== -# (seabed depth, stddev) = find_seabed(dta ptr, btm ensno, coord flag) -#====================================================================== - -# NOTE FOR YOYOS: -# - this routine only detects the BT around the depeest depth! -# - this is on purpose, because it is used for [listBT] - -# This is a PAIN: -# if the instrument is too close to the bottom, the BT_RANGE -# readings are all out; the only solution is to have a sufficiently -# large search width (around the max(depth) ensemble) so that -# a part of the up (oops, I forgot to finish this comment one year -# ago during A0304 and now I don't understand it any more :-) - -my($search_width) = 200; # # of ensembles around bottom to search -my($mode_width) = 10; # max range of bottom around mode -my($z_offset) = 10000; # shift z to ensure +ve array indices - -sub find_seabed($$$) -{ - my($d,$be,$beamCoords) = @_; - my($i,$dd,$sd,$nd); - my(@guesses); - - return undef unless ($be-$search_width >= 0 && - $be+$search_width <= $#{$d->{ENSEMBLE}}); - for ($i=$be-$search_width; $i<=$be+$search_width; $i++) { - next unless (defined($d->{ENSEMBLE}[$i]->{BT_RANGE}[0]) && - defined($d->{ENSEMBLE}[$i]->{BT_RANGE}[1]) && - defined($d->{ENSEMBLE}[$i]->{BT_RANGE}[2]) && - defined($d->{ENSEMBLE}[$i]->{BT_RANGE}[3])); - my(@BT) = $beamCoords - ? velInstrumentToEarth($d,$i, - velBeamToInstrument($d, - @{$d->{ENSEMBLE}[$i]->{BT_VELOCITY}})) - : velApplyHdgBias($d,$i,@{$d->{ENSEMBLE}[$i]->{BT_VELOCITY}}); - next unless (abs($BT[3]) < 0.05); - $d->{ENSEMBLE}[$i]->{DEPTH_BT} = - $d->{ENSEMBLE}[$i]->{DEPTH_BT} = - $d->{ENSEMBLE}[$i]->{BT_RANGE}[0]/4 + - $d->{ENSEMBLE}[$i]->{BT_RANGE}[1]/4 + - $d->{ENSEMBLE}[$i]->{BT_RANGE}[2]/4 + - $d->{ENSEMBLE}[$i]->{BT_RANGE}[3]/4; - $d->{ENSEMBLE}[$i]->{DEPTH_BT} *= -1 - if ($d->{ENSEMBLE}[$i]->{XDUCER_FACING_UP}); - $d->{ENSEMBLE}[$i]->{DEPTH_BT} += $d->{ENSEMBLE}[$i]->{DEPTH}; - $guesses[int($d->{ENSEMBLE}[$i]->{DEPTH_BT})+$z_offset]++; - $nd++; - } - return undef unless ($nd>5); - - my($mode,$nmax); - for ($i=0; $i<=$#guesses; $i++) { # find mode - $nmax=$guesses[$i],$mode=$i-$z_offset - if ($guesses[$i] > $nmax); - } - - $nd = 0; - for ($i=$be-$search_width; $i<=$be+$search_width; $i++) { - next unless defined($d->{ENSEMBLE}[$i]->{DEPTH_BT}); - if (abs($d->{ENSEMBLE}[$i]->{DEPTH_BT}-$mode) <= $mode_width) { - $dd += $d->{ENSEMBLE}[$i]->{DEPTH_BT}; - $nd++; - } else { - $d->{ENSEMBLE}[$i]->{DEPTH_BT} = undef; - } - } - return undef unless ($nd >= 2); - - $dd /= $nd; - for ($i=$be-$search_width; $i<=$be+$search_width; $i++) { - next unless defined($d->{ENSEMBLE}[$i]->{DEPTH_BT}); - $sd += ($d->{ENSEMBLE}[$i]->{DEPTH_BT}-$dd)**2; - } - - return ($dd, sqrt($sd/($nd-1))); -} - -#====================================================================== -# c = soundSpeed() -#====================================================================== - -# Taken from the RDI BroadBand primer manual. The reference given there -# is Urick (1983). - -sub soundSpeed($$$) -{ - my($salin,$temp,$depth) = @_; - return 1449.2 + 4.6*$temp -0.055*$temp**2 + 0.00029*$temp**3 + - (1.34 - 0.01*$temp) * ($salin - 35) + 0.016*$depth; -} - -1; 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: Sun May 23 16:33:33 2010 +# dlm: Sat Aug 21 23:24:42 2010 # (c) 2003 A.M. Thurnherr -# uE-Info: 527 42 NIL 0 0 72 11 2 4 NIL ofnI +# uE-Info: 123 56 NIL 0 0 72 11 2 4 NIL ofnI #====================================================================== # Extract Bottom-Track Data @@ -120,7 +120,7 @@ if ($dta{BT_MIN_EVAL_AMPLITUDE} > $opt_a); die("$0: maximum instrument BT error velocity ($dta{BT_MAX_ERROR_VELOCITY}) " . "too small for selected criterion (-e $opt_e) --- use -f to override\n") - if ($dta{BT_MAX_ERROR_VELOCITY} < $opt_e); + if (defined($dta{BT_MAX_ERROR_VELOCITY}) && $dta{BT_MAX_ERROR_VELOCITY} < $opt_e); } require $opt_F if defined($opt_F); # load filter 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: Sat May 23 16:34:30 2009 +# dlm: Sun Aug 22 22:40:32 2010 # (c) 2006 A.M. Thurnherr -# uE-Info: 269 0 NIL 0 0 72 2 2 4 NIL ofnI +# uE-Info: 45 28 NIL 0 0 72 2 2 4 NIL ofnI #====================================================================== # Split data file into per-bin time series. @@ -42,6 +42,7 @@ # - -P)itchRoll # May 22, 2009: - added -B) # May 23, 2009: - adapted to changed beampair-velocity fun name +# Aug 22, 2010: - added -R # General Notes: # - everything (e.g. beams) is numbered from 1 @@ -68,7 +69,7 @@ require "$1RDI_Coords.pl"; require "$1RDI_Utils.pl"; -die("Usage: $0 [-r)ange ] " . +die("Usage: $0 [-r)ange ] [-R)enumber ensembles] " . "[output -f)ile ] " . "[-a)ll ens (not just those with good vels)] " . "[-M)agnetic ] " . @@ -78,7 +79,7 @@ "[-%)good ] " . "[output -b)eam coordinates] [output two separate -w) estimates] " . "\n") - unless (&Getopts("4abB:d:f:M:p:r:P:S:w") && @ARGV == 1); + unless (&Getopts("4abB:d:f:M:p:r:P:RS:w") && @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); @@ -238,7 +239,8 @@ } $lastGoodBin = 0; -for ($e=0; $e<=$#{$dta{ENSEMBLE}}; $e++) { # check/transform velocities +for ($e=0; $e<=$#{$dta{ENSEMBLE}}; $e++) { # check/transform velocities + $dta{ENSEMBLE}[$e]->{NUMBER} = $e+1 if ($opt_R); # renumber ensembles next if (defined($first_ens) && $dta{ENSEMBLE}[$e]->{NUMBER} < $first_ens); diff --git a/listEns b/listEns --- a/listEns +++ b/listEns @@ -2,9 +2,9 @@ #====================================================================== # L I S T E N S # doc: Sat Jan 18 18:41:49 2003 -# dlm: Thu Jul 30 17:42:09 2009 +# dlm: Sun Aug 15 10:29:12 2010 # (c) 2003 A.M. Thurnherr -# uE-Info: 247 47 NIL 0 0 72 2 2 4 NIL ofnI +# uE-Info: 34 35 NIL 0 0 72 2 2 4 NIL ofnI #====================================================================== # Print useful info from the ensemble list or dump ensembles to @@ -31,6 +31,7 @@ # - restructured for simplicity # Mar 2, 2009: - added # of valid bin-1 vels to non-ANTS output # Jul 30, 2009: - NaN => nan +# Aug 15, 2010: - BUG: usage typo # Notes: # - -E outputs data in earth coordinates @@ -45,7 +46,7 @@ die("Usage: $0 [-A)nts] [-Q)uiet (errcheck only)] " . "[-f)ields <[name=]FIELD[,...]>] " . "[write -E)nsemples [-M)agnetic ] [min -p)ercent-good <#>]] " . - "[-r)ange ] [in-w)ater ensembles only]" . + "[-r)ange ] [in-w)ater ensembles only] " . "\n") unless (&Getopts("AE:f:M:p:Qr:w") && $#ARGV >= 0); diff --git a/splitRDI b/splitRDI new file mode 100755 --- /dev/null +++ b/splitRDI @@ -0,0 +1,64 @@ +#!/usr/bin/perl +#====================================================================== +# S P L I T R D I +# doc: Sat Aug 21 22:20:27 2010 +# dlm: Sat Aug 21 23:18:41 2010 +# (c) 2010 A.M. Thurnherr +# uE-Info: 51 44 NIL 0 0 72 2 2 4 NIL ofnI +#====================================================================== + +# split RDI files based on list of ensemble numbers (e.g. from yoyo -t) + +# HISTORY: +# Aug 21, 2010: - created + +# NOTES: +# - it is assumed that the input file begins with ensemble #1 +# - input file extension 000 is assumed +# - turning-point ensembles are written to preceding profile, +# for compatibility with [yoyo] + +$0 =~ m{(.*/)[^/]+}; +require "$1RDI_BB_Read.pl"; +use Getopt::Std; + +die("Usage: $0 " . + "[out-file -b)asename ] " . + "[out-file -n) ] " . + " \n") + unless (&getopts('b:n:') && @ARGV>=3); + +chomp($opt_b = `basename $ARGV[0] .000`) # basename + unless defined($opt_b); +$opt_n = 2 # number of digits in profile # + unless defined($opt_n); + +readHeader($ARGV[0],\%hdr); shift; # get length of ensembles +$ens_len = $hdr{ENSEMBLE_BYTES} + 2; + +$first_ens = $ARGV[0]+1; shift; # initialize loop +$last_ens = $ARGV[0]; shift; +$cnr = 0; + +do { # split data + sysseek(WBRF,($first_ens-2)*$ens_len,0) || + die("$WBRcfn: $!"); + $last_ens++ unless defined($ARGV[0]); + $nBytes = ($last_ens-$first_ens+1) * $ens_len; + sysread(WBRF,$buf,$nBytes) == $nBytes || + die("$WBRcfn: file truncated"); + + $fn = $opt_b . sprintf("_%0${opt_n}d.000",$cnr++); + open(F,">$fn") || die("$fn: $!\n"); + syswrite(F,$buf,$nBytes) == $nBytes || + die("$fn: $!\n"); + close(F); + printf(STDERR "$fn: %d ensembles ($nBytes bytes)\n",$last_ens-$first_ens+1); + + $first_ens = $last_ens+1; + $last_ens = $ARGV[0]; shift; +} while defined($last_ens); + +exit(0); + +