二次IIR滤波器
在声音信号处理中,经常需要设计滤波器,二次IIR滤波器由于其参数设计简单,运算量少,运用的最为广泛。这里介绍几种常用的滤波器的参数计算。
Peak Equalizer
Peak Equalizer的设计参数为中心频率f,最大增益db,形状参数Q。下面这个Python小程序通过这三个参数计算滤波器的参数a和b,取样频率为44100Hz。
01def peak(f, db, Q):
02 A = 10.0**(db/40.0)
03 w0= 2*math.pi*f/44100.0
04 alpha = math.sin(w0)/2/Q
05 b0 = 1+alpha*A
06 b1 = -2*math.cos(w0)
07 b2 = 1-alpha*A
08 a0 = 1+alpha/A
09 a1 = -2*cos(w0)
10 a2 = 1-alpha/A
11 b0 /= a0
12 b1 /= a0
13 b2 /= a0
14 a1 /= a0
15 a2 /= a0
16 a0 /= a0
17 return [b0,b1,b2],[a0,a1,a2]

高、低通滤波器
二次高、低通滤波器的滤波倾斜度为-12dB/oct ,只有一个参数:滤波频率f。
01def lpf(f):
02 a = [1.0, math.sqrt(2), 1.0]
03 b = [0.0, 0.0, 1.0]
04 ad = [0.0] * 3
05 bd = [0.0] * 3
06 fs = 44100.0
07 f = 2*math.pi*f
08 b[2] *= f*f
09 a[1] *= f
10 a[2] *= f*f
11 ad[0] = a[2] - a[1]*2*fs + 4*fs*fs*a[0]
12 ad[1] = a[2]*2 - 8*fs*fs*a[0]
13 ad[2] = a[2] + 2*fs*a[1] + 4*fs*fs*a[0]
14 bd[0] = b[2] - b[1]*2*fs + 4*fs*fs*b[0]
15 bd[1] = b[2]*2 - 8*fs*fs*b[0]
16 bd[2] = b[2] + b[1]*2*fs + 4*fs*fs*b[0]
17 t = ad[2]
18 return [bd[2]/t, bd[1]/t, bd[0]/t], [ad[2]/t, ad[1]/t, ad[0]/t]
01def hpf(f):
02 a = [1.0, math.sqrt(2), 1.0]
03 b = [0.0, 0.0, 1.0]
04 ad = [0.0] * 3
05 bd = [0.0] * 3
06 fs = 44100.0
07 f = 2*math.pi*f
08 b[0],b[2] = b[2],b[0]
09 a[1] *= f
10 a[2] *= f*f
11 ad[0] = a[2] - a[1]*2*fs + 4*fs*fs*a[0]
12 ad[1] = a[2]*2 - 8*fs*fs*a[0]
13 ad[2] = a[2] + 2*fs*a[1] + 4*fs*fs*a[0]
14 bd[0] = b[2] - b[1]*2*fs + 4*fs*fs*b[0]
15 bd[1] = b[2]*2 - 8*fs*fs*b[0]
16 bd[2] = b[2] + b[1]*2*fs + 4*fs*fs*b[0]
17 t = ad[2]
18 return [bd[2]/t, bd[1]/t, bd[0]/t], [ad[2]/t, ad[1]/t, ad[0]/t]


shelving滤波器
Shelving滤波器提升或者抑制高频或者低频的增益,有三个参数: 中心频率f,增益db,类型type。type为0时,对低频进行增益处理,为1时对高频进行增益处理。
01def shelving(f, db, type):
02 if type==1:
03 f = 22050 - f
04 Wb = 2.0*f/44100.0*math.pi
05 K = math.tan(Wb/2)
06 g = 10.0 ** (db/20.0)
07 V = g ** (1.0/2.0) - 1
08 a0 = 1 + 2*K + K*K
09 a1 = 2*K*K - 2
10 a2 = 1 - 2*K + K*K
11 b0 = a0
12 b1 = a1
13 b2 = a2
14 b0 = b0 + 2*V*K*(K+1)
15 b1 = b1 + 2*V*K*(2*K)
16 b2 = b2 + 2*V*K*(K-1)
17 b0 = b0 + V*V*K*K*1
18 b1 = b1 + V*V*K*K*2
19 b2 = b2 + V*V*K*K*1
20 b0 = b0 / a0
21 b1 = b1 / a0
22 b2 = b2 / a0
23 a1 = a1 / a0
24 a2 = a2 / a0
25 a0 = a0 / a0
26 if type==1:
27 b1 = -b1
28 a1 = -a1
29 return [b0,b1,b2],[a0,a1,a2]

