--- a/HISTORY Fri Aug 05 11:02:51 2016 -0400
+++ b/HISTORY Thu Mar 16 11:53:27 2017 -0400
@@ -1,9 +1,9 @@
======================================================================
H I S T O R Y
doc: Mon Oct 12 16:09:24 2015
- dlm: Wed Jun 8 22:07:57 2016
+ dlm: Thu Mar 16 11:53:07 2017
(c) 2015 A.M. Thurnherr
- uE-Info: 182 43 NIL 0 0 72 3 2 4 NIL ofnI
+ uE-Info: 260 19 NIL 0 0 72 3 2 4 NIL ofnI
======================================================================
----------------------------------------------------------------------
@@ -175,8 +175,86 @@
- updated version to V1.3beta2 [version.pl] [.hg/hgrc]
- udated ANTS_tools lib to V1.7 (beam interpolation)
- adapted [LADCP_w_ocean] to beam interpolation
+ - minor improvement to [LADCP_w_postproc]
+ - improved [plot_wprof.pl]
-...
+ Jun 1, 2016
+ - improvements to [LADCP_w_ocean]
+ - added [default_output.pl]
+ - added [plot_residuals12.pl] [plot_residuals34.pl]
+
+ Jun 2, 2016:
+ - minor improvement and bug fix in [LADCP_w_ocean]
+
+ Jun 3, 2016:
+ - minor bug fix in [LADCP_w_ocean]
+
+ Jun 6, 2016:
+ - minor improvement in [LADCP_w_ocean] [defaults.pl] [edit_data.pl]
Jun 8, 2016:
- removed plot_attitude_biases_w.pl
+ - slight improvement to [plot_attitude_residuals.pl]
+
+ Jun 11, 2016:
+ - began debugging w12 & w34 for Earth-coord data [LADCP_w_ocean]
+
+ Jul 7, 2016:
+ - major BUG in [LADCP_w_ocean] (beam-pair velocities for Earth
+ coord data)
+
+ Jul 12, 2016:
+ - docu in [defaults.pl]
+
+ Jul 29, 2016:
+ - minor plotting bug in [LADCP_w_CTD]
+
+ Jul 31, 2016:
+ - minor bug in [LADCP_w_ocean] [defaults.pl]
+
+ Aug 5, 2016:
+ - committed version found on whoosher after repair
+ - manually uploaded from ECOGIG cruise laptop:
+ - [LADCP_w_ocean] changes since Jun 11, 2016
+ - [defaults.pl] changes since Jun 2, 2016
+ - [LADCP_w_CTD] changes since May 26, 2016
+ - updated [version.pl] to require ANTSlib V6.7
+ - updated HISTORY
+
+ Aug 16, 2016:
+ - [LADCP_VKE] increased -l default 1.2e-7 based on UK2.5 SR1
+ repeat section
+
+ Dec 22, 2016:
+ - [LADCP_w_ocean] moved $opt_p to [defaults.pl]
+
+ Dec 23, 2016:
+ - [LADCP_w_ocean] minor bug
+
+ Sep 1, 2016:
+ - [LADCP_VKE] changed -l to mean epsilon, and increased value to
+ 1e-10
+
+ Mar 6, 2017:
+ - [LADCP_w_ocean] minor bug
+ - [time_lag.pl] minor bug
+
+ Mar 7, 2017:
+ - added time lines to [plot_time_lags.pl]
+
+ Mar 9, 2017:
+ - tightened timelag editing condition in [time_lag.pl]
+ - updated [HISTORY]
+
+ Mar 12, 2017:
+ - adapted to antslib V6.8 [version.pl]
+ - adapted ADCP_tools to V1.9
+ - increased -o default from 20 to 40m
+ - updated to V1.3
+
+ Mar 15, 2017:
+ - added [loadANTS.pl] for V1.3
+ - updated howto
+
+ Mar 16, 2017:
+ - published
--- a/LADCP_VKE Fri Aug 05 11:02:51 2016 -0400
+++ b/LADCP_VKE Thu Mar 16 11:53:27 2017 -0400
@@ -2,9 +2,9 @@
#======================================================================
# L A D C P _ V K E
# doc: Tue Oct 14 11:05:16 2014
-# dlm: Wed Apr 27 19:03:02 2016
+# dlm: Tue Mar 14 17:38:00 2017
# (c) 2012 A.M. Thurnherr
-# uE-Info: 99 45 NIL 0 0 72 0 2 4 NIL ofnI
+# uE-Info: 155 45 NIL 0 0 72 0 2 4 NIL ofnI
#======================================================================
$antsSummary = 'calculate VKE from LADCP-derived vertical-velocity profiles';
@@ -97,6 +97,12 @@
# Apr 16, 2016: - assertion failed on I98S#153 => modify code to allow nans in LADCP_wspec
# output
# Apr 27, 2016: - cosmetics (usage message)
+# Aug 16, 2016: - increased -l default to 1.2e-7 based on UK2.5 SR1 repeat section
+# Sep 1, 2016: - changed -l to mean epsilon, and increased value to 1e-10
+# - added %eps.minlim
+# Mar 13, 2017: - added -a)mbient <eps>
+# Mar 14, 2017: - disabled -a) by default, because -a 0 is clearly bad and
+# I have no evidence yet that -a something is better than -l 0
($ANTSLIB) = (`which ANTSlib` =~ m{^(.*)/[^/]*$});
($WCALC) = ($0 =~ m{^(.*)/[^/]*$});
@@ -119,7 +125,8 @@
#my($c) = 0.0215; # Thurnherr et al. (GRL 2015)
my($c) = 0.026; # increased by 20% for V1.2beta7
$opt_q = 3; # Equatorial band: little more than a guess based on 2015 P16N
-$opt_l = 5e-8; # low-p0 tail; visually estimated from Fig. 2 of Thurnherr et al. (GRL 2015)
+$opt_l = 5e-11; # [W/kg]; based on DIMES data sets, method works for eps >= 2e-10; 5e-11 allows for scatter in p0
+$opt_a = nan; # assume background dissipation for samples that pass the tests but have eps below -l
$opt_z = 1; # number of w_ocean samples to require
$opt_o = 0; # remove mean before calculating spectra
$opt_s = 150; # surface layer to exclude from spectra
@@ -130,22 +137,22 @@
# Usage
#----------------------------------------------------------------------
-&antsUsage('bc:de:f:g:i:k:l:mno:p:q:r:s:tuw:x:z:',0,
+&antsUsage('a:bc:de:f:g:i:k:l:mno:p:q:r:s:tuw:x:z:',0,
"[poly-o)rder <n[$opt_o]> to de-mean data; -1 to disable>] [apply cosine-t)aper]",
'[-d)own/-u)pcast-only] [exclude -b)ottom window]', # LADCP_wspec options
"[-s)urface <layer depth to exclude[${opt_s}m]>",
"[-g)ap <max depth layer to fill with interpolation[${opt_g}m]>]",
'[-w)indow <power-of-two input-records[32]>]',
"[shortwave -c)utoff <kz or lambda[${opt_c_default}m]>]", # LADCP_VKE options
- "[e-q)uatorial cutoff <latitude[${opt_q}deg]>] [-l)ow-p0 <cutoff[${opt_l}m^2/s^2/(rad/m)]>]",
+ "[e-q)uatorial cutoff <latitude[${opt_q}deg]>]",
+ "[-l)ow-eps <cutoff[${opt_l} W/kg]>] [-a)mbient <eps[${opt_a} W/kg]>]",
"[-z) ignore velocities derived from fewer than <N[$opt_z]> samples]",
'[o-m)it spectral correction] [spectral-tilt-correction -r)ange <max[0m]>]',
-# '[-l) override ADCP bin <length>] [-a) override pulse <length>]',
"[-e)ps-parameterization <constant[${c}s^-0.5]>",
- '[-k)e dissipation <file:field>]',
- '[output -f)iles in <dir>]',
- '[output -i)ndividual spectra <basename>]',
- '[output -p)lot <ps-file>]',
+ '[include microstructure -k)e dissipation <file:field> in _VKE plot]',
+ '[write output -f)iles to <directory>]',
+ '[write output filed with -i)ndividual spectra <basename>]',
+ '[output -p)lot <ps-file[#_VKE.ps]>]',
'[file...]');
#----------------------------------------------------------------------------
@@ -172,11 +179,6 @@
open2(\*FROMCLD,\*TOCLD,"LADCP_wspec $opts") || # spawn sub-process
croak("LADCP_wspec $opts: $!\n");
print(TOCLD $antsOldHeaders); # feed already gobbled header
-# if (defined($opt_l) || defined($opt_a)) { # override bin size and/or pulse length
-# &antsAddParams('ADCP_bin_length',$opt_l) if defined($opt_l);
-# &antsAddParams('ADCP_pulse_length',$opt_a) if defined($opt_a);
-# print(TOCLD $antsCurParams);
-# }
for (my($r)=0; $r<@ants_; $r++) { # feed all .wprof records to LADCP_wspec
print(TOCLD "@{$ants_[$r]}\n");
}
@@ -223,8 +225,6 @@
open2(\*FROMCLD,\*TOCLD,"LADCP_wspec $opts") || # process it
croak("LADCP_wspec $opts: $!\n");
print(TOCLD $antsOldHeaders);
-# print(TOCLD $antsCurParams)
-# if defined($opt_l) || defined($opt_a);
for (my($r)=0; $r<@ants_; $r++) {
print(TOCLD "@{$ants_[$r]}\n");
}
@@ -241,9 +241,6 @@
undef(@ants_);
} elsif (defined(fnrNoErr('pwrdens.0'))) {
-# croak("$0: -a, -l, -d, -u, -b, -w, -s meaningless when $0 used with spectral input\n")
-# if ($opt_d || $opt_u || $opt_b || defined($opt_w) ||
-# defined($opt_s) || defined($opt_g)) || defined($opt_a) || defined($opt_l);
croak("$0: -d, -u, -b, -w, -s meaningless when $0 used with spectral input\n")
if ($opt_d || $opt_u || $opt_b || defined($opt_w) || defined($opt_s) || defined($opt_g));
} else {
@@ -262,7 +259,8 @@
$n_input_files = 1 + @specbuf; # number of input files provided
&antsAddParams('LADCP_VKE::input_files.n',$n_input_files,
- 'LADCP_VKE::wsamp.min',$opt_z);
+ 'LADCP_VKE::wsamp.min',$opt_z,
+ 'LADCP_VKE::eps.minlim',$opt_l);
&antsFloatOpt(\$opt_e,$c); # default parameterization
&antsFloatOpt(\$opt_x,1); # spectral fit stddev scale factor
@@ -572,28 +570,20 @@
#
# Additional Empirical Filters:
# - latitude > 3deg guess based on Thurnherr et al., 2015, Gregg et al., 2003, 2015 GO-SHIP P16N
- # - p0 > 5x10^-8 m^s/s^2/(rad/m) based on Fig.2 in Thurnherr et al., 2015
+ # - eps >= 5e-11 W/kg based on DIMES data, cutting where errors become > factor 2
#-----------------------------------------------------------------------------------------------------
- if ($latM > $opt_q && # 1) not (too) equatorial
- $ants_[$r][$rmsf] <= 0.4 && # 2) rms spectra misfit <= 0.4 (as in Thurnherr et al., GRL 2015)
- $ants_[$r][$slpf]>=-3 && $ants_[$r][$slpf]<=-1 && # 3) slope consistent with -2 power law
- $ants_[$r][$rf] <= -0.4 && # 4) p and k_z are well correlated
- $ants_[$r][$wsf] >= $opt_z) { # 5) minimum # of samples
-# if (defined($opt_l)) { # 6) suppress low-p0 tail on -l
- $ants_[$r][$wepsf] = ($ants_[$r][$p0f] >= $opt_l)
- ? ($ants_[$r][$p0f] / $opt_e)**2 # Thurnherr et al. (GRL 2015)
- : nan; # suppress low-p0 samples
-# } else { # -l not set
-# my($lp0) = log10($ants_[$r][$p0f]);
-# if ($lp0 >= -6.5) { # Thurnherr et al. (GRL 2015)
-# $ants_[$r][$wepsf] = ($ants_[$r][$p0f] / $opt_e)**2;
-# } else { # low-p0 power-law fit
-# my($leps) = 0.56*$lp0 - 6.06;
-# $ants_[$r][$wepsf] = 10**$leps;
-# } # if log10(p0) > -6.5
-# } # if -l
- } else { # failed any of checks 1-5
+ if ($latM > $opt_q && # 1) not (too) equatorial
+ $ants_[$r][$rmsf] <= 0.4 && # 2) rms spectra misfit <= 0.4 (as in Thurnherr et al., GRL 2015)
+ $ants_[$r][$slpf]>=-3 && $ants_[$r][$slpf]<=-1 && # 3) slope consistent with -2 power law
+ $ants_[$r][$rf] <= -0.4 && # 4) p and k_z are well correlated
+ $ants_[$r][$wsf] >= $opt_z) { # 5) minimum # of samples
+ if (($ants_[$r][$p0f]/$opt_e)**2 >= $opt_l) { # level is above eps.minlim (-l)
+ $ants_[$r][$wepsf] = ($ants_[$r][$p0f] / $opt_e)**2; # => Thurnherr et al. (GRL 2015)
+ } else { # level is below eps.minlim
+ $ants_[$r][$wepsf] = $opt_a; # => set to arbitrary background value
+ }
+ } else { # failed any of checks 1-5 => nan
$ants_[$r][$wepsf] = nan;
}
--- a/LADCP_w_CTD Fri Aug 05 11:02:51 2016 -0400
+++ b/LADCP_w_CTD Thu Mar 16 11:53:27 2017 -0400
@@ -2,9 +2,9 @@
#======================================================================
# L A D C P _ W _ C T D
# doc: Mon Nov 3 17:34:19 2014
-# dlm: Thu May 26 20:42:15 2016
+# dlm: Fri Jul 29 11:20:47 2016
# (c) 2014 A.M. Thurnherr
-# uE-Info: 639 69 NIL 0 0 72 2 2 4 NIL ofnI
+# uE-Info: 69 48 NIL 0 0 72 2 2 4 NIL ofnI
#======================================================================
$antsSummary = 'pre-process SBE 9plus CTD data for LADCP_w';
@@ -62,6 +62,11 @@
# Mar 31, 2016: - changed version %PARAM
# May 26, 2016: - renamed w[.raw] to CTD_w[.raw]
# - added winch_w
+# Jul 29, 2016: - change w_CTD plot label to reflect sign convention
+
+# NOTES:
+# w_CTD is positive during the downcast to make the sign of the apparent
+# water velocity consistent with w_ocean
($ANTS) = (`which ANTSlib` =~ m{^(.*)/[^/]*$});
($WCALC) = ($0 =~ m{^(.*)/[^/]*$});
@@ -634,7 +639,7 @@
for ($r=0; $r<@w; $r++) {
printf(GMT "%f %f\n",$elapsed[$r]/60,$winch[$r]);
}
- GMT_psbasemap('-Bg60a30f5:"Elapsed Time [min]":/g1a1f0.1:"Vertical Package Velocity [ms@+-1@+]":WeSn');
+ GMT_psbasemap('-Bg60a30f5:"Elapsed Time [min]":/g1a1f0.1:"Downward Package Velocity [ms@+-1@+]":WeSn');
GMT_unitcoords();
GMT_pstext('-F+f14,Helvetica,coral+jBR'); print(GMT "0.98 0.96 dc\n");
GMT_pstext('-F+f14,Helvetica,SeaGreen+jBR'); print(GMT "0.98 0.93 uc\n");
Binary file LADCP_w_howto.pdf has changed
--- a/LADCP_w_ocean Fri Aug 05 11:02:51 2016 -0400
+++ b/LADCP_w_ocean Thu Mar 16 11:53:27 2017 -0400
@@ -2,14 +2,15 @@
#======================================================================
# L A D C P _ W _ O C E A N
# doc: Fri Dec 17 18:11:13 2010
-# dlm: Sat Jun 11 18:34:12 2016
+# dlm: Mon Mar 6 13:46:31 2017
# (c) 2010 A.M. Thurnherr
-# uE-Info: 272 65 NIL 0 0 72 2 2 4 NIL ofnI
+# uE-Info: 280 0 NIL 0 0 72 0 2 4 NIL ofnI
#======================================================================
# TODO:
# ! use instrument tilt in sidelobe editing
-# ! use instrument tilt in e.g. Sv
+# - use instrument tilt in e.g. calculation of Sv
+# - use instrument tilt in PPI editing
#
# - plots:
# - avoid over-plotting axis labels
@@ -270,6 +271,12 @@
# calculated
# Jun 6, 2016: - removed applyTiltCorrection()
# Jun 11, 2016: - BUG: w12, w34 were wrong for Earth-coord data
+# Jul 7, 2016: - BUG: finally fixed with velEarthToBPw() from [RDI_Coords.pl]
+# Jul 31, 2016: - BUG: -d did not work because it was handled in [defaults.pl]
+# Oct 16, 2016: - cosmetics
+# Dec 22, 2016: - moved $opt_p to [defaults.pl]
+# Dec 23, 2016: - BUG: -u did not set required variables to proceed
+# Mar 6, 2017: - BUG: division by zero when water-depth ~ max(CTD_depth)
# HISTORY END
# CTD REQUIREMENTS
@@ -394,7 +401,7 @@
&antsUsageError() unless (@ARGV==1 || @ARGV==2);
&antsCardOpt(\$opt_s,0); # skip # initial ensembles
-$opt_p = '+' unless defined($opt_p); # separate uc/dc time lagging
+#$opt_p = '+' unless defined($opt_p); # separate uc/dc time lagging
&antsFloatOpt(\$opt_a,1); # pressure acceleration correction
@@ -429,7 +436,6 @@
$processing_options .= " -i $opt_i" if defined($opt_i);
$processing_options .= ' -l' if defined($opt_l);
$processing_options .= ' -q' if defined($opt_q);
-$processing_options .= ' -d' if defined($opt_d);
if (defined($opt_x)) { # eval cmd-line expression to override anything
$processing_options .= " -x $opt_x";
@@ -440,6 +446,11 @@
$processing_options .= " -4";
$RDI_Coords::minValidVels = 4;
}
+
+if ($opt_d) {
+ $processing_options .= ' -d';
+ $RDI_Coords::binMapping = 'none';
+}
if (defined($opt_r)) { # residuals filters
$processing_options .= " -r $opt_r";
@@ -734,10 +745,9 @@
$LADCP{ENSEMBLE}[$ens]->{V12}[$bin] = $LADCP{ENSEMBLE}[$ens]->{V34}[$bin] = nan;
- # for 3-beam solutions, the following code sets w12 and w34 equal to w
- my($delta) = $LADCP{ENSEMBLE}[$ens]->{ERRVEL}[$bin] / sqrt(2) / tan(rad($LADCP{BEAM_ANGLE}));
- $LADCP{ENSEMBLE}[$ens]->{W12}[$bin] = $LADCP{ENSEMBLE}[$ens]->{W}[$bin] + $delta/2;
- $LADCP{ENSEMBLE}[$ens]->{W34}[$bin] = $LADCP{ENSEMBLE}[$ens]->{W}[$bin] - $delta/2;
+ # for 3-beam solutions, w12 = w34 = w (I think)
+ ($LADCP{ENSEMBLE}[$ens]->{W12}[$bin],$LADCP{ENSEMBLE}[$ens]->{W34}[$bin]) =
+ velEarthToBPw($LADCP,$ens,@{$LADCP{ENSEMBLE}[$ens]->{VELOCITY}[$bin]});
}
}
progress("\t$nvw valid velocities in bins $LADCP_firstBin-$LADCP_lastBin\n");
@@ -963,7 +973,6 @@
$sVelProf[$d] = $sVelProf[$d-1]
unless defined($sVelProf[$d]);
}
-
#-------------------
# Determine time lag
@@ -993,6 +1002,8 @@
if ($opt_u) {
progress("\tskipping time lagging (-u)\n");
+ $CTD_time_lag[0] = $CTD{TIME_LAG};
+ $CTD_tl_toEns[0] = $lastGoodEns;
} else {
#------------------------
@@ -1078,7 +1089,7 @@
if ($ens > $CTD_tl_toEns[$cli]) { # use correct lag piece
$cli++;
- die("assertion failed!\n\ttest: \$cli(=$cli) <= \$#CTD_time_lag\n")
+ die("assertion failed!\n\ttest: \$cli(=$cli) <= \$#CTD_time_lag(=$#CTD_time_lag) at ens=$ens\n")
unless ($cli <= $#CTD_time_lag);
$lag = $CTD_time_lag[$cli];
}
@@ -1271,17 +1282,23 @@
&antsAddParams('PPI_extend_upper_limit',$PPI_extend_upper_limit)
if numberp($PPI_extend_upper_limit);
progress("Editing data to remove PPI from seabed...\n");
- progress("\tConstructing depth-average soundspeed profile...\n");
- die("assertion failed") unless defined($water_depth);
- my($dz) = $water_depth - $#sVelProf; # $#sVelProf = max_depth(profile) in meters
- my($sum) = $dz * $sVelProf[$#sVelProf];
- $DASSprof[$#sVelProf] = $sum/$dz;
- for (my($d)=$#sVelProf-1; $d>=0; $d--) {
- die("assertion failed (d=$d, #sVelProf=$#sVelProf)") unless numberp($sVelProf[$d]);
- $sum += $sVelProf[$d];
+ progress("\tConstructing depth-average soundspeed profile...\n");
+ die("assertion failed") unless defined($water_depth);
+ my($dz) = $water_depth - $#sVelProf; # $#sVelProf = max_depth(profile) in meters
+ my($sum) = $dz * $sVelProf[$#sVelProf];
+ if ($dz == 0) { # water-depth <= max CTD depth
+ warning(1,"inferred water depth very close to max(CTD depth)\n");
+ $DASSprof[$#sVelProf] = $sVelProf[$#sVelProf];
+ } else {
+ die("assertion failed") unless defined($dz > 0);
+ $DASSprof[$#sVelProf] = $sum/$dz;
+ }
+ for (my($d)=$#sVelProf-1; $d>=0; $d--) {
+ die("assertion failed (d=$d, #sVelProf=$#sVelProf)") unless numberp($sVelProf[$d]);
+ $sum += $sVelProf[$d];
$dz++;
$DASSprof[$d] = $sum/$dz;
- }
+ }
($nvrm,$nerm) = editPPI($firstGoodEns,$lastGoodEns,$water_depth);
progress("\t$nvrm velocities from $nerm ensembles removed\n");
}
--- a/LADCP_wspec Fri Aug 05 11:02:51 2016 -0400
+++ b/LADCP_wspec Thu Mar 16 11:53:27 2017 -0400
@@ -2,9 +2,9 @@
#======================================================================
# L A D C P _ W S P E C
# doc: Thu Jun 11 12:02:49 2015
-# dlm: Thu Mar 31 11:08:44 2016
+# dlm: Sun Mar 12 12:01:59 2017
# (c) 2012 A.M. Thurnherr
-# uE-Info: 154 41 NIL 0 0 72 10 2 4 NIL ofnI
+# uE-Info: 48 17 NIL 0 0 72 10 2 4 NIL ofnI
#======================================================================
$antsSummary = 'calculate VKE window spectra from LADCP profiles';
@@ -31,8 +31,8 @@
# Mar 31, 2016: - changed version %PARAM
# - BUG: nspec was nan insted of 0
# - replaced wspec:: %PARAM-prefix with LADCP_wspec::
+# Mar 12, 2017: - removed ANTSBIN (which is not public)
-($ANTSBIN) = (`which list` =~ m{^(.*)/[^/]*$});
($ANTSLIB) = (`which ANTSlib` =~ m{^(.*)/[^/]*$});
($WCALC) = ($0 =~ m{^(.*)/[^/]*$});
$WCALC = '.' if ($WCALC eq '');
@@ -44,8 +44,8 @@
require "$ANTSLIB/fft.pl";
require "$ANTSLIB/lfit.pl";
require "$ANTSLIB/nrutil.pl";
-require "$ANTSBIN/.lsfit.poly";
-require "$ANTSBIN/.nminterp.linear";
+require "$ANTSLIB/.lsfit.poly";
+require "$ANTSLIB/.nminterp.linear";
&antsAddParams('LADCP_wspec::version',$VERSION);
#----------------------------------------------------------------------
--- a/defaults.pl Fri Aug 05 11:02:51 2016 -0400
+++ b/defaults.pl Thu Mar 16 11:53:27 2017 -0400
@@ -1,9 +1,9 @@
#======================================================================
# D E F A U L T S . P L
# doc: Tue Oct 11 17:11:21 2011
-# dlm: Mon Jun 6 22:07:19 2016
+# dlm: Sun Mar 12 12:53:44 2017
# (c) 2011 A.M. Thurnherr
-# uE-Info: 179 18 NIL 0 0 72 0 2 4 NIL ofnI
+# uE-Info: 459 36 NIL 0 0 72 0 2 4 NIL ofnI
#======================================================================
# HISTORY:
@@ -76,6 +76,12 @@
# May 28, 2016: - added delta-residual filter (-r)
# Jun 1, 2016: - added $plotting_level
# Jun 2, 2016: - added $tilt_correction_*
+# Jun 6, 2016: - removed $tilt_correction_*
+# Jul 12, 2016: - updated docu on acoustic backscatter
+# Jul 31, 2016: - BUG: -d to disable bin mapping did not work, because
+# this file is read before the options are parsed
+# Aug 5, 2016: - updated header
+# Dec 22, 2016: - added $opt_p
#======================================================================
# Output Log Files
@@ -117,8 +123,7 @@
# linear bin interpolation (which is better than bin mapping).
# This only has an effect for beam-coordinate data.
-$RDI_Coords::binMapping = $opt_d ? 'none' : 'linterp';
-
+# $RDI_Coords::binMapping = 'none';
# The following variables allow bias-correcting the attiude
# sensors.
@@ -325,6 +330,15 @@
&antsFloatOpt(\$opt_3,0.6);
+# Time lagging is carried out in "partial cast pieces". By default
+# the down- and upcasts are treated separately, which is accomplished
+# by setting $opt_p to '+'. In general, the $opt_p variable contains
+# a comma-separated list of elapsed splitting times (in minutes),
+# with '+' denoting the time of maximum depth.
+
+$opt_p = '+' unless defined($opt_p);
+
+
#======================================================================
# Acoustic Backscatter and Seabed Search
#======================================================================
@@ -332,11 +346,15 @@
# After applying the method of Deines (1999), an empirical correction
# for Sv is applied to the data. The following variable determines which
# bin is chosen to construct a reference profile for Sv. The bin number
-# is automatically increased if the selected bin does not contain valid
+# is automatically increased if its value is less than LADCP_firstBin (-b),
+# and also if the selected bin does not contain valid
# data, i.e. the default value of 1 ensures that the closest valid bin
# is used to construct the reference profile. The empirical correction
# causes artifacts every 100m. To disable the empirical
-# correction, undefine the following variable.
+# correction, undefine the following variable ($Sv_ref_bin = undef;)
+#
+# NOTE: Accoustic backscatter data in the reference bin are not
+# corrected beyond the method of Deines (1999).
$Sv_ref_bin = 1;
@@ -435,9 +453,12 @@
# The -o option sets the output grid resolution in meters. The following
-# sets the default value.
+# sets the default value. It was increased from 20m to 40m in Feb 2017
+# for V1.3 because this value was required for the DoMORE-2 data and
+# also improves the profiles from 2017 P18, as well as from a recent
+# GoM data set provided by J. Ochoa.
-&antsFloatOpt(\$opt_o,20);
+&antsFloatOpt(\$opt_o,40);
# The following variables limit the bins used to grid w_oean
--- /dev/null Thu Jan 01 00:00:00 1970 +0000
+++ b/loadANTS.m Thu Mar 16 11:53:27 2017 -0400
@@ -0,0 +1,111 @@
+%======================================================================
+% L O A D A N T S . M
+% doc: Thu Jul 21 22:53:21 2011
+% dlm: Thu Aug 4 10:54:24 2011
+% (c) 2011 A.M. Thurnherr
+% uE-Info: 18 71 NIL 0 0 72 2 2 4 NIL ofnI
+%======================================================================
+
+% NOTES:
+% - very restrictive subset of ANTS standard:
+% - no empty lines
+% - no in-record comments
+
+% HISTORY:
+% Jul 21, 2011: - began working on it
+% Jul 23, 2011: - made it work
+% Aug 4, 2011: - replaced "creator" by "Matlab_import"
+% - BUG: %PARAM names with . were not handled correctly
+
+function dta_struct = loadANTS(file_name)
+
+fid = fopen(file_name); % open file
+if fid < 0
+ error(sprintf('%s: not such file',file_name));
+end
+
+dta_struct = struct('Matlab_import',sprintf('loadANTS(''%s'')',file_name));
+
+l = fgetl(fid); % handle metadata
+while ischar(l) & regexp(l,'^#')
+ if regexp(l,'^#ANTS#ERROR#')
+ error(l);
+ elseif regexp(l,'^#ANTS#PARAMS#')
+ [tmp,tmp,ptk] = regexp(l,'([\w\.]+){([^}]*)}');
+ for i=1:length(ptk)
+ pname = matlabotomize(token1(i,l,ptk));
+ if sum(size(ptk{i})) < 4 % empty def
+ if isfield(dta_struct,pname)
+ dta_struct = rmfield(dta_struct,pname);
+ end
+ else
+ numval = str2double(token2(i,l,ptk));
+ if isfinite(numval) | strcmpi(token2(i,l,ptk),'nan')
+ dta_struct = setfield(dta_struct,pname,numval);
+ else
+ dta_struct = setfield(dta_struct,pname,token2(i,l,ptk));
+ end
+ end
+ end
+
+ elseif regexp(l,'^#ANTS#FIELDS#')
+ [tmp,tmp,ftk] = regexp(l,'{([^}]*)}');
+ fields = cell(1,length(ftk));
+ for i=1:length(ftk)
+ fields{i} = matlabotomize(token(i,l,ftk));
+ end
+ end % elseif
+ l = fgetl(fid);
+end % while
+
+if ischar(l) % not empty file
+ isnumeric = zeros(1,length(fields)); % determine data types
+ [tmp,tmp,tk] = regexp(l,'([^ \t]+)'); % split into tokens
+ for i=1:length(fields)
+ numval = str2double(token(i,l,tk));
+ if isfinite(numval) | strcmpi(token(i,l,tk),'nan')
+ isnumeric(i) = 1;
+ end
+ end
+end
+
+dta = cell(length(fields),1); % init data cell array
+
+while ischar(l) % loop through records
+ [tmp,tmp,tk] = regexp(l,'([^ \t]+)'); % split into tokens
+ for i=1:length(tk)
+ if isnumeric(i)
+ dta{i} = [dta{i} str2double(token(i,l,tk))];
+ else
+ dta{i} = [dta{i} {token(i,l,tk)}];
+ end
+ end
+ l = fgetl(fid);
+end
+
+for i=1:length(fields) % copy to structure
+ dta_struct = setfield(dta_struct,fields{i},dta{i});
+end
+
+fclose(fid);
+
+return
+
+%----------------------------------------------------------------------
+
+function tk = token(i,str,tki)
+ tk = str(tki{i}(1):tki{i}(2));
+ return;
+
+function tk = token1(i,str,tki)
+ tk = str(tki{i}(1,1):tki{i}(1,2));
+ return;
+
+function tk = token2(i,str,tki)
+ tk = str(tki{i}(2,1):tki{i}(2,2));
+ return
+
+function s = matlabotomize(s)
+ s = strrep(s,'.','_');
+ s = strrep(s,'3','three');
+
--- a/plot_time_lags.pl Fri Aug 05 11:02:51 2016 -0400
+++ b/plot_time_lags.pl Thu Mar 16 11:53:27 2017 -0400
@@ -1,9 +1,9 @@
#======================================================================
# P L O T _ T I M E _ L A G S . P L
# doc: Tue Jul 28 13:21:09 2015
-# dlm: Tue May 24 22:11:30 2016
+# dlm: Tue Mar 7 12:04:21 2017
# (c) 2015 A.M. Thurnherr
-# uE-Info: 59 81 NIL 0 0 72 2 2 4 NIL ofnI
+# uE-Info: 35 0 NIL 0 0 72 2 2 4 NIL ofnI
#======================================================================
# HISTORY:
@@ -12,7 +12,7 @@
# Mar 16, 2016: - adapted to gmt5
# May 18, 2016: - added version
# May 24, 2016: - fixed for partial-depth casts
-
+# mar 7, 2017: - added time lines for -p
require "$ANTS/libGMT.pl";
@@ -30,7 +30,19 @@
my($R) = "-R$xmin/$xmax/$ymin/$ymax";
GMT_begin($pfn,'-JX10',$R,'-P');
- GMT_psxy('-Sc0.1 -Gcoral');
+ GMT_psxy('-W1,grey30'); # time lines
+ for (my($x)=round($xmin,10); $x<=$xmax; $x+=10) {
+ printf(GMT "%f $ymin\n%f $ymax\n>\n",$x,$x);
+ }
+
+ GMT_psxy('-W8,yellow'); # offset runs
+ for (my($i)=0; $i<@bmo_buf; $i++) {
+ printf(GMT ">\n%f %f\n%f %f\n",
+ $fg_buf[$i]/60-0.5,$bmo_buf[$i],
+ $lg_buf[$i]/60+0.5,$bmo_buf[$i]);
+ }
+
+ GMT_psxy('-Sc0.1 -Gcoral'); # individual offsets
for (my($wi)=0; $wi<@elapsed_buf; $wi++) {
last unless ($elapsed_buf[$wi]<$LADCP{ENSEMBLE}[$LADCP_atbottom]->{ELAPSED});
printf(GMT "%f %f\n",$elapsed_buf[$wi]/60,$so_buf[$wi]);
@@ -41,14 +53,7 @@
printf(GMT "%f %f\n",$elapsed_buf[$wi]/60,$so_buf[$wi]);
}
- GMT_psxy('-W1,grey20');
- for (my($i)=0; $i<@bmo_buf; $i++) {
- printf(GMT ">\n%f %f\n%f %f\n",
- $fg_buf[$i]/60-0.5,$bmo_buf[$i],
- $lg_buf[$i]/60+0.5,$bmo_buf[$i]);
- }
-
- GMT_unitcoords(); # LABELS
+ GMT_unitcoords(); # labels
GMT_pstext('-F+f9,Helvetica,orange+jTR -N -Gwhite');
print(GMT "0.99 0.99 V$VERSION\n");
GMT_pstext('-F+f14,Helvetica,blue+jTL -N');
--- a/time_lag.pl Fri Aug 05 11:02:51 2016 -0400
+++ b/time_lag.pl Thu Mar 16 11:53:27 2017 -0400
@@ -1,9 +1,9 @@
#======================================================================
# T I M E _ L A G . P L
# doc: Fri Dec 17 21:59:07 2010
-# dlm: Mon Mar 7 18:31:34 2016
+# dlm: Tue Mar 7 09:03:20 2017
# (c) 2010 A.M. Thurnherr
-# uE-Info: 71 68 NIL 0 0 72 2 2 4 NIL ofnI
+# uE-Info: 73 63 NIL 0 0 72 2 2 4 NIL ofnI
#======================================================================
# HISTORY:
@@ -69,6 +69,8 @@
# Feb 19, 2016: - added support for -l
# - added warning
# Mar 7, 2016: - BUG: editing did not work correctly in all cases
+# Mar 6, 2017: - BUG: assertion in mad_w failed with 2017 P18 DL#206
+# Mar 9, 2017: - tightened timelag editing (good_diff: 2->1)
# DIFFICULT STATIONS:
# NBP0901#131 this requires the search-radius doubling heuristic
@@ -78,23 +80,24 @@
my($TINY) = 1e-6;
-sub mad_w($$$) # mean absolute deviation
+sub mad_w($$$) # mean absolute deviation
{
- my($fe,$le,$so) = @_; # first/last LADCP ens, CTD scan offset
- my($sad) = my($n) = 0;
+ my($fe,$le,$so) = @_; # first/last LADCP ens, CTD scan offset
my($LADCP_mean_w,$CTD_mean_w,$nsamp) = (0,0,0);
- for (my($e)=$fe; $e<=$le; $e++) { # first, calculate mean w in window
+ for (my($e)=$fe; $e<=$le; $e++) { # first, calculate mean w in window
my($s) = int(($LADCP{ENSEMBLE}[$e]->{ELAPSED} + $CTD{TIME_LAG} - $CTD{ELAPSED}[0]) / $CTD{DT} + 0.5);
- next unless ($s>=0 && $s<=$#{$CTD{ELAPSED}});
+
+# THE FOLLOWING LINE CAUSES AN ASSERTION FAILURE WITH 2017 P08 DL#206. I AM NOT SURE WHETHER MY
+# FIX SOLVES THE UNDERLYING PROBLEM OR ONLY THIS SPECIAL CASE.
+# next unless ($s>=0 && $s<=$#{$CTD{ELAPSED}});
+
+ next unless ($s>0 && $s<=$#{$CTD{ELAPSED}});
die("assertion failed\n" .
"\ttest: abs($LADCP{ENSEMBLE}[$e]->{ELAPSED} + $CTD{TIME_LAG} - $CTD{ELAPSED}[$s]) <= $CTD{DT}/2\n" .
"\te = $e, s = $s, ensemble = $LADCP{ENSEMBLE}[$e]->{NUMBER}"
) unless (abs($LADCP{ENSEMBLE}[$e]->{ELAPSED} + $CTD{TIME_LAG} - $CTD{ELAPSED}[$s]) <= $CTD{DT}/2+$TINY);
next unless numberp($LADCP{ENSEMBLE}[$e]->{REFLR_W});
- my($dw) = $LADCP{ENSEMBLE}[$e]->{REFLR_W}-$LADCP_mean_w - ($CTD{W}[$s+$so]-$CTD_mean_w);
- next unless (abs($dw) <= $w_max_lim);
-
$LADCP_mean_w += $LADCP{ENSEMBLE}[$e]->{REFLR_W};
$CTD_mean_w += $CTD{W}[$s+$so];
$nsamp++;
@@ -103,15 +106,18 @@
$LADCP_mean_w /= $nsamp;
$CTD_mean_w /= $nsamp;
- for (my($e)=$fe; $e<=$le; $e++) { # now, calculate mad
+ my($sad) = 0; # now, calculate mad
+ for (my($e)=$fe; $e<=$le; $e++) {
+ next unless numberp($LADCP{ENSEMBLE}[$e]->{REFLR_W});
my($s) = int(($LADCP{ENSEMBLE}[$e]->{ELAPSED} + $CTD{TIME_LAG} - $CTD{ELAPSED}[0]) / $CTD{DT} + 0.5);
+ next unless ($s>=0 && $s<=$#{$CTD{ELAPSED}});
+ next unless numberp($LADCP{ENSEMBLE}[$e]->{REFLR_W});
my($dw) = $LADCP{ENSEMBLE}[$e]->{REFLR_W}-$LADCP_mean_w - ($CTD{W}[$s+$so]-$CTD_mean_w);
- next unless numberp($LADCP{ENSEMBLE}[$e]->{REFLR_W});
+# print(STDERR "dw = $dw ($LADCP{ENSEMBLE}[$e]->{REFLR_W}-$LADCP_mean_w - ($CTD{W}[$s+$so]-$CTD_mean_w)\n");
next unless (abs($dw) <= $w_max_lim);
$sad += abs($dw);
- $n++;
}
- return ($n>0) ? $sad/$n : 9e99; # n == 0, e.g. in bottom gap
+ return $sad/$nsamp;
}
@@ -121,12 +127,16 @@
my($bestso) = 0; # error at first-guess offset
my($bestmad) = mad_w($fe,$le,0);
+# print(STDERR "bestLag($fe,$le,$ww,$soi)\n");
for (my($dso) = 1; $dso <= int($ww/2/$CTD{DT} + 0.5); $dso+=$soi) {
my($mad) = mad_w($fe,$le,-$dso);
+# print(STDERR "-$dso $mad\n");
$bestmad=$mad,$bestso=-$dso if ($mad < $bestmad);
$mad = mad_w($fe,$le,$dso);
+# print(STDERR " $dso $mad\n");
$bestmad=$mad,$bestso=$dso if ($mad < $bestmad);
}
+# print(STDERR "-> $bestso $bestmad\n");
return ($bestso,$bestmad);
}
@@ -279,7 +289,7 @@
#----------------------------------------------------
my(@fg,@lg);
- my($min_runlength) = 7; my($scan_runlength) = 7; my($min_good) = 4; my($good_diff) = 2;
+ my($min_runlength) = 7; my($scan_runlength) = 7; my($min_good) = 4; my($good_diff) = 1;
unless ($failed || $scan_increment>1) {
my($state) = 0;
for (my($i)=0; 1; $i++) {
--- a/version.pl Fri Aug 05 11:02:51 2016 -0400
+++ b/version.pl Thu Mar 16 11:53:27 2017 -0400
@@ -1,9 +1,9 @@
#======================================================================
# V E R S I O N . P L
# doc: Tue Oct 13 10:40:57 2015
-# dlm: Thu Jun 2 12:20:26 2016
+# dlm: Sun Mar 12 12:11:05 2017
# (c) 2015 A.M. Thurnherr
-# uE-Info: 27 110 NIL 0 0 72 0 2 4 NIL ofnI
+# uE-Info: 26 15 NIL 0 0 72 0 2 4 NIL ofnI
#======================================================================
# HISTORY:
@@ -17,12 +17,14 @@
# Apr 16, 2016: - V1.2beta8
# May 12, 2016: - V1.2
# May 19, 2016: - updated ADCP tools to V1.6
+# Aug 5, 2017: - updated ANTS lib to V6.7
+# Mar 12, 2017: - updated ANTS lit to V6.8
#$VERSION = '1.1'; # Jan 4, 2016
#$VERSION = '1.2'; # May 12, 2016
-$VERSION = '1.3beta2';
+$VERSION = '1.3';
-$antsMinLibVersion = 6.6;
-$ADCP_tools_minVersion = 1.7; # May 2016 (RDI_Coords with bin interpolation & better pitch/roll rotation)
+$antsMinLibVersion = 6.8; # Feb 2017: - .lsfit.poly, .nminterp.linear added
+$ADCP_tools_minVersion = 1.9; # Feb 2017: - round() namespace clash resolved