代码考古:PID

某比赛中外置的计算机里,提供了示例代码。部分代码藏了起来。一些“重要的代码”被丢进了 .a 静态库里。今天的主角是里面的 PID 算法。

代码考古:PID

某比赛中外置的计算机里,提供了示例代码。部分代码藏了起来。一些“重要的代码”被丢进了 .a 静态库里。今天的主角是里面的 PID 算法。

对于懂 PID 算法,而且会编程的人,以下一段内容可能是废话。

PID 就是比例、积分、微分三个单词的简写。想象一个电炉,令其恒温。有温度探头,但是只允许控制电热丝的电流。对于人类来说,温度不够就加电流;温度超了就减少电流,是多是少靠经验把控;PID 就是最简单的方案,输入误差 $e$(距离预期温度的偏差值),输出控制量 $u$(比如电流)。比例项 P 表示 $u_p = k_p  e$;积分项 I 表示 $u_i = k_i \int_{0}^{T}{e \mathrm dt}$;微分项 D 表示 $u_d = k_d \frac{\mathrm de}{\mathrm dt}$,三项相加即可。然后剩下的问题就是调节三个 $k$ 的值。

贴一份我的早期代码存档,作为离散版本 PID 的简单示范。

typedef struct {
    float K_p, K_i, K_d;
    float accm_error;
    float prev_error;
    bool initialised;
} PidContext;

float PidUpdate(PidContext *ctx, float e) {
    float delta_e, param_p, param_i, param_d;
    if(likely(ctx->initialised))
        delta_e = e - ctx->prev_error;
    else
        delta_e = 0;
    param_p = ctx->K_p * e;
    param_i = ctx->K_i * e + ctx->accm_error;
    param_d = ctx->K_d * delta_e;
    ctx->accm_error = param_i;
    ctx->prev_error = e;
    return param_p + param_i + param_d;
}

而故事的主角在这,根据提交时间戳,是 2019 年 1 月:

#ifndef PID_H
#define	PID_H

/***********************************/
typedef struct PID
{              //结构体定义
    double SetPoint;         //设定值
    double Proportion;         // Proportion 比例系数
    double  Integral;            // Integral   积分系数
    double  Derivative;          // Derivative  微分系数
    double  LastError;          // Error[-1]  前一拍误差
    double  PreError;           // Error[-2]  前两拍误差


}PID;


void PIDInit (struct PID *pp);
double PIDCal(struct PID *pp, double ThisError);




#endif

struct PID 这种写法在这份代码里的合理性(这是一份 C++ 程序)以及混乱的换行,让人不得不怀疑这是哪里抄来的代码……

考古的第一步是找到特征字符串。

这个函数名本身就很特殊,计算 calculate 通常简写成 calc 但这里却是进一步写成了 PIDCal;其次,结构体的成员名也是特征;最后,它用到了 double 类型。

百度搜索一下这些关键词看看什么结果,越古老我越欢喜。

搜索结果
https://www.arduino.cn/thread-19582-1-1.html

2016 年,好的开始。加上日期限制:

第二次搜索结果

2011 年都出来了,点开来看看。

2009,十周年庆典。

https://www.cnblogs.com/Yz81128/p/3305012.html

这个的稍微再早一些。但是值得注意的是,代码开始出现大量英文了。

于是抓起其中一句 “While the PID function works, main is just a dummy program showing” 丢进 Google。一个搜索的技巧是加双引号,可以要求完全匹配整串文字。

在 Google 里搜索全句

第一条结果很有趣。

https://www.edaboard.com/showthread.php?75438-Explanation-of-a-C-code-pp-amp-sPID-memset

不仅时间推到了 2006 年,这位发帖人提问这段代码的含义,说明作者依然另有其人。另外一个重点是,这里用的函数叫做 PIDCalc,注意最后有字母 c。

让我们加上时间限制,继续搜索:

加上时间限制

第一条第二条没什么线索,第三条是一个幻灯片。这个幻灯片中这一段代码的开头注释引起了我的注意:

第一次出现了“作者”这个字样。不难想象,其他的代码可能是将其抹去之后再发布的。因为搜索引擎对于一些古老的页面可能没有记录日期,所以我们带着这个作者名,查询“这个版本”的代码。

去掉时间限制,加上作者信息

第二条结果实际上比上次见到的 EDABoard 论坛的那个帖子还要晚一个月。只是好在有附加作者信息,这个人也是从别处找来的:

第一条搜索结果是什么呢? 看起来很像最原始的版本:

