是个啥坑？

+ K3 |' C& D- Y  m7 b9 ^一个项目中用到一个传感器测量一物理量，这里假定测量温度吧。需要判断其变化趋势，利用这个变化趋势去做一些应用。
8 h9 U) n" w( f5 X/ D5 L那么要怎么判断一个物理量的变化趋势呢？我们能自然能想到去求取该随机序列的变化率。这里涉及到一些数序定义。随机序列有很多可能的来源，最为常见是我之前在<<模数采样知多少>>中介绍的模数采样。
这样将S(t)信号转换为离散信号序列S(n)，那么对于当前时刻其斜率怎么求取呢？(这里忽略中间的过度态，仅将其看为线段相连，当然现实应用中如果有更高要求，可以做曲线拟合)


但是如果只判断，斜率极容易误判，比如下面这样的情况：

其斜率一会儿正，一会儿负，但是其总体趋势又是在增加的，所以只考察斜率显然不可取，获取需要在代码在加各种复杂的条件或者限值去判断。即使加这么多条件系统仍然可能表现的非常不健壮。

对于模拟信号2而言，趋势又在不断变化。那么怎么做才能稳定呢？先卖个关子？
- O0 T4 ^; {: N0 H& @
函数的凹凸性

凹函数

凹函数是一个定义在某个向量空间的凸集C（区间）上的实值函数f。设f为定义在区间I上的函数，若对I上的任意两点x1<x2和任意的实数t属于(0,1),总有,
# s# }* f3 y5 g( \# \/ ^' a, u
, p# ]) l2 e& H& t$ U4 i则称函数f为l上凹函数，有的书上也称为下凸函数。
; S+ l5 o/ u1 Z  x6 y  ~" h
4 Z! E: ^8 K5 I5 O如果把上述条件中的“≥”改成“>”，则叫做严格上凹函数，或叫做严格下凸函数。
上面是一维函数情况，这里来个2维函数的图，方便理解。
& n* u& r. ]) b9 r1 C; w0 O& Q' ?
0 J- K1 b. z' m, ~& X+ l6 c  D0 G- k: g


# b( R2 h% n& J8 [
- G# _+ [; v& u. I0 Q) B& b2 U) D& }/ J凸函数

设f为定义在区间I上的函数，若对I上的任意两点x1<x2和任意的实数t属于(0,1)，上面不等式变成大于等于，则在该区间为凸函数。

. A1 u9 a! e# L# a# l# j, z5 P可见，凹凸是相对的，如f(x)在某区间为凹，则-f(x)则在该区间为凸。
. S; E- E. P; b. t9 N/ ~, s性质

若一个函数在某区间二阶可导且大于0，则函数在该区间为凹函数
若一个函数在某区间二阶可导且小于0，则函数在该区间为凸函数
/ A5 N. @. }: r/ o/ [1 C2 e证明，这里就不推导了，可以利用拉格朗日中值定理可以推导出上面这个性质。
来看一下会动的图,加深一下理解:
6 e+ \0 U1 y3 k7 I
函数 切线为蓝色，曲线向上凹,绿色表示曲线是向下凹的，红色表示曲线的拐点。
sin(2x)的一阶导数为：
* R+ c& C. Q8 E0 Qsin(2x)的二阶导数为：
& E( U- E  W. I8 W
装逼结束，也可能没装对～～～

回到坑里

/ p6 A5 r* D# y8 I5 [% z3 K通过上面装逼，是否可以利用离散序列的求导数来判断传感器的变化趋势。啥？导数？又要开始表演了？

/ V" C0 V) d. T% ~0 b前面说了一阶导数是这样的：
那么二阶导数是哪样捏？

化简一下：
/ e9 m0 B9 K8 K" P- `
其中S[n]表示当前测量点，S[n-1]表示前一个测量点，S[n-2]表示前第2个测量点。
上代码
- A# Y0 C8 b# n# Z& R# R0 X- X
/ O6 D& O  O4 C7 i6 o" \7 ^4 @8 }#include <stdio.h>
% A  q9 B0 O1 Z. I! t% a8 P0 M#include <math.h>
#include <string.h>
typedef struct _T_2ND_DRV
+ @' b8 o$ i+ O8 L{
9 }" p: ^! k/ P. _5 T. G# r/ s    float xn1;
    float xn2;
}t_2ND_DRV;
* \7 D1 F) p* R; Utypedef struct _T_1ST_DRV
2 E) r) h& ]0 v$ S{
# _* A2 F; {4 h% [  x    float xn1;
0 O8 H3 O  g8 D/ V}t_1ST_DRV;
. N# x( p5 X$ N& w4 o9 K$ y2 u4 K
" U: v: _; O9 S% S7 j# S% f) cvoid init_second_derivative(t_2ND_DRV *pSndDrv)
{
    pSndDrv->xn1 = 0;
, S# k" s$ ^% O  {$ l6 F0 x/ U    pSndDrv->xn2 = 0;
, {) E/ K8 F+ q}
; l# h! z2 h" V# |
  W3 G9 D8 r: V! w7 Q& vfloat second_derivative(t_2ND_DRV *pSndDrv, float xn,float T)
