医学数据处理：mat文件的3D矩阵转dicom文件

滴~~~~~~

已于 2022-06-16 15:23:33 修改

阅读量1k

点赞数 1

文章标签： python

于 2022-06-06 20:04:52 首次发布

本文链接：https://blog.csdn.net/weixin_43730079/article/details/125152760

版权

文章目录

前言
一、dicom
二、具体代码
总结

前言

因为想把3维矩阵进行3D显示，但原本的数据存储在mat文件里，3D显示需要dicom形式的，所以需要对数据进行处理。具体内容还是参考了一些文章，比如dicom形式的文章，mat的3D矩阵转dicom的文章。

一、dicom

具体格式可以参考其他文章，我在转换时遇到了一些问题，比如modality这个tag，因为我需要得到CT图像并在小赛看看里观察，所以这个tag需要设置为‘CT’。和CT图像相关的还有一些tag，不过我用sitk.writeimage写入后没有问题就没关注。

参考：https://blog.csdn.net/sinolover/article/details/114311823?utm_medium=distribute.pc_relevant.none-task-blog-2_defaultbaidujs_baidulandingword~default-0-114311823-blog-101880158.pc_relevant_paycolumn_v3&spm=1001.2101.3001.4242.1&utm_relevant_index=3

二、具体代码

代码如下：

import numpy as np
import SimpleITK as sitk
import pydicom 
from scipy.io import loadmat
import os

def mat2dicom(folderPath):
    count_study = 0
    folderPath_new = folderPath + '/3D'#这个是我的路径相关的处理
    all_study = os.listdir(folderPath_new)
    for every_study_idx in range(len(all_study)):#读取每个增强序列对应的mat文件（3维矩阵）
        count_study +=1
        tmp_MAT_path = os.path.join(folderPath_new,all_study[every_study_idx])
        source_MAT_path = folderPath
        final_mat_path = tmp_MAT_path
        #读取三维矩阵 mat文件
        mat = loadmat(final_mat_path)#读mat文件，用loadmat来处理了
        # 可以用shape查看维度信息
        mat_t = mat.get('tensor')#tensor是我的数据的key
        # print(type(mat_t),mat_t.shape,mat_t.shape[2]) #可以看看得到数据的shape
        #定义存储dicom的路径
        every_DCM_path = source_MAT_path + '/DCM' + str(every_study_idx)
        print(every_DCM_path)#打印路径
        begin = 0#表示起始位置的参数
        flag = np.zeros((88,120))#画图需要有信息的数据，用这个全0矩阵来去掉没信息的dcm
        if not os.path.exists(every_DCM_path):
            os.makedirs(every_DCM_path)
        for i in range(mat_t.shape[0]):#一个mat中的三维矩阵有k维，这个循环就对k维循环输出每一层的dcm
            out = mat_t[i, :, :].astype('int16')#这个int16是必须要增加的，不能直接使用float，在实验的时候一改就报错，且int和int16不同
            if((out == flag).all()):#保留有信息的dcm
                continue
            InstanceNumber = begin + i
            if InstanceNumber <10:
                InstanceNumber = '0000'+str(InstanceNumber)
            elif InstanceNumber <100:
                InstanceNumber = '000'+str(InstanceNumber)
            elif InstanceNumber <200:
                InstanceNumber = '00'+str(InstanceNumber)
            else:
                print('Warning!!',InstanceNumber,"dicom数量大于200")
            
            dicom_save_file_path = every_DCM_path + '/'+InstanceNumber +'.dcm'
            if os.path.exists(dicom_save_file_path):
                continue
            sitk.WriteImage(sitk.GetImageFromArray(out),dicom_save_file_path)
            
            ######再读文件，修改dicom的tag信息使这些图片成为一个序列
            read_dicom_path = every_DCM_path+'/'+InstanceNumber+'.dcm'  #每个dcm文件的具体路径

            dcm = pydicom.read_file(read_dicom_path)
            ds = pydicom.read_file(dicom_save_file_path)
            #print(dcm)
            print(dicom_save_file_path)
            ds = modify_dicom(ds,dcm)
            #保存
            ds.save_as(dicom_save_file_path)


def modify_dicom(ds,dcm):#这个函数里注释掉比较多的内容，是因为在实验中报错的地方全注释了，应该是我的数据没有这些
     #PatientInfo
    ds.PatientID = dcm.PatientID
    ds.PatientName =  dcm.PatientName
    ds.PatientBirthDate = dcm.PatientBirthDate
    ds.PatientSex = dcm.PatientSex
    #ds.PatientAge = dcm.PatientAge
    #ds.PatientWeight = dcm.PatientWeight
    #try:
        #ds.MagneticFieldStrength = dcm.MagneticFieldStrength
        #ds.Manufacturer = dcm.Manufacturer
        #ds.InstitutionName = dcm.InstitutionName
    #except TypeError:
        #print("Error:没有InstitutionName、Manufacturer、InstitutionName tag")
    #studyInfo
    ds.StudyDate =dcm.StudyDate
    ds.StudyTime = dcm.StudyTime
    #ds.StudyDescription = dcm.StudyDescription
    ds.StudyInstanceUID = dcm.StudyInstanceUID
    ds.StudyID = ds.StudyID
    # seriesInfo
    ds.SeriesInstanceUID = dcm.SeriesInstanceUID #修改Series Instance UID
    #ds.SeriesDescription = dcm.SeriesDescription
    ds.Modality = 'CT'#这个是为了适应需要获得CT图像而定义的
    ds.SeriesNumber = dcm.SeriesNumber
    ds.InstanceNumber  = dcm.InstanceNumber
    #ds.SeriesDate = dcm.SeriesDate
    #ds.SeriesTime = dcm.SeriesTime
    #ds.SliceThickness = dcm.SliceThickness
    #ds.SliceLocation = dcm.SliceLocation
    ds.FrameOfReferenceUID = dcm.FrameOfReferenceUID
    ds[0X0020, 0X0052].value = dcm[0X0020, 0X0052].value 
    #ImageInfo
    ds.SOPClassUID  = dcm.SOPClassUID
    ds.SOPInstanceUID = dcm.SOPInstanceUID
    ds.PixelData = dcm.PixelData
    #ds.ReferencedSOPInstanceUID =dcm.ReferencedSOPInstanceUID
    #ds.ReferencedSOPClassUID =dcm.ReferencedSOPClassUID
    return ds

if __name__=="__main__":
    folderPath ='D:/code final/data/3D'#路径
    mat2dicom(folderPath)

修改代码的时候遇到的问题主要是dicom形式不熟悉、sitk.writeimage写入图像时out必须是int16、modality这个tag改为CT、去除没有信息的dcm。其中int16和modality这两个问题耗费时间比较久，int16也找了网上的其他文章，也有人遇到这种问题，modality这个是因为想在小赛看看里看到具体的CT图像，而不修改为CT在里面看到的都是全黑的。

代码参考：https://www.jianshu.com/p/4afb1cb5b5c6

总结

具体的代码还得参考实验来修改，不过对dicom形式有了一定了解。

后记：
1.发现得到的dicom不是一个序列的，看了一个病例dicom中的tag信息，发现很多UID都是一样，然后实验了一下，自定义了所有UID，结果表明，确实能够产生一个序列的dicom。
2.Frame UID和SOP UID需要不同，否则小赛看看读取的是一个序列，但是radiant读取的不是一个序列。