简介：面向医学影像处理初学者与需要在 MATLAB 中完成 DICOM 数据转换的科研人员，这份资源围绕 DCM 读取及 NII、NPZ、MAT 三种常用格式导出需求，整理了可直接调用的函数脚本与配套说明，适合医学图像预处理、神经影像数据整理等场景。压缩包共 9 个文件，约 20KB，其中 8 个为 .m 脚本或函数，涵盖 NIfTI 头信息读写、图像加载与保存、坐标变换以及三维数据提取等环节，另有 1 个 txt 说明文件用于辅助理解目录与调用关系。已有 3181 人学习下载。读者可据此搭建从 DICOM 读取到 NIfTI 转换的基础工具链，并参考 MAT、NPZ 导出思路完成数据格式适配；函数按职责拆分，便于在现有 MATLAB 工程中复用、调试和二次修改。
1. 三千个 .dcm 文件铺在一个文件夹里，怎么把它们变成能直接用的体数据
三千多个 .dcm 文件铺在一个文件夹里，文件名是一串没规律的编号，双击打不开，下游训练脚本却只认 .nii.gz 和 .npz。用 Matlab 读取 DCM 保存为 NII、NPZ 或 MAT，是这类场景里最短的一条路——不用为格式转换单独装一整套工具链，几十行代码就能把一层层切片还原成带几何信息的体数据，再按需要写成三种格式里的任意一种。这份资源做的就是这样一套转换工程：索引序列、排层序、拼仿射矩阵、写文件、回读校验，整条链路都在 Matlab 里闭环。适合手上压着一批原始 DICOM、又想直接在 Matlab 里做预处理的人；也适合把 DICOM 先转成 NIfTI 再交给别的框架训练的场景。
2. 工程骨架与选型：为什么读 DICOM 用 dicomread 而不是自己 fread
动手之前先定两件事：文件怎么摆、用什么读。这两件事定错了，后面每一行代码都在还债。DICOM 不是“一堆像素加个文件头”那么简单，它是一个允许厂商塞私有 tag 的容器，传输语法（隐式/显式 VR、大小端、有无压缩）随设备而变。自己fread解析等于从零写一个解析器，血泪经验是：做到第三个厂家的机器就崩了。
2.1 目录结构与模块边界
工程规模不大，但要按职责分开，不然改一处格式就得翻三个文件。我一般这样摆：
dcm2vol/ ├── main_convert.m % 单序列转换入口 ├── batch_convert.m % 整目录批处理入口 ├── +dcmio/ % 只负责「读」，不碰任何输出格式 │ ├── readSeries.m │ ├── sortSlices.m │ ├── buildAffine.m │ └── readPixels.m ├── +export/ % 只负责「写」，不碰 DICOM 细节 │ ├── writeNii.m │ ├── writeNpz.m │ ├── writeMat.m │ └── writeManifest.m ├── config/ │ └── default.json % 输出格式、是否反相、是否去重等开关 ├── tools/ │ └── mat2npz.py % 没有 Python 接口时的兜底转换脚本 └── tests/ └── test_roundtrip.m % 写出去再读回来比对包目录+dcmio和+export在 MATLAB 里是命名空间包，调用写成dcmio.readSeries(...)、export.writeNii(...)。这么做最直接的好处是：将来要加一个writeNrrd.m，只需要在+export下丢一个文件，+dcmio一行都不用改。
| 文件 | 输入 | 输出 | 职责边界 |
|---|---|---|---|
readSeries.m | 序列文件夹路径 | 结构体vol | 调度排序、像素、几何三步 |
sortSlices.m | dicominfo元数据数组 | 排序索引 + 层间距 | 只决定“谁在第几层” |
buildAffine.m | 参考层元数据 + 层间距 | 4×4 仿射矩阵 | 只决定“体素坐标怎么映射到病人坐标” |
writeNii.m | vol | .nii/.nii.gz | 只填 NIfTI 头 |
writeNpz.m | vol | .npz | 只处理 Matlab→numpy 的桥接 |
writeMat.m | vol | .mat | 只处理变量打包与大文件选项 |
2.2 读 DICOM 为什么交给 dicomread 和 dicomCollection
dicomread解决的是“把一个文件里的像素和元数据拿出来”这一层。它认得常见的未压缩传输语法，也对 JPEG、JPEG Lossless 这类压缩语法有支持；碰到 JPEG2000 或者某些厂商私有压缩，它会直接报错而不是给一堆错值——报错比错值好，至少能定位。真正的代价是：你得装着 Image Processing Toolbox。
dicomCollection解决的是“一个目录里到底有几个序列”这一层。一个从影像科拷出来的检查目录，常见情况是主序列、定位像、剂量报告、结构化报告混在一起，文件名毫无规律。这个函数会扫一遍目录并返回一张表，里面有SeriesInstanceUID、Modality、Filenames这些列，等于帮你先分好组。没有它的时候，大家是靠dir('*.dcm')之后自己读第一个文件的SeriesInstanceUID去过滤，能用，但每次都得重新写一遍。
2.3 环境要求与一次最小验证
工具箱层面需要 Image Processing Toolbox（dicomread、dicomwrite、niftiwrite都在里面）；走 NPZ 那条路的话，额外要一个能被 Matlab 调起来的 Python，且里面装了 numpy。先别急着写完整的转换器，用几行代码确认这条链路是通的：
% 最小验证：读一层，确认数据长得对，再谈批量 srcDir = 'D:\data\series_01'; f = fullfile(srcDir, 'IM_0001.dcm'); info = dicominfo(f); % 元数据，先肉眼看一遍关键字段 I = dicomread(info); % 直接传 info 比传路径快，省一次解析 fprintf('尺寸 %s 类型 %s\n', mat2str(size(I)), class(I)); fprintf('模态 %s 光度解释 %s\n', info.Modality, info.PhotometricInterpretation); fprintf('像素间距 %s 层厚 %g\n', mat2str(info.PixelSpacing(:)'), info.SliceThickness); imshow(I, []); % 第二个参数 [] 是关键，否则 int16 会显示成一片黑dicomread(info)这个写法值得单独说：直接传dicominfo的结果，比传文件路径少解析一次头，批量跑几千层时这个差别肉眼可见。imshow(I, [])里的[]是自动窗宽窗位，不加的话int16数据里超出 0~255 的部分会被直接截断，屏幕上就是一片黑——这个坑几乎每个人都踩过，然后怀疑是文件坏了。
跑通这一步，再往下就是把这个零散流程固化成readSeries。
3. 读取层拆解：把一个序列从文件夹还原成带几何信息的体数据
读序列这件事，难的从来不是dicomread那一行，而是三件事：层序对不对、像素值量纲对不对、几何信息全不全。这三件事任何一个出问题，写出来的 NIfTI 在别人机器上打开就是错的，而且往往要等到模型训歪了才发现。
3.1 用 dicomCollection 建索引并按序列分组
第一步是把“目录”变成“序列”。注意别把整个检查目录直接丢进去当成一个体数据：
function [fileList, seriesUID] = pickSeries(srcDir, wantUID) % srcDir : 检查目录或序列目录 % wantUID : 指定要哪个序列；留空则自动取文件数最多的那个 t = dicomCollection(srcDir); if isempty(t) error('目录里没有可识别的 DICOM：%s', srcDir); end if nargin >= 2 && ~isempty(wantUID) row = find(string(t.SeriesInstanceUID) == string(wantUID), 1); if isempty(row), error('找不到序列 %s', wantUID); end else % 文件数最多的通常是主序列，定位像和报告一般只有几帧 [~, row] = max(cellfun(@numel, t.Filenames)); end fileList = t.Filenames{row}; fileList = cellstr(string(fileList)); % 统一成 cellstr，后面好索引 seriesUID = string(t.SeriesInstanceUID(row)); endcellfun(@numel, t.Filenames)这一步是在数每个序列有几帧。定位像通常只有两三帧，剂量报告是一段文本，主序列动辄几百上千帧，取最大值这个启发式在绝大多数检查里都成立。如果某个检查有增强前后两期，两期帧数可能一样多，这时候就必须显式传wantUID——SeriesInstanceUID可以从dicomCollection返回的表里先打印出来看。
cellstr(string(fileList))是防类型不一致：有些版本返回字符串数组，有些返回字符向量元胞，后面用{}索引时两种行为不一样，统一掉能省掉一堆class判断。
3.2 排序：别信 InstanceNumber
层序错了，整卷数据就是乱的。InstanceNumber看着最像排序依据，实际最不靠谱：有的设备从 1 开始，有的从 0，有的跳号，有的多序列合并后重复。
真正稳的做法是把每层的ImagePositionPatient投影到层法向量上，按投影值排序：
function [order, dz, flipped] = sortSlices(infoCell) % infoCell : 每层的 dicominfo，cell 数组 % order : 排好序的索引 % dz : 层间距（mm） % flipped : 法向量与层序号增长方向是否相反 n = numel(infoCell); ipp = zeros(n, 3); inst = zeros(n, 1); for k = 1:n if isfield(infoCell{k}, 'ImagePositionPatient') ipp(k, :) = infoCell{k}.ImagePositionPatient(:)'; end inst(k) = infoCell{k}.InstanceNumber; end iop = infoCell{1}.ImageOrientationPatient(:)'; normal = cross(iop(1:3), iop(4:6)); % 层法向量，未归一化也能排序 if any(ipp(:) ~= 0) proj = ipp * normal(:); % 每层在法向量上的带符号位置 [projSorted, order] = sort(proj); dz = median(diff(projSorted)); if dz < 0 % 层间距算出负值说明法向量方向和实际扫描方向相反 order = flipud(order); dz = -dz; flipped = true; else flipped = false; end else % 兜底：老设备可能不写 ImagePositionPatient [~, order] = sort(inst); dz = infoCell{order(1)}.SliceThickness; flipped = false; end end几个参数的含义要拎清：ImageOrientationPatient是六个值，前三个是行方向余弦、后三个是列方向余弦，叉乘得到层法向量；ImagePositionPatient是每层左上角像素中心在病人坐标系（LPS）里的位置；proj就是这层沿着法向量走了多远。取median(diff(...))而不是直接diff的第 1 个元素，是因为偶尔会有重复层或者个别设备插了个异常值，中位数扛得住。
flipped这个返回值不是可有可无的。如果dz算出负数，说明你按法向量排序得到的顺序和真实扫描方向相反，这时候要么翻order，要么后面把体数据沿第三维flip，两者选一个，别两个都做。
3.3 像素值：RescaleSlope、PixelRepresentation 和 MONOCHROME1
dicomread给你的是存储值，不是物理值。CT 的存储值要经过RescaleSlope和RescaleIntercept才变成 HU；这一步漏了，写出来的 NIfTI 看着也挺像那么回事，但窗宽窗位全错，骨和软组织糊成一团。
function [vol, meta] = readPixels(fileList, order, ref) % fileList : 原始文件列表（未排序） % order : sortSlices 返回的排序索引 % ref : 首层的 dicominfo n = numel(order); [r, c] = size(dicomread(ref)); % 先建 double 累加，最后统一 cast，避免中途截断 vol = zeros(r, c, n, 'double'); for k = 1:n info = dicominfo(fileList{order(k)}); I = dicomread(info); % 多帧 DICOM 会返回 4 维，把帧轴挤到第 3 维 if ndims(I) == 4 I = reshape(I, size(I,1), size(I,2), []); end I = double(I(:,:,1)); % 关键一步：存储值 -> 物理值 slope = 1; intercept = 0; if isfield(info, 'RescaleSlope'), slope = info.RescaleSlope; end if isfield(info, 'RescaleIntercept'), intercept = info.RescaleIntercept; end I = I * slope + intercept; vol(:,:,k) = I; end % MONOCHROME1 是反相显示，最小值代表最亮，不翻回来就是一张负片 inverted = strcmpi(ref.PhotometricInterpretation, 'MONOCHROME1'); if inverted vol = (max(vol(:)) + min(vol(:))) - vol; end meta = struct( ... 'SeriesInstanceUID', char(ref.SeriesInstanceUID), ... 'Modality', ref.Modality, ... 'PhotometricInterpretation', ref.PhotometricInterpretation, ... 'RescaleSlope', ref.RescaleSlope, ... 'RescaleIntercept', ref.RescaleIntercept, ... 'RescaleType', 'HU', ... 'Inverted', inverted, ... 'PatientPosition', ref.PatientPosition); end为什么要先转double再算：如果文件本身是int16（PixelRepresentation = 1，CT 常见），直接I * slope会在int16里做饱和运算，-1024 以下的 HU 全被夹到 -32768，肺野那种低密度区域直接糊掉。double只是中转，最后落盘时再按需要cast回int16，内存和精度都照顾到了。
Inverted这个标志写进 meta 而不是只体现在数据里，是因为下游如果要做统一的归一化，知道“这卷是翻过的”比“这卷看起来是正的”更有用。
3.4 几何：从 ImageOrientationPatient 拼出 4×4 仿射矩阵
NIfTI 头里最核心的就是那个 4×4 仿射矩阵，它决定了体素坐标(i,j,k)怎么映射到世界坐标(x,y,z)。填错了，体数据本身没错，但配准、融合、ROI 映射全废。
function A = buildAffine(info, dz) % info : 排序后首层的 dicominfo % dz : 层间距（mm，正数） % A : 4x4，RAS+ 世界坐标 iop = info.ImageOrientationPatient(:)'; rowDir = iop(1:3) / norm(iop(1:3)); % 列索引 j 增长方向 colDir = iop(4:6) / norm(iop(4:6)); % 行索引 i 增长方向 normal = cross(rowDir, colDir); % 层索引 k 增长方向 normal = normal / norm(normal); ps = info.PixelSpacing(:)'; % <p> <a href="https://download.csdn.net/download/tianjinyikedaxue/15841878" style="color:#ec7500;font-size:14px;"> 本文还有配套的精品资源，点击获取 </a> <img alt="menu-r.4af5f7ec.gif" src="https://csdnimg.cn/release/wenkucmsfe/public/img/menu-r.4af5f7ec.gif" style="width:16px;margin-left:4px;vertical-align:text-bottom;cursor:text;"> </p>
