使用示波器播放badapple!

使用示波器播放badapple!

说实话我想做这个很久了,但是一直太懒不去查资料,就一直在鸽.高中住校生活实在是太无聊了,于是在课堂走神的时候大概构思了一下.结果回来查资料实践的时候发现核心的东西cv2都有了,复杂的算法根本不需要自己写.

原理

众所周知,示波器xy模式可以显示李萨如图像

这是一张李萨如图形的示意图:

李萨如图形示意

可以看出两路信号分别独立控制两个轴上的取值,任意时刻两个信号的值合起来就能确定一个点了.

只要我们把我们要组成的图像上的所有点都反映在这两路信号里就可以了.比如我有四个点:

1
A(114,514) B(191,981) C(514,114) D(981,191)

那么这两路信号就是

1
2
X通道:(114,191,514,981)
Y通道:(514,981,114,191)

显示出来就是以这四个点为顶点的四边形.

预处理

准备一个badapple的视频

用ffmpeg将所有帧抽出来:

1
2
mkdir raw && cd raw
ffmpeg -i ../badapple_music.mp4 -r 30 %04d.png

抽完后的目录内容:

1
2
3
4
5
6
7
8
9
10
11
12
$ tree | head
  raw
    ├── 0001.png
    ├── 0002.png
    ├── 0003.png
    ├── 0004.png
    ├── 0005.png
    ├── 0006.png
    ├── 0007.png
    ├── 0008.png
    ├── 0009.png
    ....

显示图片

要想显示视频,首先要能显示图片.

显然我们在示波器上要显示的肯定只能是目标图像内容的"边界",即内容的描边,其对于badapple这种简单的视频是非常简单的.

cv2给我们提供了一个专门寻找边界的函数:findContours().它的用法非常简单,只需要把二值化的图片和一些相关的参数传进去就搞定了.有关它的基础用法可以参考这篇文章:https://blog.csdn.net/weixin_44690935/article/details/109008946

现在我们可以写出以下的代码,它找出了图像的所有边界点坐标:

1
2
3
img = cv2.imread(src, 0) #src是图片
_, thosd = cv2.threshold(img, 127, 255, cv2.THRESH_BINARY) #原图二值化
points, _ = cv2.findContours(thosd, cv2.RETR_TREE, cv2.CHAIN_APPROX_NONE) #找出边界点

返回的这些点是排好序的,可以不用二次写个排序算法.

wav的pcm段左右声道是挨在一起的,类似[左采样1,右采样1,左采样2,右采样2,左采样2,右采样2,左采样2,右采样2,.....左采样n,右采样n]

这说明我们我们把把每个点的xy采样放在一起即可.

返回的点的数组有好几维,我们把它展开成想要的样子:

1
f_array = np.array(np.concatenate(points).flatten('C')) #展开成[x1,y1,x2,y2,x3,y3...]
1
2
3
with wave.open('result.wav',mode='wb') as f:
    f.setparams((2,2,200000,0,"NONE","NOT COMPRESSED")) #从左往右依次为双声道,量化精度两字节,采样率200k,不压缩,不压缩
    f.writeframes((f_array*40).astype(np.int16)) #40是一个音量系数,用以在合理的范围内不失真地扩大图像大小以及精度(原来都是浮点数,转为整型会丢失精度),太大会超出量化精度的范围导致极大失真

这样出来的波形就可以直接放进fl里看效果了.

输出波形

晚上给人录的视频里截的,凑活看吧()

因为我没有正经的示波器所以用fl的Wave Candy代替

显示效果

fl的这个示波器比较特殊,出来的图像会是这样,成倾斜45度角.因为fl里的这个示波器R轴和L轴与横平竖直的XY轴成45度,所以最后点组成的图像也成了45度.当然如果是正经示波器就木有这种烦恼.

根据坐标轴旋转的有关数学公式:https://zhuanlan.zhihu.com/p/339668569

可将代码做如下修改:

