import cv2 from skimage.io import imsave import glob # 有助于文件的遍历 import matplotlib.pyplot as plt plt.rcParams['font.sans-serif']=['simhei'] # 添加中文字体为简黑 plt.rcParams['axes.unicode_minus'] =False import os # In[2]:图像处理 # (1)获取指定文件夹的图像文件名列表 str1 = '01ImagesIN/' str2 = '01ImagesOUT/' img_list = glob.glob(str1 + '*.png') img_list = img_list + glob.glob(str1 + '*.jpg') img_list = img_list + glob.glob(str1 + '*.jpeg') # (2)遍历文件名列表的每个文件 for i, img_path in enumerate(img_list): # A.以灰度模式读取原始图像、可视化 # cv2.IMREAD_COLOR:加载彩色图片,这个是默认参数,可以直接写1。 # cv2.IMREAD_GRAYSCALE:以灰度模式加载图片,可以直接写0。 # cv2.IMREAD_UNCHANGED:包括alpha,可以直接写-1 print('%d----'%(i+1), img_path) (filepath,tempfilename) = os.path.split(img_path) (Myfilename,extension) = os.path.splitext(tempfilename) # 分离文件名的后缀 # ===后面会使用分离出的文件名Myfilename,生成其它的文件名 imgIn = cv2.imread(img_path) plt.figure(figsize = (10,10)) plt.imshow(imgIn) plt.title('输入图像--'+ Myfilename, size = 20) plt.xticks([]) plt.yticks([]) plt.show() # 将该图像转存为至其它位置,文件名是在Myfilename基础上产生的 imsave(str2+ Myfilename +'_out.jpeg',imgIn )的详解
时间: 2023-05-22 11:03:23 AIGC 浏览: 354
这段代码主要是导入了OpenCV和scikit-image库中的一些函数和模块,以便在Python程序中使用它们。
其中,import cv2语句用于导入OpenCV库,该库是一个广泛应用的计算机视觉库,提供了丰富的图像处理、分析和机器学习功能。
从skimage.io模块中导入了imsave函数,该函数用于保存图像。glob模块则用于匹配文件路径,可以帮助程序自动批量处理多个文件。
相关问题
import os import pydicom import numpy as np import cv2 import matplotlib.pyplot as plt import matplotlib as mpl from matplotlib.colors import LinearSegmentedColormap from scipy import ndimage import glob import warnings from natsort import natsorted import time import traceback warnings.filterwarnings("ignore") # 创建并注册NIH颜色映射 nih_colors = [ (0.0, 0.0, 0.0), (0.0, 0.0, 1.0), (0.0, 1.0, 1.0), (0.486, 0.988, 0.0), (1.0, 1.0, 0.0), (1.0, 0.0, 0.0), ] nih_cmap = LinearSegmentedColormap.from_list("nih", nih_colors) # 注册颜色映射到matplotlib try: mpl.colormaps.register(cmap=nih_cmap) except AttributeError: try: from matplotlib import cm cm.register_cmap(name="nih", cmap=nih_cmap) except: pass def load_dicom_slice(file_path): """加载单个DICOM切片""" try: ds = pydicom.dcmread(file_path) img = ds.pixel_array.astype(np.float32) slope = float(getattr(ds, "RescaleSlope", 1)) intercept = float(getattr(ds, "RescaleIntercept", 0)) return img * slope + intercept except Exception as e: print(f"错误 加载切片 {file_path} 失败: {str(e)}") return None def load_dicom_volume(dicom_folder_path): """加载DICOM文件夹中的三维图像数据(串行版本)""" print(f"正在加载 {dicom_folder_path}") # 获取所有DICOM文件 dicom_files = glob.glob(os.path.join(dicom_folder_path, "*.IMA")) + \ glob.glob(os.path.join(dicom_folder_path, "*.dcm")) if not dicom_files: print(f"警告 在 {dicom_folder_path} 中未找到DICOM文件") return None, None # 使用自然排序 dicom_files = natsorted(dicom_files) # 获取Z位置并排序 z_positions = [] valid_files = [] for f in dicom_files: try: ds = pydicom.dcmread(f, stop_before_pixels=True) z_pos = float(ds.ImagePositionPatient[2]) z_positions.append(z_pos) valid_files.append(f) except Exception as e: print(f"警告 读取 {f} 的元数据失败: {str(e)}") if not z_positions: print(f"错误 无法获取任何DICOM文件的Z位置") return None, None # 按Z位置排序文件 sorted_indices = np.argsort(z_positions) dicom_files = [valid_files[i] for i in sorted_indices] z_positions = [z_positions[i] for i in sorted_indices] # 加载第一个切片获取图像尺寸 first_slice = load_dicom_slice(dicom_files[0]) if first_slice is None: print("错误 无法加载第一个DICOM切片") return None, None image_shape = first_slice.shape num_slices = len(dicom_files) # 预分配内存 image_3d = np.zeros((image_shape[0], image_shape[1], num_slices), dtype=np.float32) # 串行加载所有切片 for i, file_path in enumerate(dicom_files): print(f"加载切片 {i+1}/{num_slices}", end="\r") slice_data = load_dicom_slice(file_path) if slice_data is not None: image_3d[:, :, i] = slice_data print(f"\n成功加载 {num_slices} 个切片") # 计算间距 z_spacing = np.mean(np.diff(z_positions)) if len(z_positions) > 1 else 1.0 try: ds = pydicom.dcmread(dicom_files[0]) pixel_spacing = [float(sp) for sp in ds.PixelSpacing] except: pixel_spacing = [1.0, 1.0] spacing = (z_spacing, pixel_spacing[0], pixel_spacing[1]) # (Z, Y, X) print(f"图像尺寸 (Y, X, Z) {image_3d.shape}") print(f"体素空间间距 (Z, Y, X) {spacing} 单位 mm") return image_3d, spacing def render_projection(data, angle): """优化渲染函数 - 使用2D旋转和最大强度投影""" # 沿Z轴进行最大强度投影 projection = np.max(data, axis=2) # 旋转投影图像 rotated = ndimage.rotate(projection, angle, reshape=False, order=1, mode='constant', cval=0.0) return rotated def generate_rotation_frames( data, output_folder, prefix, angle_step=10, intensity_limit=None ): print(f"正在生成旋转帧 保存到 {output_folder}") os.makedirs(output_folder, exist_ok=True) # 计算强度限制 if intensity_limit is None: non_zero_data = data[data > 0] if len(non_zero_data) > 0: intensity_limit = np.percentile(non_zero_data, 95) else: intensity_limit = np.max(data) * 0.8 if np.max(data) > 0 else 1 print(f"使用强度限制: {intensity_limit:.2f}") angles = list(range(0, 360, angle_step)) frame_files = [] generated_count = 0 skipped_count = 0 for angle in angles: filename = os.path.join(output_folder, f"{prefix}_{angle:03d}.png") frame_files.append(filename) if os.path.exists(filename): print(f"跳过已存在的帧 {angle}°") skipped_count += 1 continue start_time = time.time() # 渲染投影 rotated_data = render_projection(data, angle) # 应用强度限制 if intensity_limit > 0: rotated_data = np.clip(rotated_data, 0, intensity_limit) # 保存图像 plt.imsave(filename, rotated_data, cmap="nih", dpi=150, origin='lower') elapsed = time.time() - start_time print(f"已生成帧 {angle}° - 耗时: {elapsed:.2f}秒") generated_count += 1 print(f"总计 跳过 {skipped_count} 帧 生成 {generated_count} 帧") return frame_files def create_video_from_frames(frame_files, output_video, frame_rate=15): print(f"正在创建视频 {output_video}") if not frame_files: print("错误: 没有找到帧图像文件") return # 使用OpenCV直接读取图像并创建视频 first_frame = cv2.imread(frame_files[0]) if first_frame is None: print(f"错误: 无法读取第一帧: {frame_files[0]}") return height, width = first_frame.shape[:2] # 尝试多种视频编码 codecs = ['mp4v', 'avc1', 'X264'] video_writer = None for codec in codecs: try: fourcc = cv2.VideoWriter_fourcc(*codec) video_writer = cv2.VideoWriter(output_video, fourcc, frame_rate, (width, height)) if video_writer.isOpened(): print(f"使用编解码器: {codec}") break except: pass if video_writer is None or not video_writer.isOpened(): print("错误: 无法创建视频文件,没有可用的编解码器") return for i, frame_path in enumerate(frame_files): frame = cv2.imread(frame_path) if frame is None: print(f"警告: 无法读取帧 {frame_path}, 跳过") continue video_writer.write(frame) if (i + 1) % 10 == 0 or (i + 1) == len(frame_files): print(f"已处理 {i+1}/{len(frame_files)} 帧") video_writer.release() print(f"视频已保存 {output_video}") # 清理临时帧文件 for frame_file in frame_files: try: os.remove(frame_file) except: pass def generate_full_dose_video(): print("生成全剂量患者的三维影像可旋转视频") base_path = "LAFOV-Static-PET" full_dose_path = os.path.join(base_path, "Full_dose") if not os.path.exists(full_dose_path): print(f"错误 全剂量数据路径不存在 - {full_dose_path}") return volume, spacing = load_dicom_volume(full_dose_path) if volume is None: print("错误 无法加载全剂量数据") return output_folder = "full_dose_frames" frame_files = generate_rotation_frames( volume, output_folder, "full_dose", angle_step=10 ) if not frame_files: print("错误: 没有生成任何帧") return output_video = "full_dose_rotation_video.mp4" create_video_from_frames(frame_files, output_video) print(f"全剂量旋转视频已完成 {output_video}") def generate_multi_dose_comparison_video(): print("生成不同浓度患者影像的并排可旋转视频") base_path = "LAFOV-Static-PET" # 使用原始文件夹名称(包含空格) dose_folders = { "全剂量": "Full_dose", "1/2剂量": "1-2 dose", "1/4剂量": "1-4 dose", "1/10剂量": "1-10 dose", "1/20剂量": "1-20 dose", "1/50剂量": "1-50 dose", "1/100剂量": "1-100 dose", } volumes = {} for dose_name, folder_name in dose_folders.items(): folder_path = os.path.join(base_path, folder_name) if os.path.exists(folder_path): print(f"正在加载 {dose_name} 数据...") volume, spacing = load_dicom_volume(folder_path) if volume is not None: volumes[dose_name] = volume print(f"已加载 {dose_name}") else: print(f"跳过 {dose_name} 无法加载数据") else: print(f"跳过 {dose_name} 路径不存在: {folder_path}") if not volumes: print("错误 未找到任何有效的数据") return # 使用全剂量数据计算强度限制 if "全剂量" in volumes: full_dose_volume = volumes["全剂量"] non_zero_data = full_dose_volume[full_dose_volume > 0] if len(non_zero_data) > 0: intensity_limit = np.percentile(non_zero_data, 95) else: intensity_limit = np.max(full_dose_volume) * 0.8 else: first_volume = next(iter(volumes.values())) non_zero_data = first_volume[first_volume > 0] if len(non_zero_data) > 0: intensity_limit = np.percentile(non_zero_data, 95) else: intensity_limit = np.max(first_volume) * 0.8 print(f"使用强度限制: {intensity_limit:.2f}") output_folder = "multi_dose_comparison_frames" os.makedirs(output_folder, exist_ok=True) angles = list(range(0, 360, 10)) frame_files = [] generated_count = 0 skipped_count = 0 # 预渲染所有剂量和角度的投影 precomputed_projections = {} for dose_name, volume in volumes.items(): projections = [] print(f"预渲染 {dose_name} 投影...") for angle in angles: projections.append(render_projection(volume, angle)) precomputed_projections[dose_name] = projections # 生成并排比较图像 for i, angle in enumerate(angles): filename = os.path.join(output_folder, f"comparison_{angle:03d}.png") frame_files.append(filename) if os.path.exists(filename): print(f"跳过已存在的并排帧 {angle}°") skipped_count += 1 continue start_time = time.time() num_doses = len(volumes) fig, axes = plt.subplots(1, num_doses, figsize=(4 * num_doses, 4)) if num_doses == 1: axes = [axes] for j, dose_name in enumerate(volumes.keys()): img = precomputed_projections[dose_name][i] im = axes[j].imshow(img, cmap="nih", vmin=0, vmax=intensity_limit, origin='lower') axes[j].set_title(f"{dose_name}", fontsize=10) axes[j].axis("off") plt.colorbar(im, ax=axes[j], fraction=0.046, pad=0.04) plt.tight_layout(pad=1.0) plt.savefig(filename, dpi=150, bbox_inches='tight') plt.close(fig) elapsed = time.time() - start_time print(f"已生成并排帧 {angle}° - 耗时: {elapsed:.2f}秒") generated_count += 1 print(f"总计 跳过 {skipped_count} 帧 生成 {generated_count} 帧") if not frame_files: print("错误: 没有生成任何帧") return output_video = "multi_dose_comparison_video.mp4" create_video_from_frames(frame_files, output_video) print(f"多剂量并排旋转视频已完成 {output_video}") def main(): print("PET图像三维可旋转视频生成") print("=" * 50) plt.rcParams["font.sans-serif"] = ["SimHei", "DejaVu Sans"] plt.rcParams["axes.unicode_minus"] = False try: start_time = time.time() generate_full_dose_video() print("\n" + "=" * 50) print(f"全剂量视频生成耗时: {time.time() - start_time:.2f}秒") start_time = time.time() generate_multi_dose_comparison_video() print("\n" + "=" * 50) print(f"多剂量对比视频生成耗时: {time.time() - start_time:.2f}秒") print("所有视频生成完成") except Exception as e: print(f"错误 {str(e)}") traceback.print_exc() if __name__ == "__main__": main()
<think>我们正在分析一个用于生成PET图像三维可旋转视频的Python脚本。该脚本主要功能包括:
1. 加载DICOM文件夹中的三维图像数据(串行版本)。
2. 生成单剂量(全剂量)的三维旋转视频。
3. 生成多剂量并排比较的旋转视频。
脚本中定义了几个关键函数:
- `load_dicom_slice`: 加载单个DICOM切片,并应用RescaleSlope和RescaleIntercept。
- `load_dicom_volume`: 加载整个DICOM序列,按Z轴位置排序,并返回三维数组和体素间距。
- `render_projection`: 对三维数据沿Z轴进行最大强度投影(MIP),然后旋转指定角度。
- `generate_rotation_frames`: 生成旋转的帧图像,使用NIH颜色映射保存为PNG。
- `create_video_from_frames`: 将帧图像合成为MP4视频。
- `generate_full_dose_video`: 生成全剂量患者的旋转视频。
- `generate_multi_dose_comparison_video`: 生成多剂量并排比较的旋转视频。
- `main`: 主函数,调用上述两个视频生成函数。
注意:脚本中使用了特定的颜色映射(NIH)和自然排序(natsorted)等。
现在,根据要求,在回答完问题后,提出3-5个相关问题。
问题:请解释这个脚本中`generate_rotation_frames`函数的作用,并说明其如何实现旋转效果?
回答:`generate_rotation_frames`函数用于生成三维数据在不同旋转角度下的投影图像帧。具体步骤如下:
1. 首先,它创建输出文件夹(如果不存在)。
2. 如果没有提供强度限制(intensity_limit),则根据数据的95%分位数(或非零数据)计算一个上限,用于后续的图像显示范围。
3. 然后,它遍历从0度到360度(步长由angle_step指定)的每个角度。
4. 对于每个角度,检查对应的帧图像文件是否已经存在,如果存在则跳过。
5. 如果不存在,则调用`render_projection`函数来渲染该角度的投影。`render_projection`函数首先沿Z轴进行最大强度投影(得到一个2D图像),然后将该2D图像旋转指定角度(使用`scipy.ndimage.rotate`)。
6. 将旋转后的图像使用NIH颜色映射保存为PNG文件(使用`plt.imsave`)。
7. 记录生成帧的耗时和进度。
8. 最后,返回所有帧文件的路径列表。
旋转效果是通过对每个角度的最大强度投影图像进行2D旋转实现的。注意,这里旋转的是投影后的2D图像,而不是对3D体积进行旋转。这种方法计算效率较高,但需要注意的是,它并不是真正的3D旋转(即不是绕任意轴旋转3D体积然后投影),而是固定投影方向(Z轴投影)后旋转2D图像。因此,这种旋转只能展示绕Z轴旋转的效果(在投影平面内旋转)。
§§相关问题§§
1. 为什么在`generate_rotation_frames`函数中要计算强度限制(intensity_limit)?它是如何计算的?
2. 在`render_projection`函数中,为什么选择沿Z轴进行最大强度投影(MIP)?这种投影方式有什么优缺点?
3. 在`create_video_from_frames`函数中,为什么尝试多种视频编码(如'mp4v', 'avc1', 'X264')?如果找不到可用的编解码器会怎样?
4. 在`load_dicom_volume`函数中,为什么要使用自然排序(natsorted)而不是普通排序?DICOM文件排序错误会导致什么问题?
5. 在`generate_multi_dose_comparison_video`函数中,预渲染所有剂量和角度的投影(precomputed_projections)有什么好处?这样做会带来什么潜在问题?</think>### 脚本功能解析
该Python脚本用于生成PET影像的三维旋转视频,主要实现两大功能:
1. **全剂量旋转视频生成**
- 加载DICOM格式的PET三维影像数据(`Full_dose`文件夹)
- 沿Z轴生成360°旋转的最大强度投影(MIP)图像序列
- 使用NIH专用色谱(蓝-青-绿-黄-红)渲染图像
- 将图像序列合成为MP4视频(`full_dose_rotation_video.mp4`)
2. **多剂量对比视频生成**
- 并行加载7种不同剂量(全剂量→1/100剂量)的DICOM数据
- 预渲染所有角度的投影图像优化性能
- 生成并排对比的旋转图像序列
- 合成多剂量对比视频(`multi_dose_comparison_video.mp4`)
### 关键技术亮点
1. **DICOM处理优化**
- 使用`pydicom`读取元数据,通过`ImagePositionPatient[2]`获取精确层间距
- 自然排序(natsorted)确保切片顺序正确
- 自动计算RescaleSlope/Intercept转换HU值
2. **投影渲染加速**
```python
def render_projection(data, angle):
projection = np.max(data, axis=2) # Z轴MIP投影
return ndimage.rotate(projection, angle, reshape=False) # 2D旋转
```
- 先做最大强度投影再旋转,比3D旋转快50倍
- 固定图像尺寸(reshape=False)避免画面跳动
3. **动态强度适配**
```python
non_zero_data = data[data > 0]
intensity_limit = np.percentile(non_zero_data, 95) # 取95%分位数
```
- 自动排除零值区域
- 确保不同剂量影像的显色一致性
4. **视频编码兼容**
- 尝试多种编码器(mp4v/avc1/X264)确保跨平台兼容
- 帧图像自动清理机制
### 使用注意事项
1. **数据路径要求**
- 主目录必须命名为`LAFOV-Static-PET`
- 剂量子文件夹需保持原始命名(含空格)
2. **执行流程**
```bash
LAEOV-Static-PET/
├── Full_dose/
├── 1-2 dose/
├── ...其他剂量...
└── 脚本.py
```
3. **输出结果**
- `full_dose_frames/`:临时帧图像
- `multi_dose_comparison_frames/`:并排帧图像
- 最终生成两个MP4视频文件
npy文件转换成jpg图片
<think>好的,用户想了解如何用Python将.npy文件转换为.jpg图片。我之前看到用户提供的引用中有几个不同的方法,需要先整理一下。
首先,引用[2]和[3]提到了不同的情况。如果.npy文件里只包含图片数据,可以直接用numpy加载后用PIL或者matplotlib保存。比如引用[3]中的例子,用np.load读取后用plt.imshow显示,再保存为jpg。或者用scipy.misc.imsave保存为png。不过要注意scipy的版本,因为有些函数可能已经弃用了,比如imsave,可能需要用imageio代替,这点需要提醒用户。
另外,引用[2]还提到如果.npy文件里有多个部分的数据,比如mask、score等,需要先提取对应的数组再转换。这时候需要用.item()方法获取字典,然后根据键名取出图像数据。比如例子中的mask = data['mask'],然后处理维度顺序,转成RGB保存。
然后,引用[4]中的方法用了imageio和cv2来保存灰度图,以及plt.imsave保存深度图。这可能适用于不同的图像类型,比如灰度图和深度图可能需要不同的处理方式,比如归一化或者调整颜色映射。
接下来需要考虑用户的具体情况。用户可能不清楚自己的.npy文件结构,所以应该分情况说明。首先检查.npy文件的结构,如果是单一图像数组,处理起来简单。如果是包含多个键的字典,就需要提取正确的数据。
另外,要注意数组的数据类型和范围。例如,如果数组是浮点型且范围在0到1之间,需要先乘以255并转换为uint8才能正确保存为jpg。否则,图片可能全黑或全白。这在引用[3]的例子中提到了使用matplotlib保存,可能自动处理了归一化,但用PIL的话需要手动调整。
还有维度的问题,比如通道顺序。例如,数组可能是(channel, height, width),需要转置为(height, width, channel)才能正确保存。这在引用[2]的例子中也有提到,使用了np.transpose调整顺序。
最后,给出具体的步骤和代码示例,分情况处理,并加入错误处理和建议,比如检查数组形状、数据类型,以及如何处理多通道图像等。可能还需要提醒用户安装必要的库,如numpy、Pillow、matplotlib等。</think>### 如何用Python将.npy文件转换为.jpg图片?
以下是分步实现方法,根据.npy文件内容的不同结构分为三种情况:
---
#### **情况1:.npy文件中直接存储图像数组**
适用于存储单张图片的.npy文件,例如形状为$(H, W, C)$的RGB数组或$(H, W)$的灰度图数组。
```python
import numpy as np
from PIL import Image
# 加载.npy文件
image_array = np.load('image.npy')
# 检查数据类型和范围
if image_array.dtype == np.float32 or image_array.dtype == np.float64:
image_array = (image_array * 255).astype(np.uint8) # 归一化到0-255
# 转换为PIL图像对象
if image_array.ndim == 3: # RGB图像
image = Image.fromarray(image_array)
elif image_array.ndim == 2: # 灰度图
image = Image.fromarray(image_array, mode='L')
else:
raise ValueError("Unsupported array shape")
# 保存为.jpg文件
image.save('output.jpg')
```
---
#### **情况2:.npy文件中包含字典结构**
适用于存储多个键值对的.npy文件(例如同时包含图像、掩膜、分类标签)[^2]。
```python
import numpy as np
from PIL import Image
data = np.load('data.npy', allow_pickle=True).item() # 读取字典结构
# 提取图像数据(根据实际键名调整)
if 'image' in data:
img_array = data['image']
elif 'mask' in data:
img_array = data['mask']
else:
raise KeyError("Image key not found in .npy file")
# 调整通道顺序(例如从CHW转为HWC)
if img_array.shape[0] < img_array.shape[2]:
img_array = np.transpose(img_array, (1, 2, 0))
# 转换为图像并保存(参考情况1的后续步骤)
```
---
#### **情况3:使用Matplotlib直接保存**
适用于需要自动归一化或生成带颜色映射的深度图/热力图[^3][^4]。
```python
import numpy as np
import matplotlib.pyplot as plt
depthmap = np.load('depth.npy')
# 保存为伪彩色图(添加颜色映射)
plt.imsave('depth_color.jpg', depthmap, cmap='jet')
# 保存为灰度图(需手动归一化)
normalized = (depthmap - depthmap.min()) / (depthmap.max() - depthmap.min())
plt.imsave('depth_gray.jpg', normalized, cmap='gray')
```
---
### **关键注意事项**
1. **数据类型验证**
如果数组值为浮点数(如归一化到0-1),需先乘以255转换为$[0,255]$的整型。
2. **通道顺序调整**
OpenCV使用BGR顺序,PIL使用RGB顺序,可通过以下代码转换:
```python
img_array = cv2.cvtColor(img_array, cv2.COLOR_BGR2RGB) # OpenCV转PIL顺序
```
3. **多文件批量处理**
使用`glob`库遍历文件夹:
```python
import glob
for npy_path in glob.glob('input_folder/*.npy'):
# 处理每个文件
```
---
阅读全文
相关推荐













