From e8dbaec6c4d503ea352e8ebc9d0d0284a6131bc7 Mon Sep 17 00:00:00 2001 From: toni Date: Mon, 18 Dec 2017 17:02:54 +0100 Subject: [PATCH] close #6 - finished convertig matlab code to java - also added moving average and simple kalman filter - a lot of refactoring for speed --- .../{BpmEstimator.java => Estimator.java} | 4 +- .../conductorswatch/WorkerFragment.java | 16 +- java/src/main/java/Main.java | 30 ++- java/src/main/java/Plot.java | 4 +- .../AccelerometerInterpolator.java | 10 +- .../AccelerometerWindowBuffer.java | 39 ++-- .../java/bpmEstimation/AutoCorrelation.java | 4 +- .../main/java/bpmEstimation/BpmEstimator.java | 189 ++++++++++++++++++ java/src/main/java/bpmEstimation/Peaks.java | 29 ++- .../src/main/java/utilities/MovingFilter.java | 33 +++ .../src/main/java/utilities/SimpleKalman.java | 33 +++ java/src/main/java/utilities/Utils.java | 106 +++------- 12 files changed, 347 insertions(+), 150 deletions(-) rename android/ConductorsWatch/app/src/main/java/de/tonifetzer/conductorswatch/{BpmEstimator.java => Estimator.java} (97%) create mode 100644 java/src/main/java/bpmEstimation/BpmEstimator.java create mode 100644 java/src/main/java/utilities/MovingFilter.java create mode 100644 java/src/main/java/utilities/SimpleKalman.java diff --git a/android/ConductorsWatch/app/src/main/java/de/tonifetzer/conductorswatch/BpmEstimator.java b/android/ConductorsWatch/app/src/main/java/de/tonifetzer/conductorswatch/Estimator.java similarity index 97% rename from android/ConductorsWatch/app/src/main/java/de/tonifetzer/conductorswatch/BpmEstimator.java rename to android/ConductorsWatch/app/src/main/java/de/tonifetzer/conductorswatch/Estimator.java index 7f2c2ad..2c738c1 100644 --- a/android/ConductorsWatch/app/src/main/java/de/tonifetzer/conductorswatch/BpmEstimator.java +++ b/android/ConductorsWatch/app/src/main/java/de/tonifetzer/conductorswatch/Estimator.java @@ -16,7 +16,7 @@ import de.tonifetzer.conductorswatch.utilities.Utils; */ //TODO: einfügen der logik autoCorr + FindPeaks -public class BpmEstimator implements SensorEventListener { +public class Estimator implements SensorEventListener { private SensorManager mSensorManager; private Sensor mAccelerometer; @@ -51,7 +51,7 @@ public class BpmEstimator implements SensorEventListener { }; - public BpmEstimator(Context mContext){ + public Estimator(Context mContext){ this.mContext = mContext; } diff --git a/android/ConductorsWatch/app/src/main/java/de/tonifetzer/conductorswatch/WorkerFragment.java b/android/ConductorsWatch/app/src/main/java/de/tonifetzer/conductorswatch/WorkerFragment.java index 34f267a..3148786 100644 --- a/android/ConductorsWatch/app/src/main/java/de/tonifetzer/conductorswatch/WorkerFragment.java +++ b/android/ConductorsWatch/app/src/main/java/de/tonifetzer/conductorswatch/WorkerFragment.java @@ -12,8 +12,6 @@ import android.view.WindowManager; import android.widget.TextView; import java.util.Vector; -import de.tonifetzer.conductorswatch.utilities.Utils; - /** * A simple {@link Fragment} subclass. @@ -21,12 +19,12 @@ import de.tonifetzer.conductorswatch.utilities.Utils; * {@link WorkerFragment.OnFragmentInteractionListener} interface * to handle interaction events. */ -public class WorkerFragment extends Fragment implements Metronome.OnMetronomeListener, BpmEstimator.OnBpmEstimatorListener{ +public class WorkerFragment extends Fragment implements Metronome.OnMetronomeListener, Estimator.OnBpmEstimatorListener{ private OnFragmentInteractionListener mListener; private Vector mBpmList; //TODO save to file. - private BpmEstimator mBpmEstimator; + private Estimator mEstimator; private Metronome mMetronome; //private Thread mBpmThread; @@ -58,9 +56,9 @@ public class WorkerFragment extends Fragment implements Metronome.OnMetronomeLis super.onCreate(savedInstanceState); // init bpm estimator and listener - mBpmEstimator = new BpmEstimator(getContext()); - mBpmEstimator.add(this); - //mBpmThread = new Thread(mBpmEstimator, "estThread"); + mEstimator = new Estimator(getContext()); + mEstimator.add(this); + //mBpmThread = new Thread(mEstimator, "estThread"); mBpmList = new Vector(); @@ -99,7 +97,7 @@ public class WorkerFragment extends Fragment implements Metronome.OnMetronomeLis super.onStart(); // start the worker thread for bpm estimator - mBpmEstimator.start(); + mEstimator.start(); // start the worker thread for metronom mMetronomeThread.start(); @@ -110,7 +108,7 @@ public class WorkerFragment extends Fragment implements Metronome.OnMetronomeLis super.onStop(); // stop the worker thread for bpm estimator - mBpmEstimator.stop(); + mEstimator.stop(); // stop the worker thread for metronom mMetronome.stop(); diff --git a/java/src/main/java/Main.java b/java/src/main/java/Main.java index ef82c27..583ff1e 100644 --- a/java/src/main/java/Main.java +++ b/java/src/main/java/Main.java @@ -6,8 +6,7 @@ import java.io.BufferedReader; import java.io.File; import java.io.FileReader; import java.io.IOException; -import java.util.LinkedList; -import java.util.stream.IntStream; + /** * Created by toni on 04/12/17. @@ -17,17 +16,20 @@ public class Main { public static void main(String [ ] args) { File folder = new File("/home/toni/Documents/programme/dirigent/measurements/wearR"); File[] listOfFiles = folder.listFiles(); +/* Utils.ShowPNG windowRaw = new Utils.ShowPNG(); Utils.ShowPNG windowAuto = new Utils.ShowPNG(); Utils.ShowPNG windowPeaksX = new Utils.ShowPNG(); Utils.ShowPNG windowPeaksY = new Utils.ShowPNG(); Utils.ShowPNG windowPeaksZ = new Utils.ShowPNG(); +*/ // iterate trough files in measurements folder for (File file : listOfFiles) { if (file.isFile() && file.getName().contains(".csv")) { AccelerometerWindowBuffer accWindowBuffer = new AccelerometerWindowBuffer(4096, 256); + BpmEstimator bpmEstimator = new BpmEstimator(accWindowBuffer, 0, 5000); //read the file line by line try (BufferedReader br = new BufferedReader(new FileReader(file))) { @@ -47,7 +49,11 @@ public class Main { //do calculation stuff if(accWindowBuffer.isNextWindowReady()){ - AccelerometerInterpolator acInterp = new AccelerometerInterpolator(accWindowBuffer, 6); + double curBpm = bpmEstimator.estimate(); + //System.out.println("BPM: " + curBpm); +/* + + AccelerometerInterpolator acInterp = new AccelerometerInterpolator(accWindowBuffer, bpmEstimator.getSampleRate_ms()); //print raw x,y,z double[] dTs = IntStream.range(0, accWindowBuffer.getTs().length).mapToDouble(i -> accWindowBuffer.getTs()[i]).toArray(); @@ -123,20 +129,22 @@ public class Main { //fill hols improve peaks //estimate bpm between detected peaks - System.out.println("BPM-X: " + pX.getBPM(6)); - System.out.println("BPM-Y: " + pY.getBPM(6)); - System.out.println("BPM-Z: " + pZ.getBPM(6)); + //System.out.println("BPM-X: " + pX.getBPM(bpmEstimator.getSampleRate_ms())); + //System.out.println("BPM-Y: " + pY.getBPM(bpmEstimator.getSampleRate_ms())); + //System.out.println("BPM-Z: " + pZ.getBPM(bpmEstimator.getSampleRate_ms())); - //todo: statistikzeuch und mit matlab vergleichen. - //todo: todos machen. lol //todo: kleiner fenstergrößen testen. so ist doch etwas langsam auf der Uhr. - +*/ int dummyForBreakpoint = 0; } - - } // line is not visible here. + + double meanBPM = bpmEstimator.getMeanBpm(); + double medianBPM = bpmEstimator.getMedianBPM(); + //System.out.println("MEAN BPM: " + Math.round(meanBPM)); + System.out.println("MEDIAN BPM: " + Math.round(medianBPM)); + } catch (IOException e) { e.printStackTrace(); } diff --git a/java/src/main/java/Plot.java b/java/src/main/java/Plot.java index 93cce70..4b82c49 100644 --- a/java/src/main/java/Plot.java +++ b/java/src/main/java/Plot.java @@ -28,8 +28,8 @@ import javax.imageio.ImageIO; */ public class Plot { - public enum Line { NONE, SOLID, DASHED }; - public enum Marker { NONE, CIRCLE, SQUARE, DIAMOND, COLUMN, BAR }; + public enum Line { NONE, SOLID, DASHED } + public enum Marker { NONE, CIRCLE, SQUARE, DIAMOND, COLUMN, BAR } public enum AxisFormat { NUMBER, NUMBER_KGM, NUMBER_INT, TIME_HM, TIME_HMS, DATE, DATETIME_HM, DATETIME_HMS } public enum LegendFormat { NONE, TOP, RIGHT, BOTTOM } diff --git a/java/src/main/java/bpmEstimation/AccelerometerInterpolator.java b/java/src/main/java/bpmEstimation/AccelerometerInterpolator.java index 7806995..e270eea 100644 --- a/java/src/main/java/bpmEstimation/AccelerometerInterpolator.java +++ b/java/src/main/java/bpmEstimation/AccelerometerInterpolator.java @@ -1,7 +1,5 @@ package bpmEstimation; -import utilities.Utils; - import java.util.Arrays; /** @@ -16,10 +14,6 @@ public class AccelerometerInterpolator { public AccelerometerInterpolator(AccelerometerWindowBuffer ab, double sampleRate_ms){ - if(sampleRate_ms <= 0){ - sampleRate_ms = Utils.mean(Utils.diff(ab.getTs())); - } - long size = (ab.getYongest().ts - (ab.getOldest().ts - (long) sampleRate_ms)) / (long) sampleRate_ms; mTsInterp = new long[(int)size]; int j = 0; @@ -48,7 +42,7 @@ public class AccelerometerInterpolator { return mZ; } - public static double[] interpLinear(double[] x, double[] y, double[] xi) throws IllegalArgumentException { + private static double[] interpLinear(double[] x, double[] y, double[] xi) throws IllegalArgumentException { if (x.length != y.length) { throw new IllegalArgumentException("X and Y must be the same length"); } @@ -95,7 +89,7 @@ public class AccelerometerInterpolator { return yi; } - public static double[] interpLinear(long[] x, double[] y, long[] xi) throws IllegalArgumentException { + private static double[] interpLinear(long[] x, double[] y, long[] xi) throws IllegalArgumentException { double[] xd = new double[x.length]; for (int i = 0; i < x.length; i++) { diff --git a/java/src/main/java/bpmEstimation/AccelerometerWindowBuffer.java b/java/src/main/java/bpmEstimation/AccelerometerWindowBuffer.java index 86e5c36..f1d846b 100644 --- a/java/src/main/java/bpmEstimation/AccelerometerWindowBuffer.java +++ b/java/src/main/java/bpmEstimation/AccelerometerWindowBuffer.java @@ -8,23 +8,14 @@ import java.util.ArrayList; */ public class AccelerometerWindowBuffer extends ArrayList { - private int mWindowSize; - private int mOverlapSize; + private static int mWindowSize; + private static int mOverlapSize; private int mOverlapCounter; - private double[] mX; - private double[] mY; - private double[] mZ; - private long[] mTs; public AccelerometerWindowBuffer(int windowSize, int overlap){ - this.mWindowSize = windowSize; - this.mOverlapSize = overlap; + mWindowSize = windowSize; + mOverlapSize = overlap; mOverlapCounter = 1; - - mX = new double[this.mWindowSize]; - mY = new double[this.mWindowSize]; - mZ = new double[this.mWindowSize]; - mTs = new long[this.mWindowSize]; } //TODO: add exception handling. falseArgument if ad has no numeric x,y,z @@ -40,21 +31,12 @@ public class AccelerometerWindowBuffer extends ArrayList { removeRange(0, size() - mWindowSize); } - //TODO: bullshit... jedes mal wieder alle neue anlegen ist ja nur langsam - //update the double arrays. - for (int i = 0; i < size(); ++i) { - mX[i] = get(i).x; - mY[i] = get(i).y; - mZ[i] = get(i).z; - mTs[i] = get(i).ts; - } - ++mOverlapCounter; return r; } public boolean isNextWindowReady(){ - if(size() == mWindowSize && mOverlapCounter > mOverlapSize){ + if((size() == mWindowSize || size() == mWindowSize / 2 || size() == mWindowSize / 4) && mOverlapCounter > mOverlapSize){ mOverlapCounter = 1; return true; @@ -71,19 +53,22 @@ public class AccelerometerWindowBuffer extends ArrayList { } public double[] getX(){ - return mX; + return this.stream().mapToDouble(d -> d.x).toArray(); } public double[] getY(){ - return mY; + return this.stream().mapToDouble(d -> d.y).toArray(); } public double[] getZ(){ - return mZ; + return this.stream().mapToDouble(d -> d.z).toArray(); } public long[] getTs(){ - return mTs; + return this.stream().mapToLong(d -> d.ts).toArray(); } + public int getOverlapSize(){ + return mOverlapSize; + } } \ No newline at end of file diff --git a/java/src/main/java/bpmEstimation/AutoCorrelation.java b/java/src/main/java/bpmEstimation/AutoCorrelation.java index fe9c216..75b3755 100644 --- a/java/src/main/java/bpmEstimation/AutoCorrelation.java +++ b/java/src/main/java/bpmEstimation/AutoCorrelation.java @@ -9,8 +9,8 @@ import java.util.Arrays; */ public class AutoCorrelation { - int mMaxLag; - double[] mCorr; + private static int mMaxLag; + private double[] mCorr; public AutoCorrelation(double[] data, int maxLag){ diff --git a/java/src/main/java/bpmEstimation/BpmEstimator.java b/java/src/main/java/bpmEstimation/BpmEstimator.java new file mode 100644 index 0000000..21c98a8 --- /dev/null +++ b/java/src/main/java/bpmEstimation/BpmEstimator.java @@ -0,0 +1,189 @@ +package bpmEstimation; + +import utilities.MovingFilter; +import utilities.SimpleKalman; +import utilities.Utils; + +import java.util.LinkedList; + +/** + * Created by toni on 17/12/17. + */ +public class BpmEstimator { + + private AccelerometerWindowBuffer mBuffer; + private static double mSampleRate_ms; + + private LinkedList mBpmHistory_X; + private LinkedList mBpmHistory_Y; + private LinkedList mBpmHistory_Z; + + private LinkedList mBpmHistory; + private int mResetCounter; + private static int mResetLimit_ms; + + private MovingFilter mMvg; + //private SimpleKalman mKalman; + + + public BpmEstimator(AccelerometerWindowBuffer windowBuffer, double sampleRate_ms, int resetAfter_ms){ + mBuffer = windowBuffer; + mSampleRate_ms = sampleRate_ms; + + mBpmHistory_X = new LinkedList<>(); + mBpmHistory_Y = new LinkedList<>(); + mBpmHistory_Z = new LinkedList<>(); + + mBpmHistory = new LinkedList<>(); + mResetCounter = 0; + mResetLimit_ms = resetAfter_ms; + + mMvg = new MovingFilter(10); + //mKalman = new SimpleKalman(); + } + + public double estimate(){ + + double sampleRate = mSampleRate_ms; + if(sampleRate <= 0){ + sampleRate = Math.round(Utils.mean(Utils.diff(mBuffer.getTs()))); + } + + AccelerometerInterpolator interp = new AccelerometerInterpolator(mBuffer, sampleRate); + + double[] xAutoCorr = new AutoCorrelation(interp.getX(), 1024).getCorr(); + double[] yAutoCorr = new AutoCorrelation(interp.getY(), 1024).getCorr(); + double[] zAutoCorr = new AutoCorrelation(interp.getZ(), 1024).getCorr(); + + Peaks pX = new Peaks(xAutoCorr, 50, 0.1f, 0, false); + Peaks pY = new Peaks(yAutoCorr, 50, 0.1f, 0, false); + Peaks pZ = new Peaks(zAutoCorr, 50, 0.1f, 0, false); + + mBpmHistory_X.add(pX.getBPM(sampleRate)); + mBpmHistory_Y.add(pY.getBPM(sampleRate)); + mBpmHistory_Z.add(pZ.getBPM(sampleRate)); + + double estimatedBPM = getBestBpmEstimation(pX, pY, pZ); + if(estimatedBPM != -1){ + + //moving avg (lohnt dann, wenn wir viele daten haben) + //mMvg.add(estimatedBPM); + //mBpmHistory.add(mMvg.getAverage()); + + //moving median (lohnt nur bei konstantem tempo, da wir nur ein tempo damit gut halten können.) + //mMvg.add(estimatedBPM); + //mBpmHistory.add(mMvg.getMedian()); + + //kalman filter (lohnt dann, wenn wir konstantes tempo haben, mit startangabe!) + + mBpmHistory.add(estimatedBPM); + } + else { + int resetAfter = (int) Math.round(mResetLimit_ms / (mBuffer.getOverlapSize() * sampleRate)); + if(++mResetCounter > resetAfter){ + mBpmHistory.clear(); + mResetCounter = 0; + } + return -1; + } + + //last element + return mBpmHistory.getLast(); + } + + public double getMeanBpm(){ + return Utils.mean(mBpmHistory); + } + + public double getMedianBPM(){ + return Utils.median(mBpmHistory); + } + + private double getBestBpmEstimation(Peaks peaksX, Peaks peaksY, Peaks peaksZ) throws IllegalArgumentException { + + int cntNumAxis = 0; + double sumCorr = 1; //to prevent division by zero + double sumRms = 1; + double sumNumInter = 1; + + double corrMeanX = 0, corrRmsX = 0; + int corrNumInterX = 0; + if(peaksX.hasPeaks()){ + corrMeanX = Utils.geometricMean(peaksX.getPeaksValueWithoutZeroIdxAndNegativeValues()); + corrRmsX = Utils.rms(peaksX.getPeaksValueWithoutZeroIdx()); + corrNumInterX = Utils.intersectionNumber(peaksX.getData(), 0.2f); + + ++cntNumAxis; + sumCorr += corrMeanX; + sumRms += corrRmsX; + sumNumInter += corrNumInterX; + } + + double corrMeanY = 0, corrRmsY = 0; + int corrNumInterY = 0; + if(peaksY.hasPeaks()){ + corrMeanY = Utils.geometricMean(peaksY.getPeaksValueWithoutZeroIdxAndNegativeValues()); + corrRmsY = Utils.rms(peaksY.getPeaksValueWithoutZeroIdx()); + corrNumInterY = Utils.intersectionNumber(peaksY.getData(), 0.2f); + + ++cntNumAxis; + sumCorr += corrMeanY; + sumRms += corrRmsY; + sumNumInter += corrNumInterY; + } + + double corrMeanZ = 0, corrRmsZ = 0; + int corrNumInterZ = 0; + if(peaksZ.hasPeaks()){ + corrMeanZ = Utils.geometricMean(peaksZ.getPeaksValueWithoutZeroIdxAndNegativeValues()); + corrRmsZ = Utils.rms(peaksZ.getPeaksValueWithoutZeroIdx()); + corrNumInterZ = Utils.intersectionNumber(peaksZ.getData(), 0.2f); + + ++cntNumAxis; + sumCorr += corrMeanZ; + sumRms += corrRmsZ; + sumNumInter += corrNumInterZ; + } + + //no peaks, reject + if(cntNumAxis == 0){ + //throw new IllegalArgumentException("All Peaks are empty! -> Reject Estimation"); + return -1; + } + + /* + System.out.println("RMS-X: " + corrRmsX); + System.out.println("GEO-X: " + corrMeanX); + System.out.println("INTER-X: " + corrNumInterX); + + System.out.println("RMS-Y: " + corrRmsY); + System.out.println("GEO-Y: " + corrMeanY); + System.out.println("INTER-Y: " + corrNumInterY); + + System.out.println("RMS-Z: " + corrRmsZ); + System.out.println("GEO-Z: " + corrMeanZ); + System.out.println("INTER-Z: " + corrNumInterZ); + */ + + //values to low, reject + //TODO: this is a pretty simple assumption. first shot! + if(corrRmsX < 0.2 && corrRmsY < 0.2 && corrRmsZ < 0.2){ + return -1; + } + + double quantityX = ((corrMeanX / sumCorr) + (corrRmsX / sumRms) + (corrNumInterX / sumNumInter)) / cntNumAxis; + double quantityY = ((corrMeanY / sumCorr) + (corrRmsY / sumRms) + (corrNumInterY / sumNumInter)) / cntNumAxis; + double quantityZ = ((corrMeanZ / sumCorr) + (corrRmsZ / sumRms) + (corrNumInterZ / sumNumInter)) / cntNumAxis; + + //get best axis by quantity and estimate bpm + if(quantityX > quantityY && quantityX > quantityZ){ + return mBpmHistory_X.getLast(); + } + else if(quantityY > quantityZ){ + return mBpmHistory_Y.getLast(); + } + else { + return mBpmHistory_Z.getLast(); + } + } +} diff --git a/java/src/main/java/bpmEstimation/Peaks.java b/java/src/main/java/bpmEstimation/Peaks.java index cd7c6b6..6076749 100644 --- a/java/src/main/java/bpmEstimation/Peaks.java +++ b/java/src/main/java/bpmEstimation/Peaks.java @@ -74,15 +74,14 @@ public class Peaks { public double[] getPeaksValueWithoutZeroIdx() { double[] values = new double[mPeaksIdx.size() - 1]; - int i = 0; int mid = (mData.length / 2); - for(Integer idx : mPeaksIdx) { - if(!idx.equals(mid)){ - values[i] = mPeaksValue.get(idx); - ++i; + int j = 0; + for(int i = 0; i < mPeaksIdx.size(); ++i){ + if(!(mPeaksIdx.get(i) == mid)){ + values[j] = mPeaksValue.get(i); + ++j; } } - return values; } @@ -103,16 +102,15 @@ public class Peaks { public double[] getPeaksValueWithoutZeroIdxAndNegativeValues(){ double[] values = new double[mPeaksIdx.size() - 1]; - int i = 0; int mid = (mData.length / 2); - for(Integer idx : mPeaksIdx) { - double curVal = mPeaksValue.get(idx); - if(!idx.equals(mid) && curVal> 0){ - values[i] = curVal; - ++i; + int j = 0; + for(int i = 0; i < mPeaksIdx.size(); ++i){ + double curVal = mPeaksValue.get(i); + if(!(mPeaksIdx.get(i) == mid) && curVal > 0){ + values[j] = curVal; + ++j; } } - return values; } @@ -131,7 +129,8 @@ public class Peaks { if(hasPeaks()){ //todo diff and mean method for linkedlists for speed - return 60000 / Utils.mean(Utils.diff(mPeaksPos.stream().mapToDouble(i -> i * sampleRate_ms).toArray())); + //return 60000 / Utils.mean(Utils.diff(mPeaksPos.stream().mapToDouble(i -> i * sampleRate_ms).toArray())); + return 60000 / (sampleRate_ms * Utils.mean(Utils.diff(mPeaksPos.stream().mapToDouble(i -> i).toArray()))); } return -1; @@ -154,7 +153,7 @@ public class Peaks { //TODO: findPeaks method identical to Matlab... with PeakProminence private void simplePeakFinder(double[] data, int width, double threshold, double decayRate, boolean isRelative) { - int maxp = 0; + int maxp; int mid = 0; int end = data.length; double av = data[0]; diff --git a/java/src/main/java/utilities/MovingFilter.java b/java/src/main/java/utilities/MovingFilter.java new file mode 100644 index 0000000..e4f883d --- /dev/null +++ b/java/src/main/java/utilities/MovingFilter.java @@ -0,0 +1,33 @@ +package utilities; + +import java.util.LinkedList; + +/** + * Created by toni on 18/12/17. + */ +public class MovingFilter { + + private static int mSize; + private double mTotal = 0d; + private LinkedList mSamples = new LinkedList<>(); + + public MovingFilter(int size) { + this.mSize = size; + } + + public void add(double x) { + mTotal += x; + mSamples.add(x); + if(mSamples.size() > mSize){ + mTotal -= mSamples.remove(); + } + } + + public double getAverage() { + return mTotal / mSamples.size(); + } + + public double getMedian() { + return Utils.median(mSamples); + } +} diff --git a/java/src/main/java/utilities/SimpleKalman.java b/java/src/main/java/utilities/SimpleKalman.java new file mode 100644 index 0000000..54f858a --- /dev/null +++ b/java/src/main/java/utilities/SimpleKalman.java @@ -0,0 +1,33 @@ +package utilities; + +/** + * Created by toni on 18/12/17. + */ +public class SimpleKalman { + + private double mSigmaUpdate; + private double mSigmaPrediction; + private double mMu; + private double mSigma; + + SimpleKalman(double initialMu, double initialSigma, double sigmaUpdate, double sigmaPrediction){ + mSigmaUpdate = sigmaUpdate; + mSigmaPrediction = sigmaPrediction; + mMu = initialMu; + mSigma = initialSigma; + } + + public double update(double data){ + + //prediction + mMu = mMu; + mSigma += mSigmaPrediction; + + //update + double k = mSigma / (mSigma + mSigmaUpdate); + mMu = mMu + k * (data - mMu); + mSigma = (1 - k) * mSigma; + + return mMu; + } +} diff --git a/java/src/main/java/utilities/Utils.java b/java/src/main/java/utilities/Utils.java index 8464e31..68b284f 100644 --- a/java/src/main/java/utilities/Utils.java +++ b/java/src/main/java/utilities/Utils.java @@ -4,7 +4,11 @@ import bpmEstimation.Peaks; import javax.swing.*; import java.awt.*; import java.awt.image.BufferedImage; +import java.text.Collator; import java.util.Arrays; +import java.util.Collections; +import java.util.Comparator; +import java.util.LinkedList; //TODO: change from double to generic type public class Utils { @@ -36,12 +40,20 @@ public class Utils { return sum; } + public static double sum (LinkedList data){ + double sum = 0; + for (int i = 0; i < data.size(); i++) { + sum += data.get(i).doubleValue(); + } + return sum; + } + //TODO: Could be slow.. faster method? public static double rms(double[] nums) { double sum = 0.0f; for (double num : nums) sum += num * num; - return (double) Math.sqrt(sum / nums.length); + return Math.sqrt(sum / nums.length); } public static double[] diff(double[] data){ @@ -68,8 +80,24 @@ public class Utils { return sum(data) / data.length; } - public static long mean(long[] data){ - return sum(data) / data.length; + public static double mean(long[] data){ + return (double) sum(data) / (double) data.length; + } + + public static double mean(LinkedList data){ + return sum(data) / data.size(); + } + + public static double median(LinkedList data){ + data.sort(Comparator.naturalOrder()); + + double median; + if (data.size() % 2 == 0) + median = (data.get(data.size()/2) + data.get(data.size()/2 - 1))/2; + else + median = data.get(data.size()/2); + + return median; } //TODO: Could be slow.. faster method? @@ -78,7 +106,7 @@ public class Utils { for (int i = 1; i < data.length; i++) { sum *= data[i]; } - return (double) Math.pow(sum, 1.0 / data.length); + return Math.pow(sum, 1.0 / data.length); } public static int intersectionNumber(double[] signal, double border){ @@ -108,76 +136,6 @@ public class Utils { return cnt; } - public static double getBestBpmEstimation(Peaks peaksX, Peaks peaksY, Peaks peaksZ) throws IllegalArgumentException { - - int cntNumAxis = 0; - double sumCorr = 1; //to prevent division by zero - double sumRms = 1; - double sumNumInter = 1; - - double corrMeanX = 0, corrRmsX = 0; - int corrNumInterX = 0; - if(peaksX.hasPeaks()){ - corrMeanX = geometricMean(peaksX.getPeaksValueWithoutZeroIdxAndNegativeValues()); - corrRmsX = rms(peaksX.getPeaksValueWithoutZeroIdx()); - corrNumInterX = intersectionNumber(peaksX.getData(), 0.2f); - - ++cntNumAxis; - sumCorr += corrMeanX; - sumRms += corrRmsX; - sumNumInter += corrNumInterX; - } - - double corrMeanY = 0, corrRmsY = 0; - int corrNumInterY = 0; - if(peaksY.hasPeaks()){ - corrMeanY = geometricMean(peaksY.getPeaksValueWithoutZeroIdxAndNegativeValues()); - corrRmsY = rms(peaksY.getPeaksValueWithoutZeroIdx()); - corrNumInterY = intersectionNumber(peaksY.getData(), 0.2f); - - ++cntNumAxis; - sumCorr += corrMeanY; - sumRms += corrRmsY; - sumNumInter += corrNumInterY; - } - - double corrMeanZ = 0, corrRmsZ = 0; - int corrNumInterZ = 0; - if(peaksZ.hasPeaks()){ - corrMeanZ = geometricMean(peaksZ.getPeaksValueWithoutZeroIdxAndNegativeValues()); - corrRmsZ = rms(peaksZ.getPeaksValueWithoutZeroIdx()); - corrNumInterZ = intersectionNumber(peaksZ.getData(), 0.2f); - - ++cntNumAxis; - sumCorr += corrMeanZ; - sumRms += corrRmsZ; - sumNumInter += corrNumInterZ; - } - - if(cntNumAxis == 0){ - throw new IllegalArgumentException("All Peaks are empty! -> Reject Estimation"); - } - - double quantityX = ((corrMeanX / sumCorr) + (corrRmsX / sumRms) + (corrNumInterX / sumNumInter)) / cntNumAxis; - double quantityY = ((corrMeanY / sumCorr) + (corrRmsY / sumRms) + (corrNumInterY / sumNumInter)) / cntNumAxis; - double quantityZ = ((corrMeanZ / sumCorr) + (corrRmsZ / sumRms) + (corrNumInterZ / sumNumInter)) / cntNumAxis; - - //get best axis by quantity and estimate bpm - if(quantityX > quantityY && quantityX > quantityZ){ - - } - else if(quantityY > quantityZ){ - - } - else { - - } - - //check if values are empty. - - return 0.0f; - } - public static double[] removeZero(double[] array){ int j = 0; for( int i=0; i