http://ftp1.digi.com/support/documentation/rabbit_pid.c

这里混乱的换行和对齐让我想到某些古老的编辑器。说起来为什么这个文件的名字那么可爱?兔子 PID 是什么说法?

说起来这个域名的 www.digi.com 是个什么公司?点开看看,似乎是个物联网解决方案公司……似乎和 PID 确实有点关联。

继续搜索一下 Greg Young 这个人吧。

嵌入式……Rabbit 3000 Microprocessor…… Rabbit 是个处理器的商标?

查一下 Rabbit 半导体公司,发现 2006 年被 Digi International 收购了,而 Z-World 和 Rabbit 本身就是一块的。这就说得通了。

仔细阅读第一条,看到了这个开发者和这家公司的渊源:

https://books.google.co.jp/books?id=4SkqFozu80MC&pg=PA37&dq=Greg+Young+Z-World

1987 年,Norm 从圣何塞(加州)搬到了戴维斯(加州)。他的儿子上学了,戴维斯的公立学校比那些圣何塞的学校好。Z-World 租下了 111 平方的办公室,办公室里只有四个职工。
作为加州大学的所在地,戴维斯吸引着年轻的大学生。一位叫做 Greg Young 的年轻的计算机科学专业的大学生联系了 Norm。这位大学生独自编写了若干 Z80 嵌入式系统的开发工具。Norm 注意到了他的过人天赋和激情。Z-World 很快就雇用了他。Greg 写出了第一个 Dynamic C 的集成开发环境。后来,Greg 设立、带队公司的技术支持小组,并且担当设计工程师。今天,他依然在一些特殊项目中负责顾问。

那么这份代码一定不早于 1987 年,因为那个时候 Greg 还不属于 Z-World,代码里不会这么写。我们现在可以框定一个合适的时间范围 1987-2006。

回头看看 Digi 网站,我们从前端找找如何定位到刚才那个文件。

http://ftp1.digi.com/support/documentation/rabbit_pid.c 的入口点

然后用 web.archive.org 回溯网站的历史,

http://www.rabbit.com/support/downloads/downloads_feat.shtml

再往前域名发生变化,

http://www.rabbitsemiconductor.com/support/downloads/downloads_feat.shtml

看来这份代码至少也在 2003 年 6 月以前出现了。内容始终没有变。

最终,在 Z-World 的官方论坛上找到了这份代码的最原始版本。

http://www.zworld.com/support/bb/messages/14/3591.html?1023224636

文件地址位于 https://web.archive.org/web/20020627222140/http://www.zworld.com/support/bb/messages/14/pid_0_c.unk

它的开头是这样的

/***************************************************************************\
Copyright Z-World Inc. June 2002

	PID Function

This program has been written by the Technical Support Staff at Z-World in
response to several customer requests.  As such, it has NOT had the testing and
validation procedures which our "standard" software products have.  It is being
made available as a sample.  There is no warranty, implied or otherwise.

	The PID (Proportional Integral Derivative) function is used in mainly
	control applications. PIDCalc performs one iteration of the PID
	algorithm.

	While the PID function works, main is just a dummy program showing
	a typical usage.
\***************************************************************************/

typedef struct PID
{
	double	SetPoint;		// Desired Value

	double	Proportion;		// Proportional Const
	double	Integral;		// Integral Const
	double	Derivative;		// Derivative Const

	double	LastError;		// Error[-1]
	double	PrevError;		// Error[-2]
	double	SumError;		// Sums of Errors
}  PID;

/*=========================================================================*\
	Perform One PID Iteration
\*=========================================================================*/

double PIDCalc ( PID	*pp,   double NextPoint )

{  double	dError, Error;

	pp->SumError += (Error = pp->SetPoint - NextPoint);
	dError = pp->LastError - pp->PrevError;
	pp->PrevError = pp->LastError;
	pp->LastError = Error;
	return	(  pp->Proportion * Error
					+  pp->Integral * pp->SumError
					+  pp->Derivative * dError
				);
}

/*=========================================================================*\
	Initialize PID Structure
\*=========================================================================*/

void PIDInit  ( PID *pp )

{	  memset ( pp,0,sizeof(PID) );
}

注释写得很明白了,这份 PID 算法代码是 Z-World 的员工(而我们也看到,后来的版本中写明了具体作者是 Greg Young)在论坛回应客户时的代码。第一次发布的时间最可能是 2002 年 6 月,至今已有 17 年。

为什么这么简单东西也要抄一个世纪