4-使用scipy进行积分和图片处理

tags: python,数据分析,pandas,numpy,Series,DataFrame,matplotlib,pyplot,scipy

目录

---

一、使用数值积分，求解圆周率

思路：使用scipy中的integrate对函数(1 - x^2)^0.5进行积分。

1.1. 导入积分需要使用的模块

from scipy import integrate

1.2. 定义要积分的函数表达式

def circle_half(x):
    return (1 - x**2)**0.5

或者使用匿名函数：circle_half = lambda x: (1 - x**2)**0.5

1.3. 求取积分

half_pi, deviation = integrate.quad(circle_half, -1, 1)
其中上面的代码中，积分函数要传入待积分函数名、积分上下限，最后返回的结果是个元组，元组中第一个为积分值，第二个为误差。
上面求到的积分值为圆周率的一半，最后再乘以2即可得到圆周率。

二、使用傅立叶变换实现登月图片的消噪

一般图片使用傅立叶变换实现消澡的具体步骤：

1. 第一步:：读取图片数据，使用方法把彩色图片变成黑白图片
2. 第二步：使用傅里叶变换，把图片数据变成频率域数据
3. 第三步：把频率域数据中的高频点(噪点) 赋值为0
4. 第四步：进行傅里叶逆变换
5. 第五步：去虚保实

对图片的消噪处理一般针对于黑白图片，所以需要将图片的三维数据转化为二维数据，即将其中的色彩维度去掉。

2.1. 读取登月图片数据

先导入数据处理的常用模块包：

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt

使用plt模块中的imread方法读取图片数据：
moon = plt.imread('./data/moonlanding.png')
由于该图片本身就是黑白图片，所以不需要进行灰度处理。

注意使用imread()方法读取的数据是numpy.ndarray类型。

2.2. 傅立叶变换

导入傅立叶变换需要的包：
from scipy.fftpack import fft2
使用该模块将上面的图片数据转换为频域数据：
moon_fft = fft2(moon)

scipy.fftpack模块用来计算快速傅里叶变换
速度比传统傅里叶变换更快，是对之前算法的改进
图片是二维数据，注意使用fftpack的二维转变方法

2.3. 将频域中的高频点置0

模值大于某一频率(阈值)的点我们认为是高频点即噪点，取到该点，将其置0。
moon_fft[np.abs(moon_fft) > 8e2] = 0
还有另一中方法取出高频点置0，代码如下：
result = np.where(np.abs(moon_fft) > 8e3, 0, moon_fft)

上面判断高频点的阈值大小，需要根据经验设置，需要不断调整到合适的值，没有确定值。

2.4. 傅里叶逆变换

导入傅立叶逆变换需要的包：
from scipy.fftpack import ifft2
使用该模块将上面的频域数据进行傅里叶逆变换：
moon_ifft = ifft2(moon_fft)
但是逆变换之后的数据中，还是存在虚数部分，所以需要将其中 的虚数部分去掉，即需要进行去虚保实

2.5. 去虚保实

代码如下：
moon_cleaned = np.real(moon_ifft)
执行上面的操作之后得到的数据moon_cleaned就是最后去噪后的数据，可以将该数据和原始数据画图比较，如果效果不好，则需要重新设置阈值，然后再继续操作，直到得到满意的结果。

三、使用scipy进行图片处理

3.1. 导入需要使用的模块，读取图片

from scipy import ndimage
face = plt.imread('./data/cat.jpg')
plt.figure(figsize=(12,9))
plt.imshow(face)

3.2. 移动图片

face_shift = ndimage.shift(face,(200,300,0))
plt.figure(figsize=(10,8))
plt.imshow(face_shift)

其中的.shift(face,(200,300,0))中需要先传入待处理的图片数据，然后传入移动的坐标，上面的代码表明，图片以左上角为原点，下移200像素，右移300像素，颜色不变。

3.3. 旋转图片

face_rotate = ndimage.rotate(face,angle=90)
plt.imshow(face_rotate)

3.4. zoom缩放图片

face_zoom = ndimage.zoom(face,zoom=(0.5,0.8,1))
plt.imshow(face_zoom)

注意，如果图片是黑白的，则zoom中的数据也要给包含两个数据的元组。

3.5. 切割图片

face_new = face[50:400,100:640]
plt.imshow(face_new)

注意：数据切片的时候，先切片的是行数据，然后切片的是列数据。

3.6. 彩色图片变成黑白

face_gray = np.dot(face,[0.299,0.587,0.114])
plt.figure(figsize=(10,8))
plt.imshow(face_gray,cmap="gray")

3.7. 图片增加随机噪声

这里我们只针对黑白图片，我们用上面的黑白图片数据进行处理

face_noise = face_gray.copy()
face_noise = face_noise.astype('float64')
face_noise += np.random.randn(*face_noise.shape) * face_noise.std() * 0.5

首先将黑白图片数据复制一份，由于后面需要进行一定的计算，这里将其转化为float64类型的数据。然后使用函数：np.random.randn(*face_noise.shape)创建保准正态分布（平均值为0，方差为1的正态分布），由于该分布需要和图片的维度、尺寸一致，所以使用.shape的方式，取出图片的维度和尺寸，由于该属性返回的数据是元组，所以需要使用*来进行解包。将分布数据和原数据的方差相乘然后取一半作为噪声数据（取一半是因为不取的话，可能噪声太大）
使用以下代码查看原数据和增加随机噪声之后的数据对应的图片：

plt.figure(figsize=(8*2,5))
axes1 = plt.subplot(1,2,1)
axes1.imshow(face_gray,cmap="gray")
axes2 = plt.subplot(1,2,2)
axes2.imshow(face_noise,cmap="gray")

3.8. 使用高斯滤波（参数sigma）对噪声图片进行处理

face_gaussian = ndimage.gaussian_filter(face_noise, sigma=1)

上图中的sigma值需要自己不断的调节，以找到最好的效果
使用下面的代码查看原数据、加噪数据和消噪数据对应的图片：

plt.rcParams['axes.unicode_minus'] = False
plt.rcParams['font.sans-serif'] = ['SimHei']
plt.figure(figsize=(6*3,5))
axes1 = plt.subplot(1,3,1,title='原始数据')
axes1.imshow(face_gray,cmap="gray")
axes2 = plt.subplot(1,3,2,title='加噪数据')
axes2.imshow(face_noise,cmap="gray")
axes3 = plt.subplot(1,3,3,title='消噪数据')
axes3.imshow(face_gaussian,cmap="gray")

3.9. 使用中值滤波（参数size）对噪声图片进行处理

中值滤波参数size：给出在每个元素上从输入数组中取出的形状位置，定义过滤器功能的输入，需要根据情况手动不断调整。

face_median = ndimage.median_filter(face_noise, size=5)

3.10. 使用维纳滤波（参数mysize）对噪声图片进行处理

signal维纳滤波mysize：滤镜尺寸的标量。

from scipy import signal
face_wiener = signal.wiener(face_noise, mysize=7)