" Y/ D& ?& j& E0 O. e- M, r{
     float result=0.0f;
! X& {* S5 h8 U     if(T<=0)
# }2 y. X) D& L         return 0x7FBFFFFF; /*非法数据*/
' ]. m8 H+ S- V     result = (xn-2*pSndDrv->xn1+pSndDrv->xn2)/T/T;
8 v1 m8 s: ~/ S+ k     pSndDrv->xn2 = pSndDrv->xn1;
7 K6 j2 f- F: V- V  B     pSndDrv->xn1 = xn;
   
+ ~2 E& a6 |4 t' Z- I     return result;
/ H4 ]9 V+ ^  w! O* P}

void init_fisrt_derivative(t_1ST_DRV *p1stDrv)
( r2 \/ I% L7 I9 }- ?{
  Q: M2 b7 R& m) L    p1stDrv->xn1 = 0;
}
* Y) {2 o$ w( Z% p. P5 |7 w  R
6 c4 ~  g: I: D2 g) [! Nfloat fisrt_derivative(t_1ST_DRV *p1stDrv, float xn,float T)
{
     float result=0.0f;
! _, t% F) Z6 o6 M     if(T<=0)
         return 0x7FBFFFFF; /*非法数据*/
, _) e& X7 X$ v( C, F     result = (xn-p1stDrv->xn1)/T;
     p1stDrv->xn1 = xn;
   
3 ?3 l- a) E! y! K     return result;
. k6 S4 X) I3 b8 ?8 w% X}
# S! L- @) E5 t8 p- g#define PI 3.1415f
#define SAMPLE_RATE 500.0f
% G0 ]7 J! J1 V$ w5 n4 J5 C/ b#define SAMPLE_T (1/SAMPLE_RATE)
" F2 X. T$ S1 U#define SAMPLE_SIZE (100)
) D! ^2 G& ]6 R: r; b: g+ K) dint main()
{
" r- T: d: g) |2 z" {) |, j    float sim1[SAMPLE_SIZE];
    float sim2[SAMPLE_SIZE];
    float out1[SAMPLE_SIZE];
    float out2[SAMPLE_SIZE];
    t_2ND_DRV sndDrv;
5 c" ?* H$ |/ A5 Z8 A9 F  y    t_1ST_DRV frtDrv;
; s. S- u8 ?+ L5 s' f2 |1 e    init_fisrt_derivative(&frtDrv);
    init_second_derivative(&sndDrv);
8 y! t4 Q+ r( T) z1 f. h1 Z* x6 S   
    FILE *pFile=fopen("./simulationSin.csv","wt+");
    if(pFile==NULL)
0 g) U5 @/ X- k0 q( Z2 a    {
        printf("simulationSin.csv opened failed");
$ l( G: {2 b. r2 Y5 ?/ y8 k        return -1;
    }
   
    for(int i=0;i<SAMPLE_SIZE;i++)
! E% W* [6 s6 O* D0 K    {
        sim1=10*sin(2*PI*10*i/500);
    }
0 \# P/ I9 o1 Y0 {& R: ^' y    for(int i=0;i<SAMPLE_SIZE;i++)
- j* A. k: g% E% M" f1 O    {
; C& i: {3 v8 p6 |, y* {) d9 m        out1=fisrt_derivative(&frtDrv,sim1,SAMPLE_T);
        out2=second_derivative(&sndDrv,sim1,SAMPLE_T);
        fprintf(pFile,"%f,%f,%f\n",sim1,out1,out2);
    }

  p1 f4 a9 O! J/ t5 f# o+ w    fclose(pFile);
9 H' P# N& w5 o' U8 \
0 |2 c2 Z; G$ p2 l; ~- c" M    return 0;
}
* H0 i) W, G7 F7 n% k7 G5 Q: s1 n: f$ ?忽略前两个点，利用excel生成曲线：
  t5 N: {+ y) ?/ Q) W# s5 Z
从图中可看出：
一阶导数为正时，函数递增趋势；
一阶导数为负时，函数递减趋势；
二阶导数为0时，出现拐点，趋势改变；此时如果左右两侧的一阶导符号相反，则出现极值。
& I# Q) s; O7 o, }二阶导数为负时，其一阶导数也即原函数斜率规律单调减，二阶导数为正时，其一阶导数也即原函数斜率规律单调增。
" Q& W" }* i2 y( w/ n" s( n再进一步：
一阶导数与二阶导数结合起来看，就可以看出测量值变化趋势的趋势，比如在前1/4周期，此区间变换趋势为增，也即一阶导数为正，而其二阶导数为负，也可以看出递增的趋势是逐渐减小到0的。
代码优化
! E. R: O" b6 ^. {( ?# {, |( l/ _
& m1 Z# Y' [! r# L5 D: A如果只是做定性判断，上述函数，完全没必要与采样周期做除法，只需要考察其增量即可，代码可优化如下：
: w4 ~, _; f6 O# Qtypedef struct _T_2ND_DRV
{
+ \' O% t' o) W0 Y% Q. N    float xn1;
    float xn2;
* G+ b" v7 h/ T& f( B}t_2ND_DRV;
$ ]1 x+ S* h/ I$ U" [typedef struct _T_1ST_DRV
{
6 S) i7 z& t% g  N    float xn1;
  I7 E7 E9 n6 [$ j, P$ I. j) G2 W}t_1ST_DRV;

