IMU & Step Counting
Lab 4
Lab 4
In this lab, we will be adding an Inertial Measurement Unit (IMU) to our wearable prototype to implement a step counter. The IMU uses a combination of an accelerometer and a gyroscope to allow us to obtain useful calculations to detect steps, user movement, etc. Essentially, the IMU is a integral part of the final prototype that we have been building towards in Labs 1 to 3. Furthermore, we will be using the various types of signal manipulation, such as filtering, in order to shape our signals to optimize the functionality of our wearable prototype.
In the last portion of the lab, we will implement a step detection algorithm using scipy's peak finding methods. Using this algorithm, we transfer from Serial communication to BLE between the Arduino and Python. The final product will be a system that is capable of detecting steps, and displays the current step count on the Arduino's OLED display. A pause/play button will also be implemented to control the Arduino's data output to the Python script.
1.) We need to install the TimerOne and i2cDev libraries. This is done by copying the ZIP files of each library to the Arduino's libraries directory:
../Arduino/libraries
Then you can add the libraries in the Arduino IDE by going to:
Arduino >> Sketch >> Include Library >> Add .ZIP Library
2.) Install the Arduino code available below. Note that the code specifically requires the Arduino to read a "1" from Serial to start printing the values and a "0" from Serial to stop. Set up the circuit by connecting the prototype shield to the Arduino board and use Serial Monitor to write "1" to the Arduino to read the IMU values.
#include "I2Cdev.h"
//#include "MPU6050.h"
#include "MPU6050_6Axis_MotionApps20.h"
#include "Wire.h"
#include "TimerOne.h"
const int MPU_addr=0x68; // I2C address of the MPU-6050, can be changed to 0x69
MPU6050 IMU(MPU_addr);
const int interruptPin = 2;
volatile bool ipinReady = false;
int16_t ax, ay, az, tp, gx, gy, gz;
int samplePeriod = 5000; // the target sample period, 5000 microseconds, 200Hz
unsigned long startTime = 0;
unsigned long volatile currentTime = 0;
bool newRead = false;
bool sending = false;
// Timer interrupt to sample our data
void printDataISR() {
// Saves the time
currentTime = (micros()-startTime)/1e3;
newRead = true;
}
// Check the interrupt pin to see if there is data available in the MPU's buffer
void interruptPinISR() {
ipinReady = true;
}
void setup(){
// Intialize the IMU and the DMP ont he IMU
IMU.initialize();
IMU.dmpInitialize();
IMU.setDMPEnabled(true);
Wire.begin();
Wire.beginTransmission(MPU_addr);
Wire.write(MPU_addr); // PWR_MGMT_1 register
Wire.write(0); // set to zero (wakes up the MPU-6050)
Wire.endTransmission(true);
Serial.begin(9600);
while (!Serial);
// Set the timer interrupt
startTime = micros();
Timer1.initialize(samplePeriod);
Timer1.attachInterrupt(printDataISR);
// Create an interrupt for pin2, which is connected to the INT pin of the MPU6050
pinMode(interruptPin, INPUT);
attachInterrupt(digitalPinToInterrupt(interruptPin), interruptPinISR, RISING);
}
void loop(){
// toggle on/off so the buffer doesn't go crazy
// IMPORTANT: Send a '1' via Serial in order to begin reading the data, and a '0' to stop sending data
if (Serial.available()) {
int incomingByte = Serial.read();
if (incomingByte == '1')
sending = true;
if (incomingByte == '0')
sending = false;
}
if (newRead && sending && ipinReady) {
Wire.beginTransmission(MPU_addr);
Wire.write(0x3B); // starting with register 0x3B (ACCEL_XOUT_H)
Wire.endTransmission(false);
Wire.requestFrom(MPU_addr,14,true); // request a total of 14 registers
// Bit Wire.read() gets 8 bits. Every axis reading is 16 bits. Threfore we read two consecutive bytes
// from Wire, and store them next to each other as [ byte1 byte2 ]. This can be done by left bit shifting
// byte1 and ORing it with byte2
//Accelerometer (3 Axis)
ax=Wire.read()<<8|Wire.read(); // 0x3B (ACCEL_XOUT_H) & 0x3C (ACCEL_XOUT_L)
ay=Wire.read()<<8|Wire.read(); // 0x3D (ACCEL_YOUT_H) & 0x3E (ACCEL_YOUT_L)
az=Wire.read()<<8|Wire.read(); // 0x3F (ACCEL_ZOUT_H) & 0x40 (ACCEL_ZOUT_L)
//Temperature
tp=Wire.read()<<8|Wire.read(); // 0x41 (TEMP_OUT_H) & 0x42 (TEMP_OUT_L)
//Gyroscope (3 Axis)
gx=Wire.read()<<8|Wire.read(); // 0x43 (GYRO_XOUT_H) & 0x44 (GYRO_XOUT_L)
gy=Wire.read()<<8|Wire.read(); // 0x45 (GYRO_YOUT_H) & 0x46 (GYRO_YOUT_L)
gz=Wire.read()<<8|Wire.read(); // 0x47 (GYRO_ZOUT_H) & 0x48 (GYRO_ZOUT_L)
// Displays the raw value and the time it was sampled
Serial.print(currentTime);
Serial.print(" ");
Serial.print(ax);
Serial.print(" ");
Serial.print(ay);
Serial.print(" ");
Serial.print(az);
Serial.print(" ");
Serial.print(gx);
Serial.print(" ");
Serial.print(gy);
Serial.print(" ");
Serial.print(gz);
//Serial.print(" ");
//Serial.print(tp); // Feel free to disregard the temp reading
Serial.println("");
newRead = false;
ipinReady = false;
}
}
3.) We need a Python script to read the values from Arduino via Serial and live plot them in one figure with 6 subplots. This is done in the following code pasted below:
import scipy.signal as sig
import serial
from matplotlib import pyplot as plt
from matplotlib import animation
from time import sleep
import numpy as np
# import filtering as flt
# ==================== Globals ====================
N = 400 # number of samples plotted
NS = 20 # number of samples to read per iteration
sample_count = 0 # current sample count
# filter_order = 3
# LPF_cutoff = 10.0/200
# HPF_cutoff = 0.5/200
# gain = 3600
# offset = 2.5
# Vcc = 5
times, ax_values, ay_values, az_values, gx_values, gy_values, gz_values, processed = np.zeros((8, N)) # data vectors
# filter_ICs = np.zeros((2, filter_order))
# filter_coeffs = np.zeros((4, filter_order + 1))
serial_port = '/dev/cu.usbmodem1411' # the serial port to use
serial_baud = 9600 # baudrate
def grab_samples(n_samples) :
global sample_count
# Create local arrays of size(n_samples) to store new data
times, ax_values, ay_values, az_values, gx_values, gy_values, gz_values = np.zeros((7,n_samples))
i = 0
ser.write(b'1') # Tell Arduino to start writing
while i < n_samples : # while we still need samples, keep reading
try:
ser.write(b'1')
data = ser.readline().strip().decode('utf-8')
t, ax, ay, az, gx, gy, gz = data.split(' ')
t = float(t)
ax = float(ax)
ay = float(ay)
az = float(az)
gx = float(gx)
gy = float(gy)
gz = float(gz)
except : # report error if we failed
print('Invalid data: ', data)
continue
# Store the new values into appropriate local arrays
times[i] = t / 1000 # convert time to seconds from ms
ax_values[i] = ax
ay_values[i] = ay
az_values[i] = az
gx_values[i] = gx
gy_values[i] = gy
gz_values[i] = gz
i += 1
ser.write(b'0')
sample_count += n_samples # Increment sample_count by amount of n_samples
return times, ax_values, ay_values, az_values, gx_values, gy_values, gz_values # Return two arrays (size n_sample) with new values
def update_plots(i):
global times, ax_values, ay_values, az_values, gx_values, gy_values, gz_values
# global processed, filter_coeffs, filter_ICs
# shift samples left by 'NS'
times[:N-NS] = times[NS:]
# values[:N-NS] = values[NS:]
ax_values[:N-NS] = ax_values[NS:]
ay_values[:N-NS] = ay_values[NS:]
az_values[:N-NS] = az_values[NS:]
gx_values[:N-NS] = gx_values[NS:]
gy_values[:N-NS] = gy_values[NS:]
gz_values[:N-NS] = gz_values[NS:]
# grab new samples
times[N-NS:], ax_values[N-NS:], ay_values[N-NS:], az_values[N-NS:], gx_values[N-NS:], gy_values[N-NS:], gz_values[N-NS:] = grab_samples(NS)
# processed, filter_ICs = flt.process_ir(values, filter_coeffs, filter_ICs, False)
# plot
[ax.set_xlim(times[0],times[N-1]) for ax in axes]
live_plots[0].set_data(times, ax_values)
live_plots[1].set_data(times, ay_values)
live_plots[2].set_data(times, az_values)
live_plots[3].set_data(times, gx_values)
live_plots[4].set_data(times, gy_values)
live_plots[5].set_data(times, gz_values)
plt.autoscale(enable=True, axis='y', tight=None)
return live_plots
if (__name__ == "__main__"):
# Open Serial
with serial.Serial(port=serial_port, baudrate=serial_baud, timeout=1) as ser:
try:
ser.flush()
sleep(2) # tell Arduino to start sending data. NOTE!!: you should change this for your setup
# initialize the figure
fig, axes = plt.subplots(2, 3)
times, ax_values, ay_values, az_values, gx_values, gy_values, gz_values = grab_samples(N)
# # get HPF and LPF coefficients
# filter_coeffs[3, :], filter_coeffs[2, :] = sig.butter(filter_order, LPF_cutoff, btype='lowpass', analog=False)
# filter_coeffs[1, :], filter_coeffs[0, :] = sig.butter(filter_order, HPF_cutoff, btype='highpass', analog=False)
# # get HPF and LPF initial conditions
# filter_ICs[1, :] = sig.lfilter_zi( filter_coeffs[3, :], filter_coeffs[2, :] )
# filter_ICs[0, :] = sig.lfilter_zi( filter_coeffs[1, :], filter_coeffs[0, :] )
# processed, filter_ICs = flt.process_ir(values, filter_coeffs, filter_ICs, True)
live_plots = []
axes = axes.flatten()
live_plots.append(axes[0].plot(times, ax_values, lw=2)[0])
live_plots.append(axes[1].plot(times, ay_values, lw=2)[0])
live_plots.append(axes[2].plot(times, az_values, lw=2)[0])
live_plots.append(axes[3].plot(times, gx_values, lw=2)[0])
live_plots.append(axes[4].plot(times, gy_values, lw=2)[0])
live_plots.append(axes[5].plot(times, gz_values, lw=2)[0])
plt.autoscale(enable=True, axis='y', tight=None)
# initialize the y-axis limits and labels
# axes[0].set_ylim(0, 1200)
axes[0].set_title('AX')
# axes[0].set_xlabel('Time (s)')
axes[0].set_ylabel('Amplitude')
# axes[1].set_ylim(-0.0005, 0.0005)
axes[1].set_title('AY')
# axes[1].set_xlabel('Time (s)')
# axes[1].set_ylabel('Amplitude (V)')
# axes[2].set_ylim(0, 1200)
axes[2].set_title('AZ')
# axes[2].set_xlabel('Time (s)')
# axes[2].set_ylabel('Amplitude (V)')
# axes[3].set_ylim(0, 1200)
axes[3].set_title('GX')
# axes[3].set_xlabel('Time (s)')
axes[3].set_ylabel('Amplitude')
# axes[4].set_ylim(0, 1200)
axes[4].set_title('GY')
axes[4].set_xlabel('Time (s)')
# axes[4].set_ylabel('Amplitude (V)')
# axes[5].set_ylim(0, 1200)
axes[5].set_title('GZ')
# axes[5].set_xlabel('Time (s)')
# axes[5].set_ylabel('Amplitude (V)')
# set and start the animation and update at 1ms interval (if possible)
anim = animation.FuncAnimation(fig, update_plots, interval= 1)
plt.show()
finally:
ser.flush()
# sleep(2)
ser.write(b'0')
Remember to change the serial_port variable to the one corresponding to your Arduino. Once you run the Python code, you should obtain a plot similar to the one below:
Live Plots of Values from IMU via Serial
1. We will use the Python script below to implement a rough draft of the step counter:
import scipy.signal as sig
import serial
from matplotlib import pyplot as plt
from matplotlib import animation
from time import sleep
import numpy as np
# import filtering as flt
# ==================== Globals ====================
N = 400 # number of samples plotted
NS = 20 # number of samples to read per iteration
sample_count = 0 # current sample count
filter_order = 3
LPF_cutoff = 10.0/200
HPF_cutoff = 0.5/200
thresh = 800
step = 0
# gain = 3600
# offset = 2.5
# Vcc = 5
times, ax_values, ay_values, az_values, gx_values, gy_values, gz_values, processed = np.zeros((8, N)) # data vectors
ax_flt, ay_flt = np.zeros((2,N))
filter_ICs = np.zeros((2, filter_order))
filter_coeffs = np.zeros((4, filter_order + 1))
serial_port = '/dev/cu.usbmodem1411' # the serial port to use
serial_baud = 9600 # baudrate
def filterSignal(data, filter_coeffs, filter_ICs, first):
# ==================== Step 1 ====================
if (first == True) :
LPF, zi_new_LPF = sig.lfilter(filter_coeffs[3, :], filter_coeffs[2, :], data, zi = filter_ICs[1, :]*data[0])
filter_ICs[1, :] = zi_new_LPF[:]
else :
LPF, zi_new_LPF = sig.lfilter(filter_coeffs[3, :], filter_coeffs[2, :], data, zi = filter_ICs[1, :])
filter_ICs[1, :] = zi_new_LPF[:]
# ==================== Step 2 ====================
# Apply HPF filter
# if (first == True) :
# proc_data, zi_new_HPF = sig.lfilter(filter_coeffs[1, :], filter_coeffs[0, :], LPF, zi = filter_ICs[0, :]*LPF[0])
# filter_ICs[0, :] = zi_new_HPF[:]
# else :
# proc_data, zi_new_HPF = sig.lfilter(filter_coeffs[1, :], filter_coeffs[0, :], LPF, zi = filter_ICs[0, :])
# filter_ICs[0, :] = zi_new_HPF[:]
# return proc_data, filter_ICs
return LPF, filter_ICs
def grab_samples(n_samples) :
global sample_count
# Create local arrays of size(n_samples) to store new data
times, ax_values, ay_values, az_values, gx_values, gy_values, gz_values = np.zeros((7,n_samples))
i = 0
ser.write(b'1') # Tell Arduino to start writing
while i < n_samples : # while we still need samples, keep reading
try:
ser.write(b'1')
data = ser.readline().strip().decode('utf-8')
t, ax, ay, az, gx, gy, gz = data.split(' ')
t = float(t)
ax = float(ax)
ay = float(ay)
az = float(az)
gx = float(gx)
gy = float(gy)
gz = float(gz)
except : # report error if we failed
print('Invalid data: ', data)
continue
# Store the new values into appropriate local arrays
times[i] = t / 1000 # convert time to seconds from ms
ax_values[i] = ax
ay_values[i] = ay
az_values[i] = az
gx_values[i] = gx
gy_values[i] = gy
gz_values[i] = gz
i += 1
ser.write(b'0')
sample_count += n_samples # Increment sample_count by amount of n_samples
return times, ax_values, ay_values, az_values, gx_values, gy_values, gz_values # Return two arrays (size n_sample) with new values
def update_plots(i):
global times, ax_values, ay_values, az_values, gx_values, gy_values, gz_values
global ax_flt, ay_flt, filter_coeffs, filter_ICs
global thresh, step
# shift samples left by 'NS'
times[:N-NS] = times[NS:]
# values[:N-NS] = values[NS:]
ax_values[:N-NS] = ax_values[NS:]
ay_values[:N-NS] = ay_values[NS:]
az_values[:N-NS] = az_values[NS:]
gx_values[:N-NS] = gx_values[NS:]
gy_values[:N-NS] = gy_values[NS:]
gz_values[:N-NS] = gz_values[NS:]
# grab new samples
times[N-NS:], ax_values[N-NS:], ay_values[N-NS:], az_values[N-NS:], gx_values[N-NS:], gy_values[N-NS:], gz_values[N-NS:] = grab_samples(NS)
ax_flt, filter_ICs = filterSignal(ax_values, filter_coeffs, filter_ICs, False)
ay_flt, filter_ICs = filterSignal(ay_values, filter_coeffs, filter_ICs, False)
# az_flt, filter_ICs = filterSignal(a_values, filter_coeffs, filter_ICs, False)
i = 0
for i in range(370):
if (max(ay_flt[i:i+30]) - min(ay_flt[i:i+30]) > thresh):
step = step + 1
i = i + 29
print (step)
# plot
[ax.set_xlim(times[0],times[N-1]) for ax in axes]
live_plots[0].set_data(times, ax_flt)
axes[0].set_ylim(min(ax_flt), max(ax_flt))
live_plots[1].set_data(times, ay_flt)
axes[1].set_ylim(min(ay_flt), max(ay_flt))
live_plots[2].set_data(times, az_values)
axes[2].set_ylim(min(az_values), max(az_values))
live_plots[3].set_data(times, gx_values)
axes[3].set_ylim(min(gx_values), max(gx_values))
live_plots[4].set_data(times, gy_values)
axes[4].set_ylim(min(gy_values), max(gy_values))
live_plots[5].set_data(times, gz_values)
axes[5].set_ylim(min(gz_values), max(gz_values))
return live_plots
if (__name__ == "__main__"):
# Open Serial
with serial.Serial(port=serial_port, baudrate=serial_baud, timeout=1) as ser:
try:
ser.flush()
sleep(2) # tell Arduino to start sending data. NOTE!!: you should change this for your setup
# initialize the figure
fig, axes = plt.subplots(2, 3)
times, ax_values, ay_values, az_values, gx_values, gy_values, gz_values = grab_samples(N)
# get HPF and LPF coefficients
filter_coeffs[3, :], filter_coeffs[2, :] = sig.butter(filter_order, LPF_cutoff, btype='lowpass', analog=False)
filter_coeffs[1, :], filter_coeffs[0, :] = sig.butter(filter_order, HPF_cutoff, btype='highpass', analog=False)
# get HPF and LPF initial conditions
filter_ICs[1, :] = sig.lfilter_zi( filter_coeffs[3, :], filter_coeffs[2, :] )
filter_ICs[0, :] = sig.lfilter_zi( filter_coeffs[1, :], filter_coeffs[0, :] )
ax_flt, filter_ICs = filterSignal(ax_values, filter_coeffs, filter_ICs, True)
ay_flt, filter_ICs = filterSignal(ay_values, filter_coeffs, filter_ICs, True)
live_plots = []
axes = axes.flatten()
live_plots.append(axes[0].plot(times, ax_flt, lw=2)[0])
live_plots.append(axes[1].plot(times, ay_flt, lw=2)[0])
live_plots.append(axes[2].plot(times, az_values, lw=2)[0])
live_plots.append(axes[3].plot(times, gx_values, lw=2)[0])
live_plots.append(axes[4].plot(times, gy_values, lw=2)[0])
live_plots.append(axes[5].plot(times, gz_values, lw=2)[0])
# initialize the y-axis limits and labels
# axes[0].set_ylim(0, 1200)
axes[0].set_title('AX')
# axes[0].set_xlabel('Time (s)')
axes[0].set_ylabel('Amplitude')
# axes[1].set_ylim(-0.0005, 0.0005)
axes[1].set_title('AY')
# axes[1].set_xlabel('Time (s)')
# axes[1].set_ylabel('Amplitude (V)')
# axes[2].set_ylim(0, 1200)
axes[2].set_title('AZ')
# axes[2].set_xlabel('Time (s)')
# axes[2].set_ylabel('Amplitude (V)')
# axes[3].set_ylim(0, 1200)
axes[3].set_title('GX')
# axes[3].set_xlabel('Time (s)')
axes[3].set_ylabel('Amplitude')
# axes[4].set_ylim(0, 1200)
axes[4].set_title('GY')
axes[4].set_xlabel('Time (s)')
# axes[4].set_ylabel('Amplitude (V)')
# axes[5].set_ylim(0, 1200)
axes[5].set_title('GZ')
# axes[5].set_xlabel('Time (s)')
# axes[5].set_ylabel('Amplitude (V)')
# set and start the animation and update at 1ms interval (if possible)
anim = animation.FuncAnimation(fig, update_plots, interval= 1)
plt.show()
finally:
ser.flush()
# sleep(2)
ser.write(b'0')
The code applies a LPF to the accelerometer's x-axis and y-axis values and then detects increases in the y-axis values within a 30 index window that exceed a pre-determined threshold value (thresh = 800). Whenever this is detected, the global step variable is incremented.
Using this Python script, we obtain a figure similar to the one below:
Figure of Filtered AX and AY Values for Step Counter
Note that only AX and AY values have been filtered. The other plots are left for comparison to demonstrate the effect that filtering has on our values from the Arduino. (Ideally, we want our signals to be noise-free for step detection purposes.)
In our first attempt to implement a step counter, we are successful in obtaining data from the Arduino to Python via Serial communication, and we are able to filter the data using the LPF. Using the filtered signal, we then increment the step variable every time the maximum value of a new set of sampled data exceeds the minimum value of the same set by a defined threshold value.
The issue with this strategy is that it only takes into account the maximum and minimum values of a set, regardless of their context. In other words, a "step" should consist of a peak in the signal, rather than a general spike. This means we will need a method of finding the peaks in our input signal, and using the detected peaks, we will increment the step count every time we detect a peak. In the next objective, I will be discussing how we will implement this functionality.
1.) We will improve our step counter algorithm in a few ways. First, we will need to center our data signal to the 0-axis, so a high-pass filter is needed in addition to the existing low-pass filter. Next, we will detect the peaks in our signal using scipy's find_peak() method. This method returns the indices of the input signal where a peak exists. The peak is detected within a set of given parameters, such as minimum height/width, in order to prevent random noise from being labeled as a peak. The following Python code is pasted below:
lab4_final_step_count.py
#Import libraries
import scipy.signal as sig
import serial
import sys
import requests
from time import sleep
from matplotlib import pyplot as plt
from matplotlib import animation
import numpy as np
# ==================== Globals ====================
N = 400 # number of samples plotted
NS = 20 # number of samples to read per iteration
sample_count = 0 # current sample count
filter_order = 3
freq = 200
LPF_cutoff = 10.0/100
HPF_cutoff = 0.5/100
thresh = 800
step = 0
first = False
times, ay_values, processed = np.zeros((3, N)) # data vectors
ay_flt, loc = np.zeros((2,N))
filter_ICs = np.zeros((2, filter_order))
filter_coeffs = np.zeros((4, filter_order + 1))
connection_secure = False
# These variables are just defined commands for the BLE module
set_clear = "AT+CLEAR"
set_renew = "AT+RENEW"
set_reset = "AT+RESET"
set_role= "AT+ROLE1"
set_imme = "AT+IMME1"
set_adty = "AT+ADTY3"
set_AT = "AT"
# Make sure to change the MAC Address to your specific HM-10 Module's MAC Address
set_connect = "AT+CON3403DE349771"
BLE_port = '/dev/cu.usbserial' # the BLE port to use
# Initializes the HM-10 Module to our preferred settings for this specific script
def initialize_BLE (ser):
ser.reset_input_buffer()
ser.write(set_AT.encode())
sleep(0.5)
print("> AT: " + read_BLE(ser))
ser.write(set_clear.encode())
sleep(0.5)
print("> CLEAR: " + read_BLE(ser))
ser.write(set_reset.encode())
sleep(0.5)
print("> RESET: " + read_BLE(ser))
sleep(2)
ser.write(set_renew.encode())
sleep(0.5)
print("> RENEW: " + read_BLE(ser))
sleep(2)
ser.write(set_reset.encode())
sleep(0.5)
print("> RESET: " + read_BLE(ser))
sleep(2)
ser.write(set_imme.encode())
sleep(0.5)
print("> IMME: " + read_BLE(ser))
ser.write(set_adty.encode())
sleep(0.5)
print("> ADTY: " + read_BLE(ser))
ser.write(set_role.encode())
sleep(0.5)
print("> ROLE: " + read_BLE(ser))
sleep(2)
ser.write(set_reset.encode())
sleep(0.5)
print("> RESET: " + read_BLE(ser))
sleep(2)
# Reads Bluetooth input buffer and stores it into a string
def read_BLE (ser):
msg = ""
if (ser.in_waiting > 0):
msg = ser.readline(ser.in_waiting).decode('utf-8')
return msg
def filterSignal(data, filter_coeffs, filter_ICs):
global first
# ==================== Step 1 ====================
if (first == True) :
LPF, zi_new_LPF = sig.lfilter(filter_coeffs[3, :], filter_coeffs[2, :], data, zi = filter_ICs[1, :]*data[0])
filter_ICs[1, :] = zi_new_LPF[:]
else :
LPF, zi_new_LPF = sig.lfilter(filter_coeffs[3, :], filter_coeffs[2, :], data, zi = filter_ICs[1, :])
filter_ICs[1, :] = zi_new_LPF[:]
# ==================== Step 2 ====================
# Apply HPF filter
if (first == True) :
proc_data, zi_new_HPF = sig.lfilter(filter_coeffs[1, :], filter_coeffs[0, :], LPF, zi = filter_ICs[0, :]*LPF[0])
filter_ICs[0, :] = zi_new_HPF[:]
else :
proc_data, zi_new_HPF = sig.lfilter(filter_coeffs[1, :], filter_coeffs[0, :], LPF, zi = filter_ICs[0, :])
filter_ICs[0, :] = zi_new_HPF[:]
return proc_data, filter_ICs
def grab_samples(n_samples) :
global sample_count, filter_coeffs, filter_ICs, first
# Create local arrays of size(n_samples) to store new data
times, ay_values = np.zeros((2,n_samples))
ay_flt = np.zeros(n_samples)
print("Grabbing Samples...")
i = 0
while i < n_samples : # while we still need samples, keep reading
try:
ser.write(b'1')
data = ser.readline().strip().decode('utf-8')
t, ay = data.split(' ')
t = float(t)
ay = float(ay)
except : # report error if we failed
print('Invalid data: ', data)
continue
# Store the new values into appropriate local arrays
times[i] = t / 1000 # convert time to seconds from ms
ay_values[i] = ay
i += 1
ser.reset_input_buffer()
ser.write(b'0')
ay_flt, filter_ICs = filterSignal(ay_values, filter_coeffs, filter_ICs)
print("Done Grabbing Samples!")
sample_count += n_samples # Increment sample_count by amount of n_samples
return times, ay_values, ay_flt # Return two arrays (size n_sample) with new values
def update_plots(i):
global times, ay_values
global ay_flt, filter_coeffs, filter_ICs, loc
global thresh, step
print("Updating Plots...")
# shift samples left by 'NS'
times[:N-NS] = times[NS:]
ay_values[:N-NS] = ay_values[NS:]
ay_flt[:N-NS] = ay_flt[NS:]
loc[:N-NS] = loc[NS:]
# grab new samples
times[N-NS:], ay_values[N-NS:], ay_flt[N-NS:] = grab_samples(NS)
loc[N-NS:] = 0
print("Detecting Peaks")
peaks, props = sig.find_peaks(ay_flt[N-NS:], height=500, distance=10, width=(None,20))
ser.reset_input_buffer()
print("Incrementing Steps...")
for i in peaks:
loc[N-NS+i] = ay_flt[N-NS+i]
step = step + 1
ser.write(b'#')
print(step)
# plot
[ax.set_xlim(times[0],times[N-1]) for ax in axes]
live_plots[0].set_data(times, ay_values)
axes[0].set_ylim(min(ay_values), max(ay_values))
live_plots[1].set_data(times, ay_flt)
axes[1].set_ylim(min(ay_flt), max(ay_flt))
live_plots[2].set_data(times, loc)
axes[2].set_ylim(min(loc), max(loc))
return live_plots
if (__name__ == "__main__"):
# Establish connection with HM-10 BLE Module
with serial.Serial(port = BLE_port, baudrate = 9600, timeout = 1) as ser:
try:
sleep(2)
# Initialize the BLE module to our preferred settings
initialize_BLE(ser)
# Establishes a connection with the BLE module
ser.write(set_connect.encode())
sleep(0.5)
ble_in = read_BLE(ser)
print("> CONNECT: " + ble_in)
print("> Attempting to Connect...")
sleep(15)
ble_in = read_BLE(ser)
if ( ble_in.find('OK+CONNF') > -1 ):
print("Connection Failed")
connection_secure = False
else:
print("Connection Successful")
connection_secure = True
ble_in = ""
# This loop will execute until a connection has been reestablished
while (connection_secure == False):
ser.reset_input_buffer()
command = input("> 1.) Hit ENTER to try Reconnecting\n> 2.) Enter QUIT to exit:\t")
if (command == ""):
ser.write(set_connect.encode())
sleep(0.5)
ble_in = read_BLE(ser)
print("> " + ble_in)
print("> Attempting to Reconnect...")
sleep(15)
ble_in = read_BLE(ser)
print("> " + ble_in)
if ( ble_in.find('OK+CONNF') > -1 ):
print("> Reconnection Unsuccessful")
else:
print("> Reconnection Successful")
connection_secure = True
if (command == "QUIT" or command == "quit"):
sys.exit(0)
print("Get HPF LPF Coefficients...")
# get HPF and LPF coefficients
filter_coeffs[3, :], filter_coeffs[2, :] = sig.butter(filter_order, LPF_cutoff, btype='lowpass', analog=False)
filter_coeffs[1, :], filter_coeffs[0, :] = sig.butter(filter_order, HPF_cutoff, btype='highpass', analog=False)
print("Get HPF LP initial conditions....")
# get HPF and LPF initial conditions
filter_ICs[1, :] = sig.lfilter_zi( filter_coeffs[3, :], filter_coeffs[2, :] )
filter_ICs[0, :] = sig.lfilter_zi( filter_coeffs[1, :], filter_coeffs[0, :] )
print("Initializing Figure...")
# initialize the figure
fig, axes = plt.subplots(3,1)
first = True
times, ay_values, ay_flt = grab_samples(N)
first = False
peaks, props = sig.find_peaks(ay_flt, height=500, distance=10, width=(None,10))
# peaks, props = sig.find_peaks(ay_flt[N-NS:], width=(None,20), threshold=300)
for i in peaks:
loc[i] = ay_flt[i]
step = step + 1
ser.write(b'#')
live_plots = []
live_plots.append(axes[0].plot(times, ay_values, lw=2)[0])
live_plots.append(axes[1].plot(times, ay_flt, lw=2)[0])
live_plots.append(axes[2].plot(times, loc, lw=2)[0])
# initialize the y-axis limits and labels
axes[0].set_ylim(min(ay_values), max(ay_values))
axes[0].set_title('Peak Detection Figures')
axes[0].set_ylabel('AY Raw Signal')
axes[1].set_ylabel('AY Filtered')
axes[1].set_ylim(min(ay_flt), max(ay_flt))
axes[2].set_xlabel('Time (s)')
axes[2].set_ylabel('Peaks')
axes[2].set_ylim(min(loc), max(loc))
# set and start the animation and update at 1ms interval (if possible)
anim = animation.FuncAnimation(fig, update_plots, interval= 1)
plt.show()
finally:
ser.flush()
As usual, make sure to replace the serial_port to be the port that your BLE module is connected to. The Python script obtains the IMU values from Arduino via Serial, filters the signals, detects the peaks that are considered "steps", and then prints out a '#' character back to the Arduino whenever a step is detected.
2.) From the Arduino side of things, we implement the following code so that the Arduino takes in a command character from Python via BLE, and based on the command, it will either print data or increment the stepCount (if a '#' character is sent) to display on the OLED screen.
Lab4FinalStepCount.ino
#include "I2Cdev.h"
//#include "MPU6050.h"
#include "MPU6050_6Axis_MotionApps20.h"
#include "Wire.h"
#include "TimerOne.h"
#include <SPI.h>
#include <Adafruit_GFX.h>
#include <Adafruit_SSD1306.h>
#include <SoftwareSerial.h>
SoftwareSerial BTserial(8,9);
#define OLED_RESET 4
Adafruit_SSD1306 display(OLED_RESET);
#define NUMFLAKES 10
#define XPOS 0
#define YPOS 1
#define DELTAY 2
#if (SSD1306_LCDHEIGHT != 32)
#error("Height incorrect, please fix Adafruit_SSD1306.h!");
#endif
const int MPU_addr=0x68; // I2C address of the MPU-6050, can be changed to 0x69
MPU6050 IMU(MPU_addr);
const int interruptPin = 2;
volatile bool ipinReady = false;
int buttonState = 0;
int lastState = 0;
const int buttonPin = 4;
volatile int16_t ax, ay, az, tp, gx, gy, gz;
int samplePeriod = 10000; // the target sample period, 5000 microseconds, 200Hz
unsigned long startTime = 0;
unsigned long volatile currentTime = 0;
bool newRead = false;
int sampleSize = 0;
bool awake = true;
bool volatile sending = false;
int stepCount = 0;
const int timebuf_size = 10;
char timebuf[timebuf_size];
char buf[timebuf_size];
// Function to print step count to OLED
void printOLED() {
display.clearDisplay();
display.setCursor(0,0);
display.print("Step Count: ");
display.println(stepCount);
display.display();
}
// Timer interrupt to sample our data
void printDataISR() {
// Saves the time
currentTime = (micros()-startTime)/1e3;
newRead = true;
}
// Check the interrupt pin to see if there is data available in the MPU's buffer
void interruptPinISR() {
ipinReady = true;
}
void setup(){
// initalize BTserial
BTserial.begin(9600);
delay(2000);
BTserial.write("AT");
delay(200);
BTserial.write("AT+CLEAR");
delay(200);
BTserial.write("AT+RENEW");
delay(200);
BTserial.write("AT+RESET");
delay(500);
// initialize the pushbutton pin as input:
pinMode(buttonPin, INPUT_PULLUP);
buttonState = digitalRead(buttonPin);
// initialize with the I2C addr 0x3C (for the 128x32)
display.begin(SSD1306_SWITCHCAPVCC, 0x3C);
display.clearDisplay();
display.setTextSize(1);
display.setTextColor(WHITE);
display.setCursor(0,0);
printOLED();
// Intialize the IMU and the DMP on the IMU
IMU.initialize();
IMU.dmpInitialize();
IMU.setDMPEnabled(true);
Wire.begin();
Wire.beginTransmission(MPU_addr);
Wire.write(MPU_addr); // PWR_MGMT_1 register
Wire.write(0); // set to zero (wakes up the MPU-6050)
Wire.endTransmission(true);
// Set the timer interrupt
startTime = micros();
Timer1.initialize(samplePeriod);
Timer1.attachInterrupt(printDataISR);
// Create an interrupt for pin2, which is connected to the INT pin of the MPU6050
pinMode(interruptPin, INPUT);
attachInterrupt(digitalPinToInterrupt(interruptPin), interruptPinISR, RISING);
}
void loop(){
// Reads ButtonPin
lastState = buttonState;
buttonState = digitalRead(buttonPin);
if ((buttonState == 0) && (lastState == 1) )
{
awake = !awake;
sending = false;
}
// IMPORTANT: Send a '1' via Serial in order to begin reading the data, and a '0' to stop sending data
if (BTserial.available()) {
int incomingByte = BTserial.read();
if (incomingByte == '#') {
stepCount++;
printOLED();
}
if (incomingByte == '1')
sending = true;
if (incomingByte == '0')
sending = false;
}
if (newRead && sending && ipinReady && awake) {
Wire.beginTransmission(MPU_addr);
Wire.write(0x3B); // starting with register 0x3B (ACCEL_XOUT_H)
Wire.endTransmission(false);
Wire.requestFrom(MPU_addr,14,true); // request a total of 14 registers
// Bit Wire.read() gets 8 bits. Every axis reading is 16 bits. Threfore we read two consecutive bytes
// from Wire, and store them next to each other as [ byte1 byte2 ]. This can be done by left bit shifting
// byte1 and ORing it with byte2
//Accelerometer (3 Axis)
ax=Wire.read()<<8|Wire.read(); // 0x3B (ACCEL_XOUT_H) & 0x3C (ACCEL_XOUT_L)
ay=Wire.read()<<8|Wire.read(); // 0x3D (ACCEL_YOUT_H) & 0x3E (ACCEL_YOUT_L)
az=Wire.read()<<8|Wire.read(); // 0x3F (ACCEL_ZOUT_H) & 0x40 (ACCEL_ZOUT_L)
//Temperature
tp=Wire.read()<<8|Wire.read(); // 0x41 (TEMP_OUT_H) & 0x42 (TEMP_OUT_L)
//Gyroscope (3 Axis)
gx=Wire.read()<<8|Wire.read(); // 0x43 (GYRO_XOUT_H) & 0x44 (GYRO_XOUT_L)
gy=Wire.read()<<8|Wire.read(); // 0x45 (GYRO_YOUT_H) & 0x46 (GYRO_YOUT_L)
gz=Wire.read()<<8|Wire.read(); // 0x47 (GYRO_ZOUT_H) & 0x48 (GYRO_ZOUT_L)
// Parse data as a buffer to send via BLE
// Displays the raw value and the time it was sampled
dtostrf(currentTime, 9, 0, timebuf);
BTserial.print(timebuf);
BTserial.print(' ');
dtostrf(ay, 4, 2, buf);
BTserial.println(ay);
newRead = false;
ipinReady = false;
sending = false;
}
}
3.) Upload the Arduino sketch and run the Python script. You should see a live plot of the raw signal, filtered signal, and the peak detection. An example can be seen below:
Output Live Plot of Final Step Counter
The top graph plots the raw signal received from Arduino via BLE. The middle graph is the signal after being put through a low-pass and high-pass filter. Lastly, the bottom plot is a visualization of detected peaks from our step counting algorithm. Whenever a peak is detected, the bottom plot displays a spike based on the amplitude of the signal at the instance of the detected peak.
Congratulations! You have finished Lab 4. To ensure your system works properly, refer below to a video displaying the final system at work:
Progress of Events: