Skip to content

Latest commit

 

History

History

poscat

Folders and files

NameName
Last commit message
Last commit date

parent directory

..
 
 

去吧!追寻自由的电波

这题用使用默认配置的 mpv 调整速度会产生回声,进而导致听不清... google 搜索 "mpv slow playback echo" 得知解决方法为加上 --no-audio-pitch-correction 参数(或写入 mpv.conf

p😭q

这题其实比较简单,只要读懂代码并找到合适的库即可,甚至不需要理解傅立叶变换是什么(当然,理解傅立叶变换/STFT肯定对理解代码有帮助)。

代码解读

num_freqs = 32
quantize = 2
min_db = -60
max_db = 30
fft_window_size = 2048
frame_step_size = 512
window_function_type = 'hann'
red_pixel = [255, 0, 0]
white_pixel = [255, 255, 255]
y, sample_rate = librosa.load("flag.mp3")  # sample rate is 22050 Hz

这里定义了很多参数,并把音频载入了 ysample_rate。查阅 librosa 文档可知 y 代表 audio time series,sample_rate 代表 sampling rate。

spectrogram = (
    numpy.around(
        librosa.power_to_db(
            librosa.feature.melspectrogram(
                y,
                sample_rate,
                n_mels=num_freqs,
                n_fft=fft_window_size,
                hop_length=frame_step_size,
                window=window_function_type
            )
        ) / quantize
    ) * quantize
)

这里对载入的音频文件先后应用了 melspectrogrampower_to_db 以及 around 三个函数。查阅 librosa 文档可知 melspectrogram 将音频进行短时距傅里叶变换 (STFT) 后又把时频谱 (spectrogram) 的横轴 (频率) 的单位变为了 mel,并且 n_mels (即 num_freqs) 这个参数代表了有多少个不同的频率区间。而 power_to_db 则将时频谱纵轴 (振幅/声强) 的单位从功率变为了分贝。查阅 numpy 文档可知 around 将矩阵每个元素都进行了向下取整,考虑到先除 quantize = 2 后乘 2 可知这里将分贝向下 round 到了偶数。

最后输出的就是一个时频谱,其纵坐标单位为 db, 并且被向下 round 到了偶数,横坐标为 mel。

gif_data = [
    numpy.kron(
        numpy.array([
            [red_pixel if freq % 2 and round(frame[freq // 2]) > threshold else white_pixel for threshold in list(range(
                    min_db, max_db + 1, quantize))[::-1]] for freq in range(num_freqs * 2 + 1)
        ]),
        numpy.ones([quantize, quantize, 1])
    )
    for frame in spectrogram.transpose()
]

查阅 librosa 文档可知 librosa 中时频谱的第一个下标代表频率,第二个下标代表帧数,矩阵元素代表强度,这里先对 spectrogram 进行了一个转置,这样遍历时就可以方便的遍历每一帧了。观察对每一个 frame 的操作,可以发现它把每一个 frame 变成了一个二阶张量(矩阵),其中每一个元素类型都是 RGB 值,而且只有红 (255,0,0) 和白 (255,255,255) 两种情况。也就是说这段代码把时频谱的每一帧转换(映射)为了 gif 图像的一帧。

再看具体转换方法: 第一个 for 遍历了 0 到 num_freqs * 2,第二个 for 遍历了 min_db 到 max_db, 并且步长为 2,最后判断时频谱该帧在频率 freq//2 上的强度大小是否大于 threshold,如果大于则涂红色,否则涂白色。这到底在干嘛?其实很容易看出或者猜出这在画柱状图。那为什么要遍历 0 到 num_freqs * 2,判断 freq % 2 是否是 1 又有什么用?这其实是在用纯白的一列像素分隔开各个柱子。之后这个 gif 帧又和一个 2x2x1 的全 1 三阶张量进行了一个 kron 运算。查阅 numpy 文档可知 kron 为克罗内克积。查阅维基百科可知这里克罗内克积将图片的每一个像素放大到了 2x2。打开图片正式了这一猜测:柱和柱间间隔的宽度都是两个像素。注意到遍历强度时我们的步长为 2,这里经过缩放后强度(与 min_db 的差值)正好等于纵向像素高度。

自此代码解读完毕。

实现解码

TLDR:

import numpy
from librosa import db_to_power
from librosa.feature.inverse import mel_to_audio
import PIL.Image
import soundfile as sf

# run `ffmpeg -i flag.gif -vsync 0 flag%d.png` before running this program

num_freqs = 32
quantize = 2
min_db = -60
max_db = 30
fft_window_size = 2048
frame_step_size = 512
window_function_type = 'hann'
red_pixel = (255, 0, 0, 255)
white_pixel = (255, 255, 255, 255)
sample_rate = 22050

frames = 587


def get_max_db(pix, x):
    max = min_db - 1
    for y in range(92):
        pix_v = pix[x, 91-y]
        if pix_v == red_pixel:
            max = min_db + y
        else:
            break
    return max


def build_db_frame(i):
    img = PIL.Image.open("flag{}.png".format(i))
    pix = img.load()
    return [get_max_db(pix, 4*freq + 2) for freq in range(num_freqs)]


def main():
    db_frames = numpy.array([build_db_frame(i+1) for i in range(frames)])
    db_spectrogram = db_frames.transpose()
    power_spectrogram = db_to_power(db_spectrogram)
    audio = mel_to_audio(power_spectrogram, sample_rate,
                         n_fft=fft_window_size, hop_length=frame_step_size,
                         window=window_function_type)
    sf.write("flag.wav", audio, sample_rate, subtype='PCM_24')


if __name__ == '__main__':
    main()

1. 图像处理

我们先把 gif 图片分离为每一帧,可以使用 ffmpeg:

$ ffmpeg -i flag.gif -vsync 0 flag%d.png

执行后可以看到生成了 flag1.pngflag587.png 共 587 帧图像。

接下来我们要把每一帧再转换为坐标为 mel-db 的频谱,我这里使用了 pillow 库,pillow 可以直接把每一个像素读入,但需要注意 pillow 库的坐标原点在左上方。遍历所有 4n+3 列,并记录红色像素的个数,即可得到频谱。

def get_max_db(pix, x):
    max = min_db - 1
    for y in range(92):
        pix_v = pix[x, 91-y]
        if pix_v == red_pixel:
            max = min_db + y
        else:
            break
    return max


def build_db_frame(i):
    img = PIL.Image.open("flag{}.png".format(i))
    pix = img.load()
    return [get_max_db(pix, 4*freq + 3) for freq in range(num_freqs)]

def main():
    db_frames = numpy.array([build_db_frame(i+1) for i in range(frames)])

转置后得到时频谱

db_spectrogram = db_frames.transpose()

2. 解码时频谱

接下来我们需要解码时频谱。在 google 上搜索可知 librosa 提供了 mel_to_audio 这个函数将以 mel 为横坐标单位的时频谱转换为音频,而这个函数又要求强度坐标为功率,而不是分贝,再次搜索 google 可知 librosa 也提供了一个强度和分贝之间转换的函数 db_to_power,这样我们就可以解码了:

power_spectrogram = db_to_power(db_spectrogram)
audio = mel_to_audio(power_spectrogram, sample_rate,
                     n_fft=fft_window_size,
                     hop_length=frame_step_size,
                     window=window_function_type)

3. 输出音频文件

最后我们需要输出音频文件,在 google 上搜索可以找到 librosa 文档中 "Advanced I/O Use Cases" 一节,其中提到了使用 PySoundFile 可以非常方便的将音频输出为各种格式:

sf.write("flag.wav", audio, sample_rate, subtype='PCM_24')

最后打开 flag.wav 即可得到 flag。

Easy RSA

这题在理解 RSA 原理,以及一些初等数论后并不难(主要是欧拉函数的性质)。你只需找出 p, q, 求出 φ(n) = (p - 1)(q - 1) (注:φ 为欧拉函数) 后求出 e⁻¹ mod φ(n),即私钥,即可解密密文。

求 p

y 很大,因而计算 y! (mod x) 及其慢,但是注意到 x, y 很接近,而且 x 为素数(我是通过 wolframalpha 知道的)联想到可以使用 wilson's theorem 简化 y! (mod x)

求 q

难度在于求解 e 次同余方程 q^e ≡ value_q (mod n) (注意区别这里的 n 和上面提到的 n,这里的 n 指出现在 get_q 中的那个),显然,由于 q 的大小 (2²⁵⁶ 量级),暴力是不现实的。这里注意到 e 为素数(事实上这是已知最大的费马素数),推出 gcd(e, φ(n)) = 1,推出 e⁻¹ mod φ (n) 一定存在。联想到 RSA 解密过程有关的一个引理: ed ≡ 1 (mod φ(n)) ⇒ ∀x∈ℤ, x^(ed) ≡ x (mod n),我们可以通过计算 value_q 的 e⁻¹ 次幂解密得到 q。

代码

import sympy
from Crypto.Util.number import *

c = 110644875422336073350488613774418819991169603750711465190260581119043921549811353108399064284589038384540018965816137286856268590507418636799746759551009749004176545414118128330198437101472882906564195341277423007542422286760940374859966152871273887950174522820162832774361714668826122465471705166574184367478

value_q = 5591130088089053683141520294620171646179623062803708281023766040254675625012293743465254007970358536660934858789388093688621793201658889399155357407224541324547522479617669812322262372851929223461622559971534394847970366311206823328200747893961649255426063204482192349202005330622561575868946656570678176047822163692259375233925446556338917358118222905050574458037965803154233167594946713038301249145097770337253930655681648299249481985768272321820718607757023350742647019762122572886601905212830744868048802864679734428398229280780215896045509020793530842541217790352661324630048261329493088812057300480085895399922301827190211956061083460036781018660201163819104150988531352228650991733072010425499238731811243310625701946882701082178190402011133439065106720309788819

e = 65537

def get_p():
  # https://cs.stackexchange.com/questions/1495/what-is-the-most-efficient-way-to-compute-factorials-modulo-a-prime
  x = 11124440021748127159092076861405454814981575144744508857178576572929321435002942998531420985771090167262256877805902135304112271641074498386662361391760451
  y = 11124440021748127159092076861405454814981575144744508857178576572929321435002942998531420985771090167262256877805902135304112271641074498386662361391661439
  n = x - y
  l = [pow(a, -1, x) for a in range(y+1,x-1)]
  prod = 1
  for i in l:
      prod = (prod * i) % x
  p = sympy.nextprime(prod)
  return p

def get_q():
  value = [80096058210213458444437404275177554701604739094679033012396452382975889905967]
  for i in range (1, 10):
    value.insert(0,sympy.prevprime(value[0]))
  n = 1
  # φ(n)
  t = 1
  for i in value:
    n = n * i
    t = t * (i - 1)
  # e^-1 (mod φ(n))
  e_inv = pow(e, -1, t)
  q0 = pow(value_q, e_inv, n)
  q = sympy.nextprime(q0)
  return q

def main():
  p = get_p()
  q = get_q()
  n = p * q
  # φ(n)
  t = (p - 1) * (q - 1)
  # e^-1 (mod φ(n))
  d = pow(e, -1, t)
  m = pow(c, d, n)
  print(m.to_bytes((m.bit_length() + 7) // 8, byteorder='big').decode("ascii"))

if __name__ == '__main__':
  main()

执行代码即可打印出 flag