* o6 x. a  d- \, }4 X' |5 v2 \void init_second_derivative(t_2ND_DRV *pSndDrv)
{
- A: j( E6 S$ C    pSndDrv->xn1 = 0;
    pSndDrv->xn2 = 0;
% D  j1 Q$ X) `/ N4 w0 P# H}

4 [1 }0 d: {7 V" m: Y' j9 kfloat second_derivative(t_2ND_DRV *pSndDrv, float xn)
! k& b" K4 A% U2 D" R5 g1 y7 n* ~{
     float result=0.0f;
9 m% O9 |: Q9 H$ l     result = xn-2*pSndDrv->xn1+pSndDrv->xn2;
     pSndDrv->xn2 = pSndDrv->xn1;
     pSndDrv->xn1 = xn;
3 p1 F9 m% V% S  j   
     return result;
}
* {3 T" _6 ?6 j8 u: E
void init_fisrt_derivative(t_1ST_DRV *p1stDrv)
3 F* x( J, ]# W' B{
8 l9 ]* O% h3 u5 I1 x% h" Z    p1stDrv->xn1 = 0;
& \& n0 u- L) p! U. @6 ~& g}

$ g! f8 _; E! B: A3 Qfloat fisrt_derivative(t_1ST_DRV *p1stDrv, float xn)
. X% u1 p4 `' [# W{
     float result=0.0f;
     result = xn-p1stDrv->xn1;
: t% F% Q( M9 G3 D% P6 R     p1stDrv->xn1 = xn;
8 W* X: U0 R7 L' [0 j; F   
     return result;
( i& I* m. h6 t* W" |}
0 @" M7 U4 q: S" r意外收获

% |2 l. r% b8 L* _2 p6 v4 A这里意外引入一个可能很多人没注意的知识点NaN，在计算中，NaN代表非数字，是数字数据类型的成员，可以将其解释为不确定的或无法表示的值，尤其是在浮点运算中。1985年，IEEE 754浮点标准引入了NaN的系统使用，并表示了其他无限量（如无穷大）。
前述函数返回0x7FBFFFFF,也就是表示无穷大。
* q$ Q9 U8 e9 q; e不同的操作系统和编程语言可能具有NaN的不同字符串表示形式：
nan
NaN
NaN%
NAN
( [+ e; A, d+ X1 Q% K NaNQ
NaNS
qNaN
+ D* Q0 d9 K9 D1 {  P$ _ sNaN
1.#SNAN
1.#QNAN
5 w# a$ R; R. c9 g1 u2 n, H -1.#IND
实际上，由于编码的NaN具有符号，因此通常也可以在NaN的字符串表示中找到它们，例如：
-NaN
; F3 {9 J; M: G7 u  NaN12345
3 L4 z' ^+ o- a: }2 _ -sNaN12300
: v. d- ?: Z1 M% U -NaN(s1234)
5 X1 U/ F. `0 r6 e9 o工程应用

这里给出我的建议方案：

将传感器信号经由电路处理，模数采样，在进入前级数字滤波器，滤除不必要的噪声，在进行一阶/二阶求导。对于一阶和二阶求导再做一级移动平均滤波，最后在按照上面描述进行判别变化趋势，基本比较健壮了。实际移动均值滤波长度不宜选择过长，否则响应就比较滞后了。不能对传感器的变化趋势做出实时的判别。加了后级均值滤波器，则会消除由于波形忽上忽下的随机噪声干扰影响，使得系统判别更为健壮，实际滤波器长度需根据不同的场合进行调试优化。或者也可以选择别的IIR/FIR滤波器形式实现。
7 ]6 P- T1 R9 Q6 _7 z总结一下
! q  O3 h4 Q% Y6 a% X  G. U$ g1 x
做为嵌入式er编程，有时候有必要去看看数学书，了解一下数学原理的背后故事，可能会给你带来意想不到的作用哦。

求导怎么计算？太高深了

高手啊，在测试校准中要用到这些知识吧？
期望能够分享更多的知识！

这是在干啥？

太高深 了