源程序高亮显示代码公布
经过一段时间的调试,高亮度显示Python源代码的程序似乎已经没有什么BUG了,于是就给网站添加了一个功能: 用这段显示本网站的代码,不过由于安全方面的考虑,目前只能浏览如下3个源文件,今后将用这种方式发布一些Python写的小程序。
colorit.py : 高亮度显示Python源代码
showpy.py : 调用colorit.py显示网站的源代码
sudokuSolver.py : 数独游戏求程序
Python源程序的高亮度化
最近博客中的Python源程序越贴越多,但是由于不支持源程序的高亮度显示,因此阅读大段的程序比较困难。
很早就想在博客中加入源程序的高亮度显示功能。也找到了不少Python的库可以完美地完成这些功能。可是有一个问题:因为服务器速度慢、内存少,为了一个高亮度显示而载入一大堆库实在是负担太大。
于是我就着手自己写一个最最简单的源程序高亮度化的Python小程序。这个问题初看上去似乎比较困难,可是仔细思考之后才发现原来很简单。不到1个小时这个程序就完成了,而且也只有九十多行,处理速度也很快, 因此我就不考虑缓存结果,而是每次显示的时候都重新计算一遍。
下面是一个高亮度化的例子,可以看看效果:
※程序还是有一些小BUG,所以暂时不公开 (源程序已经公布 2007/9/22)
class Dummy(object):
def __init__(self):
self.n = 1000
def go(self):
for i in range(1000):
print "number %s:" % i
用回溯法解迷宫谜题
按照下面的规则,找到从开始方格(红色)到目标方格(绿色)的路径。
| 规则
|
例如上讲这个题目的解答为:
2 | 6 | 4 | 1 | 2 | 2 | |||||||||||||||||||||||||||||||||||||
| 3 | |||||||||||||||||||||||||||||||||||||||||
3 | ||||||||||||||||||||||||||||||||||||||||||
2 | ||||||||||||||||||||||||||||||||||||||||||
1 | ||||||||||||||||||||||||||||||||||||||||||
3 | ||||||||||||||||||||||||||||||||||||||||||
5 |
下面是一个10×10的谜题,这个题目有些难度,手工解决花了不少时间,于是就想用Python写一个程序解解看。
4 | 1 | 9 | 4 | 6 | 5 | 3 | 3 | 4 | 8 | |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
| 4 | |||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
6 | ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
7 | ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
4 | ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
4 | ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
3 | ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
3 | ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
7 | ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
1 | ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
8 |
解这种题目最简单的方法就是回溯法,回溯法最简单的实现就是用递归。这里先给出完整的程序,然后再简单介绍一下编程思路。
01from sets import Set
02
03class Maze:
04 def __init__(self, N):
05 self.maze = Set()
06 self.sx = [0]*N
07 self.sy = [0]*N
08 self.N = N
09
10 def add(self, pos):
11 self.maze.add(pos)
12 self.sx[pos[0]] += 1
13 self.sy[pos[1]] += 1
14
15 def remove(self, pos):
16 self.maze.remove(pos)
17 self.sx[pos[0]] -= 1
18 self.sy[pos[1]] -= 1
19
20 def check(self, pos, X, Y):
21 x,y = pos
22 if self.sx[x] > X[x]: return False
23 if self.sy[y] > Y[y]: return False
24 if self.sx[x] == X[x]:
25 for i in xrange(x):
26 if self.sx[i]!=X[i]: return False
27 if self.sy[y] == Y[y]:
28 for i in xrange(y):
29 if self.sy[i]!=Y[i]: return False
30 return True
31
32 def printout(self):
33 for i in xrange(self.N):
34 for j in xrange(self.N):
35 if (j,i) in self.maze: print "*",
36 else: print ".",
37 print
38
39def near_pos(pos, N):
40 x,y=pos
41 if x-1>=0: yield x-1,y
42 if x+1<=N-1:yield x+1,y
43 if y-1>=0: yield x,y-1
44 if y+1<=N-1:yield x,y+1
45
46def maze_step(maze, pos, target, N):
47 maze.add(pos)
48 if maze.check(pos,X,Y):
49 if pos == target:
50 maze.printout()
51 return
52 for npos in near_pos(pos, N):
53 if npos not in maze.maze:
54 maze_step(maze, npos, target, N)
55 maze.remove(pos)
56
57X = [4,1,9,4,6,5,3,3,4,8]
58Y = [4,6,7,4,4,3,3,7,1,8]
59maze_step(Maze(10),(0,0),(9,9),10)
Maze类的maze属性是一个集合(Set),用来记录当前走过的所有的方格的坐标。sx和sy属性分别用来记录maze中每列以及每行的方格数。add和remove方法用来添加和删除某个方块,我们看到添加和删除方块的时候,同时更新sx和sy的相应坐标的值。
check方法用来检测当前的走法是否符合规则。pos参数为最后一个添加进去的方格坐标,也就是需要检查的对象。我们先看看这个方格所在的行和列的和是否满足条件:
if self.sx[x] > X[x]: return False
if self.sy[y] > Y[y]: return False
由于入口在左上,而出口在右下,因此一旦某一行或者列的方格数等于了规则所给出的方格数,那么它左边的列或者上面的行就必须也符合规定的数目,因为如果要走回去改变方格数的话,那么就必须再经过一次已经等于规定数目的行或者列,这样就不能满足条件了。这个条件很关键,如果去除掉这个条件检测的话,那么等几个小时可能都出不来结果。
if self.sx[x] == X[x]:
for i in xrange(x):
if self.sx[i]!=X[i]: return False
if self.sy[y] == Y[y]:
for i in xrange(y):
if self.sy[i]!=Y[i]: return False
maze_step完成递归回溯功能,注意在回溯(也就是函数返回)之前,需要通过执行remove方法删除掉对maze的修改。
最后这个程序的输出为:
* . * * . * . . * *
* . * . . * * * * *
* * * . . . . . . *
. . * * * . . . . *
. . * . * . . . . *
. . * . * . . . . *
. . * . * * * * * *
. . * . . . . . . .
. . * * * * * * * *
方便起见,程序没有记录行走路径,只输出路径所经历过的方块,如果需要记录行走路径的话,那就需要在Maze类中再添加一个path属性,它是一个列表,在add和remove中分别对其进行添加和删除的操作即可。
爸爸回家看海月
今年8月份是我最高兴的一个月, 因为爸爸妈妈都回家来看我了。特别是爸爸,我已经有一年多没有和他见面了。在这一年多的时间里,我们只能靠网络聊天相互熟悉对方的声音,靠照片相互熟悉对方的相貌。因此当我看见爸爸风尘仆仆地走进来的时候,我就躲在妈妈身后,打量这个既陌生又熟悉的爸爸。
妈妈问我:“海月,这个人是谁?”
我扭扭捏捏,挤出一句“爸爸”,就逃到别的房间去了。
我能一眼就认出爸爸来,令他很欣慰,不过接下来我的表现,就让他有些着急了。我想,这么久你都没有来看我,绝对不能马上和你靠近乎,所以我对他的态度若即若离。先是一见到爸爸就假装哭,令他手足无措,只好拿好吃的来哄我。我看爸爸的态度很诚恳,就答应让他亲我了。
不过亲归亲,爸爸想抱我,我可是不让的。以至于好几次我们全家出去玩都是妈妈抱着我,把她累了个够呛。我越不让爸爸抱我,他就越想抱我,想尽各种办法逗我开心。这样自然而然地,我们父女之间的感情就一下子好了起来。爸爸也就在不知不觉中将我抱了起来。
爸爸抱我和别人不一样,他让我跳蹦极(我轻轻一跳,爸爸就把我高高地举起来),玩飞机(把我举在胸前,满屋子乱跑),这样的活动令我很兴奋和高兴,爸爸也是好久没有这么开心地运动了。