Dec_17_2010/mergeCTD+LADCP.LADCPtimeSeries
author A.M. Thurnherr <athurnherr@yahoo.com>
Thu, 16 Apr 2015 10:31:09 +0000
changeset 26 b89d4b01fcc5
parent 0 3365828b1004
permissions -rw-r--r--
after implementing new schemes for plotting and active files

#======================================================================
#                    M E R G E C T D + L A D C P . L A D C P T I M E S E R I E S 
#                    doc: Sun May 23 16:40:53 2010
#                    dlm: Sat Oct 30 22:23:33 2010
#                    (c) 2010 A.M. Thurnherr
#                    uE-Info: 11 60 NIL 0 0 72 0 2 4 NIL ofnI
#======================================================================

# HISTORY:
#	May 23, 2010: - created from [perl-tools/RDI_Utils.pl]
#	Oct 20, 2010: - disabled max_gap profile-restarting code

# NOTES:
#	- resulting DEPTH field based on integrated w without any sound speed correction
#	- single-ping ensembles assumed, i.e. no percent-good tests applied
#	- specified bin numbers are 1-relative

sub ref_lr_w($$$$)										# calc ref-layer vert vels
{
	my($dta,$ens,$rl_b0,$rl_b1) = @_;
	my($i,@n,@bn,@v,@vel,@bv,@w);

	for ($i=$rl_b0-1; $i<=$rl_b1-1; $i++) {
		if (defined($dta->{ENSEMBLE}[$ens]->{W}[$i])) {							# valid w
			$vel[2] += $dta->{ENSEMBLE}[$ens]->{W}[$i]; $n[2]++;
			$vel[3] += $dta->{ENSEMBLE}[$ens]->{ERRVEL}[$i], $n[3]++ if defined($dta->{ENSEMBLE}[$ens]->{ERRVEL}[$i]);
			push(@w,$dta->{ENSEMBLE}[$ens]->{W}[$i]); 							# for stderr test
		}
	}

	my($w) = $n[2] ? $vel[2]/$n[2] : undef;				# w uncertainty
	my($sumsq) = 0;
	for ($i=0; $i<=$#w; $i++) {
		$sumsq += ($w-$w[$i])**2;
	}
	my($stderr) = $n[2]>=2 ? sqrt($sumsq)/($n[2]-1) : undef;
#	The following stderr test introduces a huge gap near the bottom of
#	the profiles. Without it, they seem more reasonable.
#	next if ($stderr > 0.05);

	if (defined($w)) {									# valid w
		$dta->{ENSEMBLE}[$ens]->{REFLR_W} = $w;
		$dta->{ENSEMBLE}[$ens]->{REFLR_W_ERR} = $stderr;
	}
}

#======================================================================
# ($firstgood,$lastgood,$atbottom,$w_gap_time) =
#	calcLADCPts($dta,$lr_b0,$lr_b1,$min_corr,$max_e,$max_gap);
#======================================================================

sub calcLADCPts($$$)
{
	my($dta,$rl_b0,$rl_b1,$max_gap) = @_;
	my($firstgood,$lastgood,$atbottom,$w_gap_time,$maxz);

	for (my($z)=0,my($e)=0; $e<=$#{$dta->{ENSEMBLE}}; $e++) {

		ref_lr_w($dta,$e,$rl_b0,$rl_b1);
	
		if (defined($firstgood)) {
			$dta->{ENSEMBLE}[$e]->{ELAPSED_TIME} =			# time since start
				$dta->{ENSEMBLE}[$e]->{UNIX_TIME} -
				$dta->{ENSEMBLE}[$firstgood]->{UNIX_TIME};
		} else {
			if (defined($dta->{ENSEMBLE}[$e]->{REFLR_W})) {		# start of prof.
				$firstgood = $lastgood = $e;		    
				$dta->{ENSEMBLE}[$e]->{ELAPSED_TIME} = 0;
			}
			next;
		}
	
		#--------------------------------------------------
		# within profile: both $firstgood and $lastgood set
		#--------------------------------------------------
	
		if (!defined($dta->{ENSEMBLE}[$e]->{REFLR_W})) {			# gap
			$w_gap_time += $dta->{ENSEMBLE}[$e]->{UNIX_TIME} -
						   $dta->{ENSEMBLE}[$e-1]->{UNIX_TIME};
			next;
		}
	
		my($dt) = $dta->{ENSEMBLE}[$e]->{UNIX_TIME} -		# time step since
				  $dta->{ENSEMBLE}[$lastgood]->{UNIX_TIME}; # ... last good ens
	
#		if ($dt > $max_gap) {
#			printf(STDERR "\n\nWARNING: long gap (%ds). Profile restarted at ensemble #$dta->{ENSEMBLE}[$e]->{NUMBER}\n",$dt);
#			$firstgood = $lastgood = $e;
#			undef($atbottom); undef($maxz);
#			$z = 0;
#			$dta->{ENSEMBLE}[$e]->{ELAPSED_TIME} = 0;
#			$w_gap_time = 0;
#			next;
#		}
	
		#-----------------------------------
		# The current ensemble has a valid w
		#-----------------------------------
	
		$z += $dta->{ENSEMBLE}[$lastgood]->{REFLR_W} * $dt;			# integrate
	
		$atbottom = $e, $maxz = $z if ($z > $maxz); 
		$lastgood = $e;
	}
	
	return ($firstgood,$lastgood,$atbottom,$w_gap_time);
}

1;