diff --git a/LICENSE b/LICENSE index 8dada3e..8653bf3 100644 --- a/LICENSE +++ b/LICENSE @@ -186,7 +186,7 @@ same "printed page" as the copyright notice for easier identification within third-party archives. - Copyright {yyyy} {name of copyright owner} + Copyright 2017 LucaCerina Licensed under the Apache License, Version 2.0 (the "License"); you may not use this file except in compliance with the License. diff --git a/README.md b/README.md new file mode 100644 index 0000000..630c0c5 --- /dev/null +++ b/README.md @@ -0,0 +1,15 @@ +README.md rev. 09 May 2017 by Luca Cerina. +Copyright (c) 2017 Luca Cerina. +Distributed under the Apache 2.0 License in the accompanying file LICENSE. + +# Automatic Multiscale-based Peak Detection (AMPD) + +ampdLib implements automatic multiscale-based peak detection (AMPD) algorithm +as in An Efficient Algorithm for Automatic Peak Detection in Noisy Periodic and +Quasi-Periodic Signals, by Felix Scholkmann, Jens Boss and Martin Wolf, +Algorithms 2012, 5, 588-603. + +### Python required dependencies +- Python > 3 +- Numpy + diff --git a/ampdLib.py b/ampdLib.py new file mode 100644 index 0000000..86ceefc --- /dev/null +++ b/ampdLib.py @@ -0,0 +1,98 @@ +# -*- coding: utf-8 -*- +__author__ = "Luca Cerina" +__copyright__ = "Copyright 2017, Luca Cerina" +__email__ = "luca.cerina@polimi.it" + + +import numpy as np + +# AMPD function +def ampd(sigInput): +"""Find the peaks in the signal with the AMPD algorithm + + Parameters + ---------- + sigInput: ndarray + The 1D signal given as input to the algorithm + + Returns + ------- + pks: ndarray + The ordered array of peaks found in sigInput +""" + + sigInput = sigInput.reshape(len(sigInput), 1) + + # Create preprocessing linear fit + sigTime = np.arange(0, len(sigInput)) + + fitPoly = np.polyfit(sigTime, sigInput, 1) + sigFit = np.polyval(fitPoly, sigTime) + + # Detrend + dtrSignal = sigInput - sigFit + + N = len(dtrSignal) + L = int(np.ceil(N / 2.0)) - 1 + + # Generate random matrix + LSM = np.random.uniform(1.0, 2.0, size = (L,N)) # uniform + alpha = 1 + + # Local minima extraction + for k in np.arange(1, L): + for i in range(k+1, N-k): + if((sigInput[i-1, 0] > sigInput[i-k-1, 0]) + & (sigInput[i-1, 0] > sigInput[i+k-1, 0])): + LSM[k-1, i-1] = 0 + + # Find minima + G = np.sum(LSM, 1) + l = np.where(G == G.min())[0] + + LSM = LSM[0:l, :] + + S = np.std(LSM, 0) + + pks = np.flatnonzero(S == 0) + return pks + +# Fast AMPD +def ampdFast(sigInput, order): +"""A slightly faster version of AMPD which divides the signal in 'order' windows + + Parameters + ---------- + sigInput: ndarray + The 1D signal given as input to the algorithm + order: int + The number of windows in which sigInput is divided + + Returns + ------- + pks: ndarray + The ordered array of peaks found in sigInput +""" + + # Check if order is valid (perfectly separable) + if(len(sigInput)%order != 0): + print("AMPD: Invalid order, decreasing order") + while(len(sigInput)%order != 0): + order -= 1 + print("AMPD: Using order " + str(order)) + + N = int(len(sigInput) / order / 2) + + # Loop function calls + for i in range(0, len(sigInput)-N, N): + print("\t sector: " + str(i) + "|" + str((i+2*N-1))) + pksTemp = ampd(sigInput[i:(i+2*N-1)]) + if(i == 0): + pks = pksTemp + else: + pks = np.concatenate((pks, pksTemp+i)) + + # Keep only unique values + pks = np.unique(pks) + + return pks + \ No newline at end of file