[Research] Window peak detection algorithm for step counting

1. Introduction

Location-based services become of more importance in people's lives with the advancement of the internet of things (IoT). How to provide continuous, indoor and outdoor location service attracts more and more people' attention. Pedestrian Dead-Reckoning (PDR) system is one of the most important positioning and navigation system because it only needs the IMU sensors, which is more convenient than other techniques that either need infrastructure or GPS.

For a typical PDR system, the a continuous walk detection algorithm and a robust estimate of the number of steps taken are two vital parts. Apart from these, the heading and the matric are also important. Here, firstly I will introduce some start-of-arts works about the walk detection and step counting. Then, I will show you one popular step counting algorithm, peak detection method.

2. State-of-art works

1) Algorithms classifications

• Thresholding techniques: As with any thresholding technique, difficulties arise in selecting the optimal threshold value, which can vary between users, surfaces, and shoes. For an unconstrained smartphone, the freedom of movement might additionally affect the property of the accelerometer signal and mistakenly trigger the threshold.
• Detecting a stride (two steps): Typical stride frequencies are around 1–2 Hz, a range that few activities other than walking exhibit.
• Peak detection [25, 48] or zero-crossing counting on low-pass filtered accelerometer signals can be used to find specific steps. 
• The short-term Fourier transform (STFT) can be applied to evaluate the frequency content of successive windows of data. However, this windowing introduces resolution issues. 
• Wavelet transforms (CWT/DWT), which repeatedly correlate a ‘mother’ wavelet with the signal, temporally compressing or dilating as appropriate, are more suitable as they can capture sudden changes in the acceleration, but are known for being computationally more expensive.
• auto-correlation to detect the period directly in the time domain. 
• dynamic time warping (DTW) 

2) Feature

Commonly used features include mean, variance or standard deviation, energy, entropy, correlation between axes, and FFT coefficients and Features obtained by dimensionality reduction (e.g., principal component analysis).

3) Classification methods

Most studies use these features for training classifiers such as neural networks, support vector machines (SVMs) or Gaussian mixture models (GMMs), or give them as an input to non-parametric classification methods such as k-nearest neighbour. These techniques can be effective (particularly the Bayesian ones), but are supervised, so offline training is required—often for each individual and each carrying position.

Unsupervised techniques are used much less (notably, self-organizing maps and k-means clustering) with the exception of hidden Markov models (HMMs), which can be used both for classification and in an unsupervised fashion. HMMs are particularly interesting as they are trained on a sequence of features, which tends to improve the accuracy as temporal correlations get captured in the model. 

3. Window peak detection algorithm

1). Window peak detection model

A cycle of human walking is composed of two phases: standing and walking phase. The step detection means a recognition of the walking phase. The walking phase is divided into a swing phase and a heel-touch-down phase. Each phase is shown in figure 1.

Figure 1. The walking phase

In the 1st swing phase, the foot is located on behind of gravity center of the human body. And the foot is located on front of gravity centre of human body in the 2nd swing phase. The foot accelerated during swing phase. The acceleration is composed of vertical and horizontal components as shown in figure 2, where a , h , g means horizontal acceleration, vertical acceleration and gravity force, respectively.

Figure 2. components of the acceleration

Then we deduce the relationship between those three components. In figure 3 and 4, θ(t) is inclination angle of the leg at time t and it has two situations, at the left the vertical direction and at the right of the vertical direction. In the 1st swing phase, people lift their leg and step forward. And in the 2nd swing phase, people step forward continuously. We have equations as shown below:

Figure 3. 1st swing phase

Figure 4.2nd swing phase

Formula. 1

We adopted common assumptions that a typical inclination of leg was within the limit of 30 degree ~ 50 degree and a , h have a range of 0.8 ~ 2.3g and 0.6 ~2.0g.
In most PDR system, a vertical impact is measured to detect a step. A step is counted when the measured vertical impact is larger than the given threshold. A vertical impact results in the peak of the vertical acceleration. Figure 5 and 6 show the signal pattern of one walking cycle.

Figure.5 Vertical acceleration signal pattern in walking phase

Figure.6 horizontal acceleration signal pattern in walking phase

From fig,5 and 6, we can see at the moment of the leg touching the ground after two swing phases, the peak occurs on the vertical acceleration direction. That is very important as the principle of this window peak detection method is using this peak to count steps as shown below:

Figure.7 peak detection example

2). Step detection method

From the aforesaid, the process we count steps is to detect the peak of the vertical acceleration. The distance between two peaks is one step. In practical, vertical acceleration is not easy to get as the smartphone is not firmly attached to a human body. Generally, we directly use the mode of the three-axis accelerometer.

# -*- coding: utf-8 -*-
"""
Created on Mon Jun 25 20:56:54 2018

@author: mzx bwu
"""
from read_log_file_IPIN2018 import read_log_file_IPIN2018
from step_counter import detect_peak
import pandas as pd
import matplotlib.pyplot as plt
import numpy as np
#import math
import os
import natsort

# =============================================================================
# there are two kinds of parameters: active parameters and functional parameters
# active parameters are those with direct affection with the step counting 
value_num = 4    #the size of threshold list
time_interval = 0.25   #the timeinterval between two peaks
threshold_value = 2.5     #the threshold for the difference between valley and peak
initialValue=1.8         #the pre-set value for the difference between valley and peak
                         #it is used to decide whether to update the threshold_value
