Merge branch 'JC159_post_cruise' into 'master'

Jc159 post cruise

See merge request !1

source/mextras/ctd_apply_celltm.m

+function cond = ctd_apply_celltm(time,temp,cond);
+
+% dy040 bak 21 dec 2015
+% code up celltm algorithm
+
+alpha = 0.03; % sbe default values are alpha = 0.03; 1/beta = 7;
+beta = 1/7;
+
+num = numel(time);
+
+ctm = zeros(size(cond));
+
+kfirst = min(find(isfinite(time+temp)));
+
+if isempty(kfirst)
+    msg = ['No finite time+temp data found for celltm calculation in ctd_apply_celltm'];
+    fprintf(2,'%s\n',msg);
+    cond = cond+nan;
+    return
+end
+
+timelast = time(kfirst);
+templast = temp(kfirst);
+ctmlast = ctm(kfirst);
+
+for kl = kfirst+1:num
+    if ~isfinite(time(kl)+temp(kl))
+        ctm(kl) = ctm(kl-1); % keep ctm fixed if not possible to update it due to missing temperature
+        ctmlast = ctm(kl);
+        continue
+    end
+    dtime = time(kl)-timelast;
+    dtemp = temp(kl)-templast;
+    
+    a = 2 * alpha/(dtime*beta + 2);
+    b = 1 - (2 * a/alpha);
+    dcdt = 0.1*(1+0.006*(temp(kl)-20));
+    ctm(kl) = -1 * b * ctmlast + a*dcdt*dtemp;
+    timelast = time(kl);
+    templast = temp(kl);
+    ctmlast = ctm(kl);
+end
+
+
+ctm = 10 * ctm;
+
+cond = cond+ctm;
+
+return

source/mextras/mcarter.m

@@ -38,7 +38,7 @@ function dep = mcarter(lat,lon,uncdep)
 
 m_common
 MEXEC_A.Mprog = 'mcarter';
-if ~MEXEC_G.quiet; m_proghd; end
+m_proghd
 
 % fn = [MCARTER_DIRECTORY '/carter'];
 % d = load(fn); % load mat file.
@@ -119,4 +119,4 @@ dep.carter_area = carea;
 dep.default_names = {'cordep' 'carterarea'};
 dep.default_units = {'metres' 'carterarea'};
 
-return
\ No newline at end of file
+return

source/mextras/mcontr_plotlines.m

@@ -42,7 +42,6 @@ while cdim > 0
     cdim = size(posall,2);
 end
 hold on
-% keyboard
 % plot each contour line
 % first switch off lines in base contour plot
 set(hplot,'linestyle','none');
@@ -150,7 +149,6 @@ for kcont = 1:numcontour
     end
 
 end
-% keyboard
 % % % % % % % clabelall = get(hclabel);
 % % % % % % % clall = clabelall;
 % % % % % % % nstring = length(clabelall);

source/mextras/mcoxyhyst.m

-function [oxygen_out]=mcoxyhyst(oxygen_sbe,time,press,H1,H2,H3)
-% function [oxygen_out]=mcoxyhyst(oxygen_sbe,time,press,H1,H2,H3)
+function [oxygen_out C D]=mcoxyhyst(oxygen_sbe,time,press,H1in,H2in,H3in,sensor)
+% [oxygen_out]=mcoxyhyst(oxygen_sbe,time,press,H1,H2,H3)
+% [oxygen_out]=mcoxyhyst(oxygen_sbe,time,press,H1,H2,H3,sensor)
+% [oxygen_out C D]=mcoxyhyst(oxygen_sbe,time,press,H1,H2,H3)
 %
 % gdm on di346
 % function to apply an adjustment for hysteresis in the oxygen sensor
@@ -11,26 +13,63 @@ function [oxygen_out]=mcoxyhyst(oxygen_sbe,time,press,H1,H2,H3)
 % H1     -0.033     [-0.02 to -0.05]
 % H2      5000      
 % H3      1450      [1200  to  2000]
