time_series.pl
changeset 0 3365828b1004
child 2 a077ea2a9f36
equal deleted inserted replaced
-1:000000000000 0:3365828b1004
       
     1 #======================================================================
       
     2 #                    C A L C _ L A D C P _ T I S . P L 
       
     3 #                    doc: Sun May 23 16:40:53 2010
       
     4 #                    dlm: Tue Dec 21 00:34:15 2010
       
     5 #                    (c) 2010 A.M. Thurnherr
       
     6 #                    uE-Info: 87 22 NIL 0 0 72 2 2 4 NIL ofnI
       
     7 #======================================================================
       
     8 
       
     9 # HISTORY:
       
    10 #	May 23, 2010: - created from [perl-tools/RDI_Utils.pl]
       
    11 #	Oct 20, 2010: - disabled max_gap profile-restarting code
       
    12 #	Dec 17, 2010: - re-added {DEPTH} field to ensembles
       
    13 #	Dec 18, 2010: - max gap re-enabled
       
    14 #	Dec 20, 2010: - cosmetics
       
    15 
       
    16 # NOTES:
       
    17 #	- resulting DEPTH field based on integrated w without any sound speed correction
       
    18 #	- single-ping ensembles assumed, i.e. no percent-good tests applied
       
    19 #	- specified bin numbers are 1-relative
       
    20 
       
    21 sub ref_lr_w($$$$)										# calc ref-layer vert vels
       
    22 {
       
    23 	my($dta,$ens,$rl_b0,$rl_b1) = @_;
       
    24 	my($i,@n,@bn,@v,@vel,@bv,@w);
       
    25 
       
    26 	for ($i=$rl_b0-1; $i<=$rl_b1-1; $i++) {
       
    27 		if (defined($dta->{ENSEMBLE}[$ens]->{W}[$i])) {							# valid w
       
    28 			$vel[2] += $dta->{ENSEMBLE}[$ens]->{W}[$i]; $n[2]++;
       
    29 			$vel[3] += $dta->{ENSEMBLE}[$ens]->{ERRVEL}[$i], $n[3]++ if defined($dta->{ENSEMBLE}[$ens]->{ERRVEL}[$i]);
       
    30 			push(@w,$dta->{ENSEMBLE}[$ens]->{W}[$i]); 							# for stderr test
       
    31 		}
       
    32 	}
       
    33 
       
    34 	my($w) = $n[2] ? $vel[2]/$n[2] : undef;				# w uncertainty
       
    35 	my($sumsq) = 0;
       
    36 	for ($i=0; $i<=$#w; $i++) {
       
    37 		$sumsq += ($w-$w[$i])**2;
       
    38 	}
       
    39 	my($stderr) = $n[2]>=2 ? sqrt($sumsq)/($n[2]-1) : undef;
       
    40 
       
    41 	if (defined($w)) {									# valid w
       
    42 		$dta->{ENSEMBLE}[$ens]->{REFLR_W} = $w;
       
    43 		$dta->{ENSEMBLE}[$ens]->{REFLR_W_ERR} = $stderr;
       
    44 	}
       
    45 }
       
    46 
       
    47 #======================================================================
       
    48 # ($firstgood,$lastgood,$atbottom,$w_gap_time) =
       
    49 #	calcLADCPts($dta,$lr_b0,$lr_b1,$min_corr,$max_e,$max_gap);
       
    50 #======================================================================
       
    51 
       
    52 sub calcLADCPts($$$$)
       
    53 {
       
    54 	my($dta,$rl_b0,$rl_b1,$max_gap) = @_;
       
    55 	my($firstgood,$lastgood,$atbottom,$w_gap_time,$max_depth);
       
    56 
       
    57 	for (my($depth)=0,my($e)=0; $e<=$#{$dta->{ENSEMBLE}}; $e++) {
       
    58 
       
    59 		ref_lr_w($dta,$e,$rl_b0,$rl_b1);
       
    60 	
       
    61 		if (defined($firstgood)) {
       
    62 			$dta->{ENSEMBLE}[$e]->{ELAPSED} =				# time since start
       
    63 				$dta->{ENSEMBLE}[$e]->{UNIX_TIME} -
       
    64 				$dta->{ENSEMBLE}[$firstgood]->{UNIX_TIME};
       
    65 		} else {
       
    66 			if (defined($dta->{ENSEMBLE}[$e]->{REFLR_W})) {		# start of prof.
       
    67 				$firstgood = $lastgood = $e;		    
       
    68 				$dta->{ENSEMBLE}[$e]->{ELAPSED} = 0;
       
    69 			}
       
    70 			next;
       
    71 		}
       
    72 	
       
    73 		unless (defined($dta->{ENSEMBLE}[$e]->{REFLR_W})) {				# gap
       
    74 			$w_gap_time += $dta->{ENSEMBLE}[$e]->{UNIX_TIME} -
       
    75 						   $dta->{ENSEMBLE}[$e-1]->{UNIX_TIME};
       
    76 			next;
       
    77 		}
       
    78 	
       
    79 		my($dt) = $dta->{ENSEMBLE}[$e]->{UNIX_TIME} -		# time step since
       
    80 				  $dta->{ENSEMBLE}[$lastgood]->{UNIX_TIME}; # ... last good ens
       
    81 	
       
    82 		if ($dt > $max_gap) {
       
    83 			if ($max_depth>50 && $depth<0.5*$max_depth) {
       
    84 				warning(1,"long gap (%ds). Profile ended at ensemble #$dta->{ENSEMBLE}[$e]->{NUMBER}\n",$dt);
       
    85 				last;				
       
    86             }
       
    87 			warning(1,"long gap (%ds). Profile restarted at ensemble #$dta->{ENSEMBLE}[$e]->{NUMBER}\n",$dt);
       
    88 			$firstgood = $lastgood = $e;
       
    89 			undef($atbottom); undef($max_depth);
       
    90 			$depth = 0;
       
    91 			$dta->{ENSEMBLE}[$e]->{ELAPSED} = 0;
       
    92 			$w_gap_time = 0;
       
    93 			next;
       
    94 		}
       
    95 	
       
    96 		$depth += $dta->{ENSEMBLE}[$lastgood]->{REFLR_W} * $dt;			# integrate
       
    97 		$dta->{ENSEMBLE}[$e]->{DEPTH} = $depth;
       
    98 	
       
    99 		$atbottom = $e, $max_depth = $depth if ($depth > $max_depth); 
       
   100 		$lastgood = $e;
       
   101 	}
       
   102 	
       
   103 	return ($firstgood,$lastgood,$atbottom,$w_gap_time);
       
   104 }
       
   105 
       
   106 1;