From bb389ad7df57e4ae11b314d77f5230cb33de20d9 Mon Sep 17 00:00:00 2001
From: PStar User Account <pstar@192.168.122.1>
Date: Tue, 8 Oct 2019 10:36:01 +0000
Subject: [PATCH] mcoxyhyst_reverse that deals with mising data

---
 source/mextras/mcoxyhyst_reverse.m | 54 +++++++++++++++++++-----------
 1 file changed, 34 insertions(+), 20 deletions(-)
 mode change 100755 => 100644 source/mextras/mcoxyhyst_reverse.m

diff --git a/source/mextras/mcoxyhyst_reverse.m b/source/mextras/mcoxyhyst_reverse.m
old mode 100755
new mode 100644
index 8eaafb9..a9d97e5
--- a/source/mextras/mcoxyhyst_reverse.m
+++ b/source/mextras/mcoxyhyst_reverse.m
@@ -1,21 +1,35 @@
-function [oxygen_rev]=mcoxyhyst_reverse(oxygen_sbe,time,press,H1,H2,H3)
-% function [oxygen_rev]=mcoxyhyst_reverse(oxygen_sbe,time,press,H1,H2,H3)
-%
-% gdm on di346
-% function to reverse the adjustment for hysteresis in the oxygen sensor
-% applied by SeaBird 
-% 
-% the default values and ranges or the constants are:
-% H1     -0.033     [-0.02 to -0.05]
-% H2      5000      
-% H3      1450      [1200  to  2000]
-
-oxygen_rev=oxygen_sbe;
-
-for k=2:length(time)
-    D=1+H1*(exp(press(k)/H2)-1);
-    C=exp(-1*(time(k)-time(k-1))/H3);
-
-    oxygen_rev(k)=D*(oxygen_sbe(k)-C*oxygen_sbe(k-1))+C*oxygen_rev(k-1);
-    
+function [oxygen_rev]=mcoxyhyst_reverse(oxygen_sbe,time,press,H1,H2,H3)
+% function [oxygen_rev]=mcoxyhyst_reverse(oxygen_sbe,time,press,H1,H2,H3)
+%
+% gdm on di346
+% function to reverse the adjustment for hysteresis in the oxygen sensor
+% applied by SeaBird 
+% 
+% the default values and ranges or the constants are:
+% H1     -0.033     [-0.02 to -0.05]
+% H2      5000      
+% H3      1450      [1200  to  2000]
+
+% DAS chnaged to deal with absent data but note this could give slightly different results from original data as 
+% the sbe software may deal with absent data differently
+
+oxygen_rev=oxygen_sbe;
+
+kfirst = min(find(isfinite(oxygen_sbe)));
+klastgood = kfirst; % keep track of most recent good cycle
+
+for k=2:length(time)
+    D=1+H1*(exp(press(k)/H2)-1);
+    C=exp(-1*(time(k)-time(k-1))/H3);
+
+    oxygen_rev(k)=D*(oxygen_sbe(k)-C*oxygen_sbe(k-1))+C*oxygen_rev(k-1);
+    if isnan(oxygen_sbe(k)+press(k))
+        oxygen_rev(k) = nan; %already the case because of initialisation of oxygen_out
+    else
+        if press(k) < 0; press(k) = 0; end
+        D=1+H1*(exp(press(k)/H2)-1);
+        C=exp(-1*(time(k)-time(klastgood))/H3);
+		oxygen_rev(k)=D*(oxygen_sbe(k)-C*oxygen_sbe(klastgood))+C*oxygen_rev(klastgood);
+        klastgood = k;
+    end    
 end;
\ No newline at end of file
-- 
GitLab