NumPy 音頻和圖像處理

# 來源：NumPy Cookbook 2e Ch5

將圖像加載進(jìn)內(nèi)存

import numpy as np 
import matplotlib.pyplot as plt

# 首先生成一個(gè) 512x512 的圖像
# 在里面畫 30 個(gè)正方形
N = 512 
NSQUARES = 30

# 初始化
img = np.zeros((N, N), np.uint8) 
# 正方形的中心是 0 ~ N 的隨機(jī)數(shù)
centers = np.random.random_integers(0, N, size=(NSQUARES, 2))
# 正方形的邊長是 0 ~ N/9 的隨機(jī)數(shù)
radii = np.random.randint(0, N/9, size=NSQUARES) 
# 顏色是 100 ~ 255 的隨機(jī)數(shù)
colors = np.random.randint(100, 255, size=NSQUARES)

# 生成正方形
for i in xrange(NSQUARES):
    # 為每個(gè)正方形生成 x 和 y 坐標(biāo)
    xindices = range(centers[i][0] - radii[i], centers[i][0]  + radii[i])   
    xindices = np.clip(xindices, 0, N - 1)   
    yindices = range(centers[i][1] - radii[i], centers[i][1]  + radii[i])   
    
    # clip 過濾范圍之外的值
    # 相當(dāng)于 yindices = yindices[(0 < yindices) & (yindices < N - 1)]
    yindices = np.clip(yindices, 0, N - 1)
    if len(xindices) == 0 or len(yindices) == 0:
        continue
    # 將 x 和 y 坐標(biāo)轉(zhuǎn)換成網(wǎng)格
    # 如果不轉(zhuǎn)換成網(wǎng)格姿鸿，只會給對角線著色
    coordinates = np.meshgrid(xindices, yindices)     
    img[coordinates] = colors[i]
   
# tofile 以二進(jìn)制保存數(shù)組的內(nèi)容，沒有形狀和類型信息 
img.tofile('random_squares.raw') 
# np.memmap 以二進(jìn)制加載數(shù)組贞言，如果類型不是 uint8冠胯，則需要執(zhí)行
# 如果數(shù)組不是一維，還需要指定形狀
img_memmap = np.memmap('random_squares.raw', shape=img.shape)

# 顯示圖像（會自動(dòng)將灰度圖映射為偽彩色）
plt.imshow(img_memmap) 
plt.axis('off') 
plt.show()

組合圖像

import numpy as np import 
matplotlib.pyplot as plt 
from scipy.misc import lena

ITERATIONS = 10 
lena = lena() 
SIZE = lena.shape[0] 
MAX_COLOR = 255. 
x_min, x_max = -2.5, 1 
y_min, y_max = -1, 1

# 數(shù)組初始化
x, y = np.meshgrid(np.linspace(x_min, x_max, SIZE),
                   np.linspace(y_min, y_max, SIZE)) 
c = x + 1j * y 
z = c.copy() 
fractal = np.zeros(z.shape, dtype=np.uint8) + MAX_COLOR 

# 生成 mandelbrot 圖像 
for n in range(ITERATIONS):
    mask = np.abs(z) <= 4
    z[mask] = z[mask] ** 2 +  c[mask]
    fractal[(fractal == MAX_COLOR) & (-mask)] = (MAX_COLOR - 1) * n / ITERATIONS

# 繪制 mandelbrot 圖像 
plt.subplot(211) 
plt.imshow(fractal) 
plt.title('Mandelbrot') 
plt.axis('off')

# 將 mandelbrot 和 lena 組合起來
plt.subplot(212) 
# choose 的作用是躲查，如果 fractal 的元素小于 lena 的對應(yīng)元素
# 就選擇 fractal它浅，否則選擇 lena
# 相當(dāng)于 np.fmin(fractal, lena)
plt.imshow(np.choose(fractal < lena, [fractal, lena])) 
plt.axis('off') 
plt.title('Mandelbrot + Lena')
plt.show()

使圖像變模糊

import numpy as np 
import matplotlib.pyplot as plt 
from random import choice 
import scipy 
import scipy.ndimage

# Initialization 
NFIGURES = 5 
k = np.random.random_integers(1, 5, NFIGURES) 
a = np.random.random_integers(1, 5, NFIGURES)
colors = ['b', 'g', 'r', 'c', 'm', 'y', 'k']

# 繪制原始 的 lena 圖像
lena = scipy.misc.lena() 
plt.subplot(211) 
plt.imshow(lena) 
plt.axis('off')

# 繪制模糊的 lena 圖像
plt.subplot(212) 
# 使用 sigma=4 的高斯過濾器
blurred = scipy.ndimage.gaussian_filter(lena, sigma=4)
plt.imshow(blurred) 
plt.axis('off')

# 在極坐標(biāo)中繪圖
# 極坐標(biāo)無視 subplot
theta = np.linspace(0, k[0] * np.pi, 200) 
plt.polar(theta, np.sqrt(theta), choice(colors))

for i in xrange(1, NFIGURES):
    theta = np.linspace(0, k[i] * np.pi, 200)   
    plt.polar(theta, a[i] * np.cos(k[i] * theta), choice(colors))
plt.axis('off')
plt.show()

重復(fù)聲音片段

import scipy.io.wavfile 
import matplotlib.pyplot as plt 
import urllib2 
import numpy as np

# 下載音頻文件
response = urllib2.urlopen('http://www.thesoundarchive.com/ austinpowers/smashingbaby.wav') 
print(response.info()) 

# 將文件寫到磁盤
WAV_FILE = 'smashingbaby.wav' 
filehandle = open(WAV_FILE, 'w') 
filehandle.write(response.read()) 
filehandle.close() 

