# HG changeset patch
# User Chris
# Date 1482953017 0
# Wed Dec 28 19:23:37 2016 +0000
# Node ID 86e41f9d73c83cbb244b04b40b51efebbfc8a1ac
# Parent ad7e5d12d69ad6bd4e72114a599a95f026584a93
Added statelevels function code.
diff --git a/inst/statelevels.m b/inst/statelevels.m
new file mode 100644
--- /dev/null
+++ b/inst/statelevels.m
@@ -0,0 +1,80 @@
+## Copyright (C) 2016 Chris Adams
+##
+## This program is free software; you can redistribute it and/or modify it
+## under the terms of the GNU General Public License as published by
+## the Free Software Foundation; either version 3 of the License, or
+## (at your option) any later version.
+##
+## This program is distributed in the hope that it will be useful, but
+## WITHOUT ANY WARRANTY; without even the implied warranty of
+## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
+## GNU General Public License for more details.
+##
+## You should have received a copy of the GNU General Public License
+## along with this program. If not, see .
+
+## -*- texinfo -*-
+## @deftypefn {} {@var{retval} =} statelevels (@var{input1}, @var{input2})
+##
+## @seealso{}
+## @end deftypefn
+
+#Uses histogram method.
+#1 Take a histogram
+#2 find ilow and ihigh, min and max non zero hist values
+#3 make two subhist
+#4 take mean of each for lower and upper
+
+function S = statelevels (X)
+ #TODO
+ #Accept a dynamic number of histogram bins
+
+ sz = size(X);
+ if(sz(1) !=1 | !(sz(2) > 1))
+ error('dutycycle expects a row vector');
+ end
+
+ #Default params
+ #Number of bins to use for estimation
+ nBins = 100;
+
+ #Get a histogram
+ [nn,xx] = hist(X,nBins);
+
+ #find ilow,ihigh, hist indexes with nonzero counts
+ for ilow=1:nBins
+ if(nn(ilow) > 0)
+ break;
+ end
+ end
+
+ for ihigh=nBins:-1:1
+ if(nn(ihigh) > 0)
+ break;
+ end
+ end
+
+
+ #Get the two subhistograms
+ sub1_nn = nn(ilow:ceil(0.5*(ihigh-ilow)));
+ sub1_xx = xx(ilow:ceil(0.5*(ihigh-ilow)));
+ sub2_nn = nn(ilow+ceil((0.5*(ihigh-ilow))):ihigh);
+ sub2_xx = xx(ilow+ceil((0.5*(ihigh-ilow))):ihigh);
+
+ #Get the state levels based off the mean
+ #First calculate the bin centres
+ #It's the difference bewteen bin edges divide 2
+ #Tack a zero on front so we can add to bin edges
+ sub1_xxCentres = [0 diff(sub1_xx)/2] + sub1_xx;
+ #Now times centres by count and sum and divide by number of observations
+ S1 = sum(sub1_xxCentres .* sub1_nn,2)/sum(sub1_nn,2);
+
+ #Do the same for the higher level
+ sub2_xxCentres = [0 diff(sub2_xx)/2] + sub2_xx;
+ #Now times centres by count and sum and divide by number of observations
+ S2 = sum(sub2_xxCentres .* sub2_nn,2)/sum(sub2_nn,2);
+
+ #return
+ S = [S1 S2];
+
+endfunction