+%
+% overhaul by bak on dy040, following GDM on dy039. Need to allow pars to
+% vary with depth. GDM allowed choice of two. This script allows pars to be
+% defined as arbitrary functions of depth.
+% H1 H2 H3 are now arrays. Default is single, passed-in
+% value, uniform with pressure. Option is now to vary with depth. Should be
+% fully backwards compatible.
+
+m_common % on and after dy040, identify cruise from m_setup
+scriptname = 'mcoxyhyst';
+mcruise = MEXEC_G.MSCRIPT_CRUISE_STRING;
+
+if nargin == 6
+   % sensor not defined in call, assume sensor = 1, which is backwards compatible
+   sensor = 1;
+end
 
 oxygen_out=oxygen_sbe;
 
+D = nan(size(oxygen_sbe));
+C = D;
+
 kfirst = min(find(isfinite(oxygen_sbe)));
 klastgood = kfirst; % keep track of most recent good cycle
 
+oxygen_out(1) = oxygen_sbe(1)/1;
+
+%default for cruises before dy040, and when nothing fancy is required
+zz = zeros(size(D));
+H1 = H1in + zz;
+H2 = H2in + zz;
+H3 = H3in + zz;
+
+% now do fancy things on cruises that need it
+%i.e. depth-varying parameters
+get_cropt
+
+D(kfirst)=1+H1(kfirst)*(exp(press(kfirst)/H2(kfirst))-1);
+
 for k=kfirst+1:length(time)
     % bak: 23 jan 2010 need to be able to step over nans
     %     oxygen_out(k)=((oxygen_sbe(k)+(oxygen_out(k-1)*C*D))-(oxygen_sbe(k-1)*C))/D;
     % bak: 29 feb 2012 there are some nans in press after cleaning up raw
     %     file on jc069_064. raw file has spikes due to noisy telemetry
     %     through slip rings
