背景
有些产品可能需要用到旋钮来调节音量,比如对讲机,收音机等。如何实现呢?常见的一种实现方案是使用电位器接一个ADC输入,通过ADC检测电压的大小变化来调节音量大小,软件实现上需要定时检测ADC的电压变化,检测周期长短影响音量调节体验,设置的周期长,反应速度慢,会出现突然声音大或者突然声音小的情况,检测时间太频繁,系统又会频繁唤醒,无法进入低功耗状态。
另一种方式是是使用旋转编码器接几个GPIO,通过旋转编码器可以检测到正转反转,通过其输出信号可以得到旋转角度。
电位器实物图:
旋转编码器实物图:
在某项目中使用ADC+可调电位器的方式实现音量调节功能,电路如下图:
在使用ADC时发现一个问题,从ADC读取到的电压值并不是实际电压,从ADC读取的电压值范围是0~184mv,测量得到的实际电压为0~1668mv,所以需要对ADC的原始数据做处理,通过软件拟合得到真实电压。遇到问题,解决问题,学习使用python处理曲线拟合的几种方法。
python 调试环境介绍
对我的python开发调试环境做个简单介绍,操作系统是win 10,安装的python 是3.8,同时安装了 jupyter lab。
平时一些简单的python脚本编写、算法验证通常都是在 jupyter lab环境上编码调试。
为什么使用jupyter lab呢?
python对我来说是一门辅助性的工具,通常编写的代码量不大,jupyter lab 上可以运行代码,绘图显示方便,同时编写笔记,所以jupyter lab 非常适合我的使用需求,而不需要重量级的IDE开发环境。
下面是一个来自jupyter lab官网的运行界面图,除了支持python,还支持c++、R、Julia等其他编程语言。
python 线性拟合方法
线性拟合需要用到的python包有 numpy、matplotlib 和 scipy,通过pip安装即可,然后导入需要用到的python包
1
2
3
4
5
6
|
import numpy as np
import matplotlib as mpl
import matplotlib.pyplot as plt
from numpy import polyfit, poly1d
from scipy.linalg import lstsq
from scipy.stats import linregress
|
导入从adc读取到的原始数据和实际测量得到的电量数据,数据量不大,定义为列表,然后绘制ADC原始数据和实际电压的曲线图。
1
2
3
4
5
6
7
8
|
# ADC 采集的原始数据
adc = [ 192, 333, 372, 399, 421, 435, 450, 453, 459, 462, 456, 436, 419, 386, 340, 199, 184]
# 实际测量的电压数据
vo = [702, 1212, 1343, 1442, 1525, 1572, 1623, 1638, 1652, 1665, 1647, 1580, 1512, 1395, 1232, 729, 679]
# 绘图
p = plt.plot(adc, vo, 'rx')
|
从上图看,原始数据和实际电压数据基本上保持线性关系,这是比较幸运的情况,从图上看应该是一条直线,通过直线拟合的方法来得到实际电压没有问题。
怎么获取这条直线的参数呢?
方法1,numpy 提供的多项式拟合函数polyfit计算。
线性拟合为一阶多项式,一阶多项式$$y = a_1 x + a_0$$拟合成功返回两个系数 $$[a_1, a_0]$$。
1
2
3
|
coeff1 = polyfit(adc, vo, 1)
print(coeff1)
# [ 3.56483415 20.43063565]
|
使用poly1d生成多项式函数。
1
2
3
|
f = poly1d(coeff1)
print(f)
# 3.565 x + 20.43
|
然后画出拟合曲线图。
1
2
|
p = plt.plot(adc, f(adc))
p = plt.plot(adc, vo, 'rx')
|
方法2,使用scipy提供的线性回归函数linregress求解
1
2
3
4
|
# 求解系数
slope, intercept, r_value, p_value, stderr = linregress(adc, vo)
print(slope, intercept)
# 3.5648341453932266 20.430635650877775
|
绘制拟合曲线图
1
2
3
4
5
6
|
f = poly1d(coeff)
print(f)
# 3.565 x + 20.43
p = plt.plot(adc, vo, 'rx')
p = plt.plot(adc, f(adc))
|
方法3,求最小二乘解
当使用一个 N-1 阶的多项式拟合这 M 个点时,有这样的关系存在$$XC = Y$$ ,
$$\left[ \begin{matrix} x_0^{N-1} & \dots & x_0 & 1 \\ x_1^{N-1} & \dots & x_1 & 1 \\ \dots & \dots & \dots & \dots \\ x_M^{N-1} & \dots & x_M & 1 \end{matrix}\right] \left[ \begin{matrix} C_{N-1} \\ \dots \\ C_1 \\ C_0 \end{matrix} \right] = \left[ \begin{matrix} y_0 \\ y_1 \\ \dots \\ y_M \end{matrix} \right]$$
使用一阶多项式求解,即N=2,先将adc扩展为X,vo转换为Y:
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
|
x = np.array(adc)
Y= np.array(vo).reshape(-1,1)
X = np.hstack((x[:,np.newaxis], np.ones((x.shape[-1],1))))
print(X[0:8])
print(Y[0:8])
'''
[[192. 1.]
[333. 1.]
[372. 1.]
[399. 1.]
[421. 1.]
[435. 1.]
[450. 1.]
[453. 1.]]
[[ 702]
[1212]
[1343]
[1442]
[1525]
[1572]
[1623]
[1638]]
'''
|
求解:
1
2
3
4
5
6
7
8
9
10
|
C, resid, rank, s = lstsq(X, Y)
C, resid, rank, s
"""
(array([[ 3.56483415],
[20.43063565]]),
array([137.64565162]),
2,
array([1.59853058e+03, 9.95236176e-01]))
"""
|
绘图:
1
2
3
4
5
6
|
f = poly1d(C.reshape(1,-1)[0])
print(f)
# 3.565 x + 20.43
p = plt.plot(adc, vo, 'rx')
p = plt.plot(adc, f(adc))
|
非线性拟合
-
多项式拟合正弦函数
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
|
# 定义sin函数
x = np.linspace(-3 * np.pi,3 * np.pi,100)
y = np.sin(x)
# 多项式拟合,从一阶多项式到9阶多项式
y1 = poly1d(polyfit(x,y,1))
y3 = poly1d(polyfit(x,y,3))
y5 = poly1d(polyfit(x,y,5))
y7 = poly1d(polyfit(x,y,7))
y9 = poly1d(polyfit(x,y,9))
# 绘图
p = plt.plot(x, y, 'k-')
p = plt.plot(x, y1(x))
p = plt.plot(x, y3(x))
p = plt.plot(x, y5(x))
p = plt.plot(x, y7(x))
p = plt.plot(x, y9(x),"g--")
a = plt.axis([-3 * np.pi, 3 * np.pi, -1.25, 1.25])
|
-
拟合自定义曲线
$$y = a e^{-b sin( f x + \phi)}$$
定义非线性函数:
1
2
3
4
|
def function(x, a , b, f, phi):
"""a function of x with four parameters"""
result = a * np.exp(-b * np.sin(f * x + phi))
return result
|
1
2
3
4
5
6
7
8
9
10
11
|
from scipy.optimize import curve_fit
from scipy.stats import norm
# 绘制原始曲线和加入噪声后曲线
x = np.linspace(0, 2 * np.pi, 50)
actual_parameters = [3, 2, 1.25, np.pi / 4]
y = function(x, *actual_parameters)
y_noisy = y+ 1.8 * norm.rvs(size=len(x))
p = plt.plot(x, y, 'k-')
p = plt.plot(x, y_noisy, 'gx')
|
使用curve_fit求解
1
2
3
4
5
6
7
|
p_est, err_est = curve_fit(function, x, y_noisy)
# 第一个返回的是函数的参数,第二个返回值为各个参数的协方差矩阵
# 绘制结果曲线
p = plt.plot(x, y_noisy, "gx")
p = plt.plot(x, y, 'k')
p = plt.plot(x, function(x, *p_est), "r")
|
参考资料
- 04.04 curve fitting (http://lijin-thu.github.io/04.%20scipy/04.04%20curve%20fitting.html)