# 使用 SciPy 讀取音頻文件
sample_rate, data = scipy.io.wavfile.read(WAV_FILE) 
print("Data type", data.dtype, "Shape", data.shape)
# ('Data type', dtype('uint8'), 'Shape', (43584L,))

# 繪制原始音頻文件
plt.subplot(2, 1, 1)
plt.title("Original") 
plt.plot(data)

# 繪制重復(fù)后的音頻文件
plt.subplot(2, 1, 2)
# tile 用于重復(fù)數(shù)組
repeated = np.tile(data, 3)
plt.title("Repeated") 
plt.plot(repeated) 

# 保存重復(fù)后的音頻文件
scipy.io.wavfile.write("repeated_yababy.wav", sample_rate, repeated)
plt.show()

生成聲音

# 聲音可以表示為某個(gè)振幅、頻率和初相的正弦波
# 如果我們把鋼琴上的鍵編為 1 ~ 88镣煮，
# 那么它的頻率就是 440 * 2 ** ((n - 49) / 12)
# 其中 n 是鍵的編號

import scipy.io.wavfile 
import numpy as np
import matplotlib.pyplot as plt

RATE = 44100 
DTYPE = np.int16

# 生成正弦波 
def generate(freq, amp, duration, phi): 
    t = np.linspace(0, duration, duration * RATE) 
    data = np.sin(2 * np.pi * freq * t + phi) * amp
    
    return data.astype(DTYPE)

# 初始化
# 彈奏 89 個(gè)音符
NTONES = 89 
# 振幅是 200 ~ 2000
amps = 2000. * np.random.random((NTONES,)) + 200. 
# 時(shí)長是 0.01 ~ 0.2
durations = 0.19 * np.random.random((NTONES,)) + 0.01 
# 鍵從 88 個(gè)中任取
keys = np.random.random_integers(1, 88, NTONES) 
# 頻率使用上面的公式生成
freqs = 440.0 * 2 ** ((keys - 49.)/12.) 
# 初相是 0 ~ 2 * pi
phi = 2 * np.pi * np.random.random((NTONES,))

tone = np.array([], dtype=DTYPE)

for i in xrange(NTONES):   
    # 對于每個(gè)音符生成正弦波
    newtone = generate(freqs[i], amp=amps[i],  duration=durations[i], phi=phi[i])   
    # 附加到音頻后面
    tone = np.concatenate((tone, newtone))

# 保存文件
scipy.io.wavfile.write('generated_tone.wav', RATE, tone)

# 繪制音頻數(shù)據(jù)
plt.plot(np.linspace(0, len(tone)/RATE, len(tone)), tone) 
plt.show()

設(shè)計(jì)音頻濾波器

import scipy.io.wavfile 
import matplotlib.pyplot as plt 
import urllib2 
import numpy as np

# 下載音頻文件
response = urllib2.urlopen('http://www.thesoundarchive.com/ austinpowers/smashingbaby.wav') 
print(response.info()) 

# 將文件寫到磁盤
WAV_FILE = 'smashingbaby.wav' 
filehandle = open(WAV_FILE, 'w') 
filehandle.write(response.read()) 
filehandle.close() 

# 使用 SciPy 讀取音頻文件
sample_rate, data = scipy.io.wavfile.read(WAV_FILE) 
print("Data type", data.dtype, "Shape", data.shape)
# ('Data type', dtype('uint8'), 'Shape', (43584L,))

# 繪制原始音頻文件
plt.subplot(2, 1, 1)
plt.title("Original") 
plt.plot(data)

# 設(shè)計(jì)濾波器姐霍，iirdesign 設(shè)計(jì)無限脈沖響應(yīng)濾波器
# 參數(shù)依次是 0 ~ 1 的正則化頻率、
# 最大損失典唇、最低衰減和濾波類型
b,a = scipy.signal.iirdesign(wp=0.2, ws=0.1, gstop=60, gpass=1, ftype='butter')

# 傳入剛才的返回值镊折，使用 lfilter 函數(shù)來調(diào)用濾波器
filtered = scipy.signal.lfilter(b, a, data)

# 繪制濾波后的音頻
plt.subplot(2, 1, 2) 
plt.title("Filtered") 
plt.plot(filtered)

# 保存濾波后的音頻
scipy.io.wavfile.write('filtered.wav', sample_rate, filtered. astype(data.dtype))
plt.show()

Sobel 過濾器的邊界檢測

# Sobel 過濾器用于提取圖像的邊界
# 也就是將圖像轉(zhuǎn)換成線框圖風(fēng)格
import scipy 
import scipy.ndimage 
import matplotlib.pyplot as plt

# 導(dǎo)入 Lena
lena = scipy.misc.lena()

# 繪制 Lena（左上方）
plt.subplot(221) 
plt.imshow(lena) 
plt.title('Original') 
plt.axis('off')


# Sobel X 過濾器過濾后的圖像（右上方）
sobelx = scipy.ndimage.sobel(lena, axis=0, mode='constant')
plt.subplot(222) 
plt.imshow(sobelx) 
plt.title('Sobel X') 
plt.axis('off')

# Sobel Y 過濾器過濾的圖像（左下方） 
sobely = scipy.ndimage.sobel(lena, axis=1, mode='constant')
plt.subplot(223) 
plt.imshow(sobely) 
plt.title('Sobel Y') 
plt.axis('off')

# 默認(rèn)的 Sobel 過濾器（右下方）
default = scipy.ndimage.sobel(lena)
plt.subplot(224) 
plt.imshow(default) 
plt.title('Default Filter') 
plt.axis('off')
plt.show()