-    %     therefore oxygen is nan if eitehr oxygen_sbe or press is nan
-    if isnan(oxygen_sbe(k)+press(k))
+    %     therefore oxygen is nan if either oxygen_sbe or press is nan
+    
+    if ~isfinite(oxygen_sbe(k)+press(k))
         oxygen_out(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_out(k)=((oxygen_sbe(k)+(oxygen_out(klastgood)*C*D))-(oxygen_sbe(klastgood)*C))/D;
+        D(k)=1+H1(k)*(exp(press(k)/H2(k))-1);
+        C(k)=exp(-1*(time(k)-time(klastgood))/H3(k));
+
+        oxygen_out(k)=((oxygen_sbe(k)+(oxygen_out(klastgood)*C(k)*D(k)))-(oxygen_sbe(klastgood)*C(k)))/D(k);
         klastgood = k;
+        
     end
-end;
+end
+return

source/mextras/mfallrate.m

@@ -21,7 +21,7 @@ m_margslocal
 m_varargs
 
 MEXEC_A.Mprog = 'mcalc';
-if ~MEXEC_G.quiet; m_proghd; end
+m_proghd
 
 
 fprintf(MEXEC_A.Mfidterm,'%s','Enter name of input disc file  ')

source/mextras/mposspd.m

@@ -33,7 +33,7 @@ m_margslocal
 m_varargs
 
 MEXEC_A.Mprog = 'mposspd';
-if ~MEXEC_G.quiet; m_proghd; end
+m_proghd
 
 fprintf(MEXEC_A.Mfidterm,'%s','Input file name     ')
 fn_in = m_getfilename;
@@ -244,4 +244,4 @@ hist.filename = ncfile_ot.name;
 MEXEC_A.Mhistory_ot{1} = hist;
 m_write_history;
 
-return
\ No newline at end of file
+return

source/mextras/mposspd_polar.m

@@ -34,7 +34,7 @@ m_margslocal
 m_varargs
 
 MEXEC_A.Mprog = 'mposspd';
-if ~MEXEC_G.quiet; m_proghd; end
+m_proghd
 
 fprintf(MEXEC_A.Mfidterm,'%s','Input file name     ')
 fn_in = m_getfilename;
@@ -243,4 +243,4 @@ hist.filename = ncfile_ot.name;
 MEXEC_A.Mhistory_ot{1} = hist;
 m_write_history;
 
-return
\ No newline at end of file
+return

source/mscs/msload.m

@@ -58,7 +58,7 @@ nv = length(vars);
 if ~exist('varlist','var'); varlist = '/'; end
 th.fldnam = vars;
 th.noflds = nv; % create a structure equivalent to the mstar headers to parse for var names
-% keyboard
+
 varnums = m_getvlist(varlist,th);
 % time always seems to be last in the scs list; put it first if it is
 % in the load list.
@@ -75,7 +75,6 @@ end
 dc1 = nan+ones(nf,1); dc2 = dc1; totdc = 0;
 m = 'Counting data cycles';
 if ~MEXEC_G.quiet; fprintf(MEXEC_A.Mfidterm,'%s\n',m); end
-end
 
 for kf = 1:nf
     fn = fnames{kf};

source/mscs/msposinfo.m

@@ -25,7 +25,7 @@ function [lat lon] = msposinfo(dn1,navstream)
 
 m_common
 
-if ~exist('navstream','var'); navstream = MEXEC_G.uway_default_navstream; end  % use default nav stream name
+if ~exist('navstream','var'); navstream = MEXEC_G.default_navstream; end  % use default nav stream name
 
 instream = navstream; % mexec stream short name
 tstream = msresolve_stream(instream);

source/mscs/mtranslate_varnames.m

@@ -34,8 +34,8 @@ tunder(strfind(tunder,'-')) = '_';
 tunder(strfind(tunder,'.')) = '_';
 
 root_template = mgetdir('M_TEMPLATES');
-fntemplatein = [root_template '/' MEXEC_G.Mshipdatasystem '_' MEXEC_G.MSCRIPT_CRUISE_STRING '_renamelist_' tunder '.csv'];
-fntemplateot = [root_template '/' MEXEC_G.Mshipdatasystem '_' MEXEC_G.MSCRIPT_CRUISE_STRING '_renamelist_' tunder '_out.csv'];
+fntemplatein = [root_template '/' MEXEC_G.Mshipdatasystem '_renamelist_' tunder '.csv'];
+fntemplateot = [root_template '/' MEXEC_G.Mshipdatasystem '_renamelist_' tunder '_out.csv'];
 
 % clean up the translation table if needed. This code was lifted from
 % another script but should not be needed for this application.
@@ -48,6 +48,7 @@ for kline = 1:length(cellall)
     snamesin{kline} = m_remove_outside_spaces(cellrow{1});
     snamesot{kline} = m_remove_outside_spaces(cellrow{2});
     sunits{kline} = m_remove_outside_spaces(cellrow{3});
+    if length(sunits{kline})==0; sunits{kline} = ' '; end
 end
 snamesin = snamesin(:);
 snamesot = snamesot(:);
@@ -77,6 +78,7 @@ for k = 1:numvar
     end
 end
 snames_units = snames_units(:);
+
 %--------------------------------
 % 2009-01-26 07:48:13
 % mheadr

source/mscs/scs_to_mstar2

@@ -9,7 +9,7 @@ m_margslocal
 m_varargs
 
 Mprog = 'techsas_to_mstar2';
-if ~exist('vmode','var') | ~strcmp(vmode,'q'); m_proghd; end
+m_proghd
 
 
 instream = m_getinput('Type techsas stream name or mexec short name eg SBE-SBE45.TSG or posmvpos : ','s');
@@ -126,4 +126,4 @@ histin.version = [];
 histin.mstar_site = [];
 Mhistory_in{1} = histin;
 m_write_history;
-return
\ No newline at end of file
+return

source/mscs/scs_to_mstar2.m

 function scs_to_mstar2
-% function techsas_to_mstar(tstream,dn1,dn2)
-% function ncfile = techsas_to_mstar(techsas_in,ncfile,dataname)
+% function scs_to_mstar2(tstream,dn1,dn2)
+% function ncfile = scs_to_mstar2(scs_in,ncfile,dataname)
 
-% load techsas file into mstar file
+% load scs file into mstar file
 %
 % 8 Sep 2009: SCS version of original techsas script, for JR195
 % The searched directory is MEXEC_G.uway_root, which for example can be
@@ -10,13 +10,14 @@ function scs_to_mstar2
 % The var names and units are taken from ascii file
 % seatex-gga.TPL
 % for example.
+%
 
 m_common
 m_margslocal
 m_varargs
 
 MEXEC_A.Mprog = 'scs_to_mstar2';
-m_proghd;
+m_proghd
 
 
 instream = m_getinput('Type scs stream name or mexec short name eg gyro_s or gyro : ','s');
@@ -71,6 +72,7 @@ if isempty(tdata.time) % no data cycles found
 end
 
 ncfile.name = mstar_fn;
+
 ncfile = m_openot(ncfile); %check it is not an open mstar file
 % techsas_in.name = techsas_fn;
 
@@ -127,7 +129,6 @@ m_add_comment(ncfile,tstream);
 m_add_comment(ncfile,['at ' nowstring]);
 m_add_comment(ncfile,['Time converted from matlab day number to seconds after mstar time origin']);
     
-
 m_finis(ncfile);
 
 h = m_read_header(ncfile);
@@ -145,4 +146,5 @@ histin.version = [];
 histin.mstar_site = [];
 MEXEC_A.Mhistory_in{1} = histin;
 m_write_history;
+
 return

source/msubs/m_add_variable_name.m

@@ -24,7 +24,6 @@ ncfile.metadata = metadata;
 dimnames = m_unpack_dimnames(ncfile);
 varnames = m_unpack_varnames(ncfile);
 
-
 %test to see if a variable already exists with correct dimensions; if not, add it.
 
 kmatch = strmatch(string_var_name,varnames,'exact');
@@ -76,8 +75,7 @@ if ~isempty(kmatch) % It exists; check dimensions match;
     return
 end
 
- % Variable doesn't exist so we will add it; First, check each requested dimension exists in file
-
+% Variable doesn't exist so we will add it; First, check each requested dimension exists in file
 for k = 1:length(dims_array)
 
     kmatch = strmatch(dims_array{k},dimnames,'exact');
@@ -106,6 +104,7 @@ end
 v.Name = string_var_name;
 v.Dimension = dims_array;
 v.Nctype = datatype;
+%keyboard
 nc_addvar(ncfile.name,v);
 
 % % metadata = nc_info(ncfile.name); %refresh metadata

source/msubs/m_write_variable.m

@@ -41,7 +41,6 @@ for k = 1:length(krmatch)
     collength(k) = metadata.Dimension(kcmatch(k)).Length;
 end
 
-
 if length(krmatch) ~= length(kcmatch)
     error('m_write_variable weird mismatch of dimension names - investigate further')
 end

source/mtechsas/mtchoosefiles.m

@@ -35,7 +35,18 @@ for kn = 1:numf
     %if t1 > dn2; continue; end % t1 not in range
     %if t2 < dn1; continue; end % t2 not in range
     % if we get here file is useful
-    if n>0 & t1<dn2 & t2>dn1
+%     if n>0 & t1<dn2 & t2>dn1 % Bug fix by bak jc159 17 april 2018.
+%     Amazing that this bug has lasted so long. Previously the test in the
+%     next statement did not allow for equality. This crashed mtlast on sim
+%     or em120, when t2 == dn1 == dn2 == the last cycle in the file.
+%     In most cases, eg nav files, with files updating at 1Hz, then dn1=dn2 is
+%     determined as the last time in the file in mtlast, and by the time 
+%     the code arrives here in mtchoosefiles, t2 is greater than dn1. 
+%     For files updating at intervals of a few seconds, the old code would
+%     work sometimes and fail sometimes depending on whetehr the file updated
+%     between mtlast and the sub call to mtchoosefiles. Allowing t2 >= dn1 will select a
+%     file to keep when dn1 is the last time in the file.
+    if n>0 & t1<=dn2 & t2>=dn1
         keep(kn) = 1;
     end
 end