import argparse import logging import re from multiprocessing import Process, Queue from pathlib import Path import numpy as np from skimage import exposure, filters from modules.config import logger from modules.volume import volume_loading_func, volume_saving_func def normalize_intensity( np_volume: np.ndarray, relative_path: Path, logger: logging.Logger ): logger.info(f"[processing start] {relative_path}") nstack = len(np_volume) stack: np.ndarray = np_volume[nstack // 2 - 16 : nstack // 2 + 16] hist_y, hist_x = exposure.histogram(stack[stack > 0]) thr = filters.threshold_otsu(stack[stack > 0]) peak_air = np.argmax(hist_y[hist_x < thr]) + hist_x[0] peak_soil = np.argmax(hist_y[hist_x > thr]) + (thr - hist_x[0]) + hist_x[0] np_volume = np_volume.astype(np.int64) for i in range(len(np_volume)): np_volume[i] = ( (np_volume[i] - peak_air).clip(0) / (peak_soil - peak_air) * 256 / 2 ) logger.info(f"[processing end] {relative_path}") return exposure.rescale_intensity( np_volume, in_range=(0, 255), out_range=(0, 255) ).astype(np.uint8) 请详细解释每一行的代码意思
时间: 2023-06-02 22:04:39 浏览: 109
这是一段 Python 代码,主要是实现了对 3D 数据进行灰度值归一化预处理的功能,下面逐行解释:
```python
import argparse
import logging
import re
from multiprocessing import Process, Queue
from pathlib import Path
import numpy as np
from skimage import exposure, filters
from modules.config import logger
from modules.volume import volume_loading_func, volume_saving_func
```
首先是导入一些 Python 模块和工具。其中 argparse、logging、re、Path、numpy、skimage 等模块都是 Python 内置模块,而 multiprocessing 和 skimage 的 exposure、filters 模块需要额外安装。另外,这段代码还从自定义模块 modules 中导入了两个函数,分别是 volume_loading_func 和 volume_saving_func。
```python
def normalize_intensity( np_volume: np.ndarray, relative_path: Path, logger: logging.Logger ):
```
定义了一个名为 normalize_intensity 的函数,它接受三个参数:np_volume,表示一个 numpy 数组,存储了 3D 数据的原始灰度值;relative_path,表示 3D 数据的存储路径,类型为 Path;logger,表示一个 logging.Logger 对象,用于写日志。
```python
logger.info(f"[processing start] {relative_path}")
```
使用 logger 对象写日志,输出一个带有相对路径的日志信息,表示灰度归一化开始。
```python
nstack = len(np_volume)
stack: np.ndarray = np_volume[nstack // 2 - 16 : nstack // 2 + 16]
```
获取 np_volume 数组的长度 nstack,并将 np_volume 数组中的中央一部分数据赋值给 stack 变量,这里只保留了中央的 32 张切片。
```python
hist_y, hist_x = exposure.histogram(stack[stack > 0])
thr = filters.threshold_otsu(stack[stack > 0])
peak_air = np.argmax(hist_y[hist_x < thr])
peak_soil = np.argmax(hist_y[hist_x > thr])
```
利用 skimage 的 exposure 模块,对 stack 中大于 0 的像素值进行直方图统计,并使用 filters 模块的 threshold_otsu 函数获取阈值 thr。然后,寻找 hist_y 中小于 thr 的最大值所在的下标 peak_air,和 hist_y 中大于 thr 的最大值所在的下标 peak_soil。
```python
np_volume = np_volume.astype(np.int64)
for i in range(len(np_volume)):
np_volume[i] = ((np_volume[i] - peak_air).clip(0) / (peak_soil - peak_air) * 256 / 2)
```
将 np_volume 数组转化为 int64 类型,并对其中的每一个像素值进行灰度归一化处理。首先将每个像素值减去 peak_air,然后将结果限制在 [0, max_soil - peak_air] 范围内,最后将结果缩放到 [0, 128] 范围内。
```python
logger.info(f"[processing end] {relative_path}")
```
使用 logger 对象写日志,输出一个带有相对路径的日志信息,表示灰度归一化处理结束。
```python
return exposure.rescale_intensity(
np_volume, in_range=(0, 255), out_range=(0, 255)
).astype(np.uint8)
```
将归一化的 numpy 数组 np_volume 传入 skimage 的 exposure 模块的 rescale_intensity 函数,将灰度值缩放到 [0, 255] 范围内,并将结果转化为 uint8 数据类型,最终返回归一化处理后的 3D 数据。
阅读全文