From 8cf7726af8ef16cd7be68dfe10013fe8c732f376 Mon Sep 17 00:00:00 2001
From: Jens Ahrensfeld <jens@jayfield.org>
Date: Thu, 23 May 2019 16:41:45 +0000
Subject: [PATCH] - added finding peak width

git-svn-id: http://moon:8086/svn/software/trunk/projects/GnuRadio@445 b431acfa-c32f-4a4a-93f1-934dc6c82436
---
 matlab/peak_detect.m | 96 ++++++++++++++++++++++++++++++--------------
 1 file changed, 65 insertions(+), 31 deletions(-)

diff --git a/matlab/peak_detect.m b/matlab/peak_detect.m
index 0dd9006..f2f7c6a 100644
--- a/matlab/peak_detect.m
+++ b/matlab/peak_detect.m
@@ -38,10 +38,12 @@ pa = 0.3;
 d_max = 30;
 
 am = 0.04;
-as = 0.9;
-Ad = 15;
+as = 0.2;
+Ad = 20;
 Cm = 0;
-pbh = 5;
+
+pb_min = 5
+pb_max = 37
 
 s = zeros(1, Nh);
 a = zeros(1, Nh);
@@ -83,11 +85,12 @@ for k=1:numFrames
 	% Init 
 	if k == 1
 		Xm = Xm + mean(X);
+		Xs = Xs + mean(X);
 	endif
 	
 	% Peak detection
 	Xs = (1-as)*Xs + as*X;
-	gate_peak = Xs > (Ad+Xm);
+	gate_peak = (Xs-Xm) >= Ad;
 	peak_count = gate_peak.*peak_count + gate_peak;
 
 	% Output
@@ -133,35 +136,73 @@ for k=1:numFrames
 	peak_boxes_rel = zeros(1, Nh);
 	peak_boxes_abs = zeros(1, Nh);
 	for x = peak_cand_sorted,
-			pbox_start = max(1, x-pbh);
-			pbox_end = min(Nh, x+pbh);
-			
-			% Find intersection
-			box_height_rel = Xs(x) - Xm(x);
-			box_height_abs = X(x);
-			does_intersect = 0;
-			for n=pbox_start:pbox_end
-				if peak_boxes_rel(n) > box_height_rel
-					does_intersect = 1;
-					break
+		
+		% Find peak widths
+		left = x;
+		right = x;
+		nom = X(x) - 10;
+		left_found = 0;
+		right_found = 0;
+		while(~(left_found & right_found))
+			% Determine peak width left from center (x)
+			if X(left) > nom
+				if left > 1
+					left = left - 1;
+				else
+					left_found = 1;
 				end
-			endfor
-			if does_intersect == 0
-				peak_boxes_rel(pbox_start:pbox_end) = box_height_rel * ones(1, pbox_end - pbox_start + 1);
-				peak_boxes_abs(pbox_start:pbox_end) = box_height_abs * ones(1, pbox_end - pbox_start + 1);
-				peak_final = [peak_final x];
+			else
+				left_found = 1;
+			endif
+			% Determine peak width right from center (x)
+			if X(right) > nom
+				if right < Nh
+					right = right + 1;
+				else
+					right_found = 1;
+				end
+			else
+				right_found = 1;
+			endif
+			pbh_left = x - left;
+			pbh_right = right - x;
+			
+			% Take larger peak width
+			pbh = max(pbh_left, pbh_right);
+			
+			% Ensure peak box is at least pb_min
+			pbh = max(pbh, pb_min);
+		end
+		
+		% Use peak width for peak box
+		pbox_start = max(1, x-pbh);
+		pbox_end = min(Nh, x+pbh);
+		
+		% Find intersection
+		box_height_rel = X(x) - Xm(x);
+		box_height_abs = X(x);
+		does_intersect = 0;
+		for n=pbox_start:pbox_end
+			if peak_boxes_rel(n) > box_height_rel
+				does_intersect = 1;
+				break
 			end
+		endfor
+		if does_intersect == 0
+			peak_boxes_rel(pbox_start:pbox_end) = box_height_rel * ones(1, pbox_end - pbox_start + 1);
+			peak_boxes_abs(pbox_start:pbox_end) = box_height_abs * ones(1, pbox_end - pbox_start + 1);
+			peak_final = [peak_final x];
+		end
 	endfor
 
 	% Update Xm at no-peak bins
-	gate_upd = peak_boxes_rel > 0;
 	Xm = (1-am)*Xm + am*X;
 
 	% Fill peak gaps with guessed values
 	last = Xm(1);
 	for n=2:Nh,
-		if gate_upd(n) == 0
-			last = 0.5*last + 0.5*Xm(n);
+		if (Xm(n) < (last+2)) or (peak_boxes_rel(n) == 0);
+			last = (1-am)*last + am*Xm(n);
 		else
 			Xm(n) = last;	
 		end
@@ -175,15 +216,8 @@ for k=1:numFrames
 	subplot(2, 1, 1)
 	plot(1:Nh, X, peak_final, X(peak_final), 'ro', 1:Nh, peak_box_values, 'm-'); grid; axis([0, Nh, -60, 40]);
 	subplot(2, 1, 2)
-	plot(1:Nh, Xs, 1:Nh, Xm); grid; axis([0, Nh, -60, 40]);
+	plot(1:Nh, Xs-Xm, 1:Nh, Xm); grid; axis([0, Nh, -60, 40]);
 	pause(1/30);
 end
 
-function s = matchedFilter(N, B)
-n2 = (N-1)/2
-n = (0:N-1)-n2
-s = sinc(n/B)
-
-endfunction
-	
 endfunction