1
2
3
4
5
f_array = np.array([])
points = np.array(np.concatenate(points).flatten('C')) #展开成[x1,y1,x2,y2,x3,y3...]
points = points.reshape(points.size//2,2) #变成 [[x1,y1],[x2,y2],[x3,y3],...,[x114,y114]]
for x,y in points:
    f_array = np.concatenate((f_array, np.array([-(y*SIN_45+x*COS_45),-(y*COS_45-x*SIN_45)])))

这样再出来的图像就是正的了,不过引入三角函数会让处理速度降低,毕竟是浮点数计算_.不过可以通过提前把COS_45和SIN_45提前算好节省一部分计算量加快一点点.

变换后的图像

为了让图像稍微居中一下,在原坐标上加点偏移(凭感觉试,懒得推):

1
points = points.reshape(points.size//2,2)-256 #-256是图像在示波器中的左右上下偏移

代码:

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
27
28
29
30
31
32
33
import cv2
import numpy as np
import wave
import os
FRAME_RATE=100000

SIN_45=np.sin(np.pi/4)
COS_45=np.cos(np.pi/4)

def procimg(src:str):
    img = cv2.imread(src, 0)
    _, thosd = cv2.threshold(img, 127, 255, cv2.THRESH_BINARY)
    points, _ = cv2.findContours(
        thosd, cv2.RETR_TREE, cv2.CHAIN_APPROX_NONE)
    if len(points) == 0:
        return np.array([])
    f_array = np.array([])
    points = np.array(np.concatenate(points).flatten('C')) #展开成[x1,y1,x2,y2,x3,y3...]
    points = points.reshape(points.size//2,2)-256 #变成 [[x1,y1],[x2,y2],[x3,y3],...,[x114,y114]],-256是图像在示波器中的左右上下偏移
    for x,y in points:
        f_array = np.concatenate(
            (f_array, np.array([-(y*SIN_45+x*COS_45),-(y*COS_45-x*SIN_45)])))  #旋转
    return f_array

f_array=np.array([])

# for i in sorted(os.listdir('raw')):
for i in ['0240.png']: #这是个测试用的图
    print(i)
    f_array=np.concatenate((f_array,procimg('raw/'+i)))
with wave.open('result.wav',mode='wb') as f:
    f.setparams((2,2,FRAME_RATE,0,"NONE","NOT COMPRESSED"))
    f.writeframes((f_array*40).astype(np.int16))

处理视频

一般

显然,只要把每一帧的输出按顺序拼起来就可以达到目的了,那么用上文32行注释代码替代33行调试用代码,其余保持不变即可:

1
2
for i in sorted(os.listdir('raw')):
#for i in ['0240.png']: #这是个测试用的图

更加平滑的帧速

不难发现,仅仅简单堆叠每帧图像会导致帧速不统一,画面复杂/细节多的帧会明显慢于细节简单的帧.

我们可以通过把每一帧所包含的点的数量"展开"到一个统一的数,这样每一帧画面所经过的时间自然而然就相等了.

因此我们可以将procimg函数改为如下:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
def procimg(src:str):
    img = cv2.imread(src, 0)
    _, thosd = cv2.threshold(img, 127, 255, cv2.THRESH_BINARY)
    points, _ = cv2.findContours(
        thosd, cv2.RETR_TREE, cv2.CHAIN_APPROX_NONE)
    if len(points) == 0:
        return np.array([0]*AVG_FRAMESIZE*2)
    f_array = np.array([])
    points = np.array(np.concatenate(points).flatten('C')) #展开成[x1,y1,x2,y2,x3,y3...]
    points = (points.reshape(points.size//2,2)-256).tolist() #变成 [[x1,y1],[x2,y2],[x3,y3],...,[x114,y114]],-256是图像在示波器中的左右上下偏移
    perc=1
    choosed=[]
    if (len_of_points:=len(points))<AVG_FRAMESIZE:
        perc=math.floor(AVG_FRAMESIZE/len_of_points)
        choosed=random.sample(points,(AVG_FRAMESIZE%len_of_points))
    # print(len_of_points,perc,len(choosed))
    for x,y in points:
        ch=0
        if [x,y] in choosed:
            ch=1
        f_array = np.concatenate(
            (f_array, np.array([-(y*SIN_45+x*COS_45),-(y*COS_45-x*SIN_45)]*(perc+ch))))  #旋转
    # print(len(f_array))
    return f_array

其中AVG_FRAMESIZE为常量,意为每帧至少需包含多少对点;这里随便酌情取作7000.理论越大越好,不过过大会导致处理速度降低和生成文件变大.

我们还可以让程序自己计算生成音频的采样率,以让最终的结果与原版badapple保持同样的时长:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
VID_FPS=30
f_array=np.array([])
frames=sorted(os.listdir('raw'))
# frames=['1430.png']
VID_FRAMES=len(frames) #多少帧
VID_DUR_TIME=VID_FRAMES/VID_FPS #多少秒
for i in frames:
    print(i)
    f_array=np.concatenate((f_array,procimg('raw/'+i)))
SAMPLE_CNT=len(f_array)/2 #采样点按左右声道两个为一对取总对数计
AUDIO_FRAME_RATE=SAMPLE_CNT//VID_DUR_TIME
print(SAMPLE_CNT,AUDIO_FRAME_RATE)
with wave.open('badapple.wav',mode='wb') as f:
    f.setparams((2,2,AUDIO_FRAME_RATE,0,"NONE","NOT COMPRESSED"))
    f.writeframes((f_array*40).astype(np.int16))

嗯。。。这个东西如果你试一下就会发现慢到怀疑人生。。巨大多内存复制。

性能优化

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
27
28
29
30
31
32
33
...

    if (len_of_points:=len(points))<AVG_FRAMESIZE:
        perc=math.floor(AVG_FRAMESIZE/len_of_points)
        random.seed(114514)   #使用固定seed使得抽样结果固定，最终效果为相同输入相同输出。wav的hash为一个固定值。
        choosed=set([str(i) for i in random.sample(points,(AVG_FRAMESIZE%len_of_points))]) #使用set,后续查询时间复杂度O(1)
    
    cap = (len_of_points*perc+len(choosed))*2
    f_array = [None]*cap  #预分配
    cur_point = 0
    for x,y in points:
        ch=0
        if str([x,y]) in choosed: #愉快地O(1)
            ch=1
        f_array[cur_point:cur_point+(perc+ch)*2] = [-(y*SIN_45+x*COS_45),-(y*COS_45-x*SIN_45)]*(perc+ch)  #旋转
        cur_point+=(perc+ch)*2
...

frames=sorted(os.listdir('raw'))
# frames=['1430.png']
VID_FRAMES=len(frames)
f_array=[None]*VID_FRAMES*AVG_FRAMESIZE*2 #预分配内存，减少内存复制开销

...

    frame_len=(frame_data:=procimg('raw/'+i)).__len__()
    f_array[cur_len:cur_len+frame_len]=frame_data.tolist() #按索引直接插入
    print(f'set f_array[{cur_len}:{cur_len+frame_len}]')
    cur_len+=frame_len

...

    f.writeframes((np.array(f_array)*40).astype(np.int16))

改用hash表来加速，使用预分配好的内存空间降低内存开销，其实还有优化空间，不过已经够了。

最终代码:

无性能优化

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
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
import cv2
import numpy as np
import wave
import os
import math
import random

AVG_FRAMESIZE=7000 #每帧多少对坐标(采样)

VID_FPS=30


SIN_45=np.sin(np.pi/4)
COS_45=np.cos(np.pi/4)

def procimg(src:str):
    img = cv2.imread(src, 0)
    _, thosd = cv2.threshold(img, 127, 255, cv2.THRESH_BINARY)
    points, _ = cv2.findContours(
        thosd, cv2.RETR_TREE, cv2.CHAIN_APPROX_NONE)
    if len(points) == 0:
        return np.array([0]*AVG_FRAMESIZE*2)
    f_array = np.array([])
    points = np.array(np.concatenate(points).flatten('C')) #展开成[x1,y1,x2,y2,x3,y3...]
    points = (points.reshape(points.size//2,2)-256).tolist() #变成 [[x1,y1],[x2,y2],[x3,y3],...,[x114,y114]],-256是图像在示波器中的左右上下偏移
    perc=1
    choosed=[]
    if (len_of_points:=len(points))<AVG_FRAMESIZE:
        perc=math.floor(AVG_FRAMESIZE/len_of_points)
        choosed=random.sample(points,(AVG_FRAMESIZE%len_of_points))
    # print(len_of_points,perc,len(choosed))
    for x,y in points:
        ch=0
        if [x,y] in choosed:
            ch=1
        f_array = np.concatenate(
            (f_array, np.array([-(y*SIN_45+x*COS_45),-(y*COS_45-x*SIN_45)]*(perc+ch))))  #旋转
    # print(len(f_array))
    return f_array

f_array=np.array([])
frames=sorted(os.listdir('raw'))
# frames=['1430.png']
VID_FRAMES=len(frames)
VID_DUR_TIME=VID_FRAMES/VID_FPS
for i in frames:
    print(i)
    f_array=np.concatenate((f_array,procimg('raw/'+i)))
SAMPLE_CNT=len(f_array)/2
AUDIO_FRAME_RATE=SAMPLE_CNT//VID_DUR_TIME
print(SAMPLE_CNT,AUDIO_FRAME_RATE)
with wave.open('badapple.wav',mode='wb') as f:
    f.setparams((2,2,AUDIO_FRAME_RATE,0,"NONE","NOT COMPRESSED"))
    f.writeframes((f_array*40).astype(np.int16))

执行时间：18m-9.6s = 1070.4s（使用ohmyzsh的timer插件）

FpS=65731070.4≈6.14F_{pS}=\frac{6573}{1070.4}\approx6.14

嗯。慢得惊天地泣鬼神。Latex写公式显得上档次

性能优化

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
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
import cv2
import numpy as np
import wave
import os
import math
import random

AVG_FRAMESIZE=7000 #每帧多少对坐标(采样)

VID_FPS=30


SIN_45=np.sin(np.pi/4)
COS_45=np.cos(np.pi/4)

def procimg(src:str):
    img = cv2.imread(src, 0)
    _, thosd = cv2.threshold(img, 127, 255, cv2.THRESH_BINARY)
    points, _ = cv2.findContours(
        thosd, cv2.RETR_TREE, cv2.CHAIN_APPROX_NONE)
    if len(points) == 0:
        return np.array([0]*AVG_FRAMESIZE*2)
    points = np.array(np.concatenate(points).flatten('C')) #展开成[x1,y1,x2,y2,x3,y3...]
    points = (points.reshape(points.size//2,2)-256).tolist() #变成 [[x1,y1],[x2,y2],[x3,y3],...,[x114,y114]],-256是图像在示波器中的左右上下偏移
    perc=1
    choosed=[]
    if (len_of_points:=len(points))<AVG_FRAMESIZE:
        perc=math.floor(AVG_FRAMESIZE/len_of_points)
        random.seed(114514)
        choosed=set([str(i) for i in random.sample(points,(AVG_FRAMESIZE%len_of_points))])
    cap = (len_of_points*perc+len(choosed))*2
    f_array = [None]*cap #预分配
    # print(cap)
    cur_point = 0
    for x,y in points:
        ch=0
        if str([x,y]) in choosed:
            ch=1
        f_array[cur_point:cur_point+(perc+ch)*2] = [-(y*SIN_45+x*COS_45),-(y*COS_45-x*SIN_45)]*(perc+ch)  #旋转
        
        cur_point+=(perc+ch)*2
    # print(len(f_array))
    return f_array

frames=sorted(os.listdir('raw'))
# frames=['1430.png']
VID_FRAMES=len(frames)
f_array=[None]*VID_FRAMES*AVG_FRAMESIZE*2 #预分配内存
VID_DUR_TIME=VID_FRAMES/VID_FPS

print(VID_FRAMES*AVG_FRAMESIZE*2)
cur_len=0  #now pointer
for i in frames:
    print(i)
    # procimg('raw/'+i)
    frame_len=(frame_data:=procimg('raw/'+i)).__len__()
    f_array[cur_len:cur_len+frame_len]=frame_data
    print(f'set f_array[{cur_len}:{cur_len+frame_len}]')
    cur_len+=frame_len
SAMPLE_CNT=len(f_array)/2
AUDIO_FRAME_RATE=SAMPLE_CNT//VID_DUR_TIME
print(SAMPLE_CNT,AUDIO_FRAME_RATE)
with wave.open('badapple.wav',mode='wb') as f:
    f.setparams((2,2,AUDIO_FRAME_RATE,0,"NONE","NOT COMPRESSED"))
    f.writeframes((np.array(f_array)*40).astype(np.int16))

执行时间： 2m-27.9s = 92.1s（使用ohmyzsh的timer插件）

FpS=657392.1≈71.36F_{pS}=\frac{6573}{92.1}\approx71.36

效率提升了10倍甚至⑨倍！

效果

帧速不统一的

帧速统一的

https://www.bilibili.com/video/BV1au411a75Y/

#python #东方 #示波器 #opencv
使用示波器播放badapple!
https://www.hakurei.org.cn/2023/06/21/lissajous-badapple/
作者
zjkimin
发布于
2023年6月21日
更新于
2023年8月25日
许可协议