peak_threshold=20  #it is used to deal with a sudden increase of the acceleration even if it doesn't increase twice.                        
#and the parameters for computing the average_value.

# functional parameters are those just initiate the parameters. those parameters are functional ones
is_direction_up = False
continue_up_count = 0
continue_up_former_count=0
valley_of_wave = 0 
peak_of_wave = 0
last_status = False
tmp_count=0
tmp_value = []
old_measurement = 0
time_last_peak=0

def detect_peak(new_measurement, old_measurement):
    global is_direction_up,continue_up_count,continue_up_former_count, valley_of_wave, peak_of_wave 
    last_status = is_direction_up   #default is decreasing
    if new_measurement >= old_measurement:
        is_direction_up = True
        continue_up_count += 1
    else:
        continue_up_former_count = continue_up_count
        continue_up_count = 0
        is_direction_up = False
# =============================================================================
# The peak detection      
# there are four conditions:     
# The status of current point is decreasing => isDirectionUp is false      
# The last status is increasing =>lastStatus is true      
# Until this peak point, the times of increasing continuously are at least two 
# If run across the sudden increase of accelerometer, the peak must be beyond the peak_threshold (generally set 20)
# =============================================================================
    if last_status and (not is_direction_up) and (continue_up_former_count >= 2 or old_measurement >= peak_threshold): 
        peak_of_wave = old_measurement
        return  True
    elif (not last_status) and is_direction_up: #last status is direction down, and now is up,so it is valley
        valley_of_wave = old_measurement      
        return False
    else:
        return False

# =============================================================================
# the threshold unpdating of the difference between the peak and the valley 
# choose past value_num thresholds and get the mean
# =============================================================================
def peak_vally_thread(value):
    global tmp_count
    tmp_thread = threshold_value
    if(tmp_count < value_num):
        tmp_value.append(value)
        tmp_count += 1
    else:
        tmp_thread = average_value(tmp_value)
        for i in range(1,value_num):
            tmp_value[i-1] = tmp_value[i]
        tmp_value[value_num-1]=value
    return tmp_thread
        
# =============================================================================
# the gradient of threshold
# Aim is to set the threshold in a certain period.
# but now I don't know how to train these parameters.
# =============================================================================
def average_value(value_list):
    ave = np.mean(value_list)
    if ave >= 5 :
        ave = 5.1
    elif 4 <= ave < 5:
        ave = 4.1
    elif 3 <= ave < 4:
        ave = 3.1
    elif 2 <= ave < 3:
        ave = 2.1
    else:
        ave = 1.1
    return ave

if __name__ == '__main__':    #function interface
    file_dir = '../Data/ipin2018track3-logfiles/Logfiles/01 - Training'
    files = os.listdir(file_dir) 
    files = natsort.natsorted(files) 
    i=1
    file_name = file_dir+'/'+files[i]
    #print(files[i])
    
    # list order 0 posi, 1 acce 2 gyro 3 magn 4 pres 5 ligh 6 prox 7 soun 8 ahrs 9 gnss 10 rfid 11 wifi 12 imul 13 imux 14 ble4    
    raw_total_data_dic = read_log_file_IPIN2018(file_name)
    b=pd.DataFrame(raw_total_data_dic['acce'])
    #plt.plot(b[0],b[2])
    #plt.plot(b[0],b[3])
    #plt.plot(b[0],b[4])
    plt.plot(b[0],np.sqrt(b[2]**2+b[3]**2+b[4]**2))
    acce_mag = np.sqrt(b[2]**2+b[3]**2+b[4]**2)

    my_count = 0    #record the steps
    for i in range(1484,2964):#len(acce_mag)             
        if old_measurement == 0:
             old_measurement = acce_mag[i]
             time_last_peak = b[0][i]    
        else:
            if detect_peak(acce_mag[i], old_measurement):
                time_now = b[0][i]
    #                 print(time_now-time_last_peak)
    # =============================================================================
    #              The process of step counting
#                    First detect the valley and use this valley to find the next peak
#                    when find the peak, obtain the difference between the peak and the valley.
# condition 1: the time difference between two close peaks should be beyond the time_interval 
                   # condition 2: The difference bwteen the peak and the valley should beyond the threshold_value               
    # =============================================================================
                if ((time_now-time_last_peak) >= time_interval)&((peak_of_wave-valley_of_wave)>=threshold_value):   
                    my_count +=1
                    #print(my_count)
                    time_last_peak = b[0][i]
    # =============================================================================
                   #Only the difference bwteen the valley and the peak is beyond the initiation value, we update the threshold value. Or it will lower down the threshold.
    # =============================================================================                
                if ((time_now-time_last_peak) >= time_interval)&((peak_of_wave-valley_of_wave)>=initialValue):   
                    time_last_peak = b[0][i]
                    threshold_value = peak_vally_thread(peak_of_wave-valley_of_wave)  
        old_measurement = acce_mag[i]

Reference

Paper http://user.ceng.metu.edu.tr/~tcan/se542_f1314/Schedule/week13_papers/WD_SC_evaluationPaper.pdf
Opensource dataset: https://www.repository.cam.ac.uk/handle/1810/266947