從numpy數組創建pydicom文件

Question

我試圖從標準大小（512 x 512或256 x 256）numpy陣列創建一個mew dicom映像。看起來這應該是直截了當的，我已經調整了我的代碼http://code.google.com/p/pydicom/source/browse/source/dicom/examples/write_new.py，它似乎執行相同的過程，但是當我保存文件時，我無法在RadiAnt或MicroDicom中查看它。

import dicom, dicom.UID 
from dicom.dataset import Dataset, FileDataset 

def write_dicom(pixel_array,filename): 

    file_meta = Dataset() 
    ds = FileDataset(filename, {},file_meta = file_meta,preamble="\0"*128) 
    ds.PixelData = pixel_array.tostring() 
    ds.save_as(filename) 
    return 

if __name__ == "__main__": 
    import numpy as np 
    pixel_array = np.tile(np.arange(256).reshape(16,16), (16,16)) * 4 
    write_dicom(pixel_array,'pretty.dcm') 

Answer 1

在上面的例子中工作，但會導致許多工具抱怨DICOMs，他們甚至不能在所有使用ITK/SimpleITK作爲堆棧讀。我發現從numpy製作DICOM的最好方法是使用SimpleITK工具並逐片生成DICOM。一個基本的例子（https://github.com/zivy/SimpleITK/blob/8e94451e4c0e90bcc6a1ffdd7bc3d56c81f58d80/Examples/DicomSeriesReadModifyWrite/DicomSeriesReadModifySeriesWrite.py）示出如何在一個堆棧加載，執行變換，然後重新保存該文件，但可以很容易地通過使用最終在輸出圖像的標籤

import SimpleITK as sitk 
filtered_image = sitk.GetImageFromArray(my_numpy_array) 

數量進行修改是相當大的所以手動創建它們都很乏味。此外，SimpleITK支持8,16,32位圖像以及RGB，因此比在pydicom中更容易。

(0008, 0008) Image Type       CS: ['DERIVED', 'SECONDARY'] 
(0008, 0016) SOP Class UID      UI: Secondary Capture Image Storage 
(0008, 0018) SOP Instance UID     UI: 1.2.826.0.1.3680043.2.1125.1.35596048796922805578234000521866725 
(0008, 0020) Study Date       DA: '20170803' 
(0008, 0021) Series Date       DA: '20170803' 
(0008, 0023) Content Date      DA: 0 
(0008, 0030) Study Time       TM: '080429.171808' 
(0008, 0031) Series Time       TM: '080429' 
(0008, 0033) Content Time      TM: 0 
(0008, 0050) Accession Number     SH: '' 
(0008, 0060) Modality       CS: 'OT' 
(0008, 0064) Conversion Type      CS: 'WSD' 
(0008, 0090) Referring Physician's Name   PN: '' 
(0010, 0010) Patient's Name      PN: '' 
(0010, 0020) Patient ID       LO: '' 
(0010, 0030) Patient's Birth Date    DA: '' 
(0010, 0040) Patient's Sex      CS: '' 
(0018, 2010) Nominal Scanned Pixel Spacing  DS: ['1', '3'] 
(0020, 000d) Study Instance UID     UI: 1.2.826.0.1.3680043.2.1125.1.33389357207068897066210100430826006 
(0020, 000e) Series Instance UID     UI: 1.2.826.0.1.3680043.2.1125.1.51488923827429438625199681257282809 
(0020, 0010) Study ID       SH: '' 
(0020, 0011) Series Number      IS: '' 
(0020, 0013) Instance Number      IS: '' 
(0020, 0020) Patient Orientation     CS: '' 
(0020, 0052) Frame of Reference UID    UI: 1.2.826.0.1.3680043.2.1125.1.35696880630664441938326682384062489 
(0028, 0002) Samples per Pixel     US: 1 
(0028, 0004) Photometric Interpretation   CS: 'MONOCHROME2' 
(0028, 0010) Rows        US: 40 
(0028, 0011) Columns        US: 50 
(0028, 0100) Bits Allocated      US: 32 
(0028, 0101) Bits Stored       US: 32 
(0028, 0102) High Bit       US: 31 
(0028, 0103) Pixel Representation    US: 1 
(0028, 1052) Rescale Intercept     DS: "0" 
(0028, 1053) Rescale Slope      DS: "1" 
(0028, 1054) Rescale Type      LO: 'US' 
(7fe0, 0010) Pixel Data       OW: Array of 8000 bytes 

Answer 2

DICOM是一個非常複雜的格式。有很多方言，相容性是一個運氣問題。您也可以嘗試nibabel，或許它的方言對RadiAnt或MicroDicom更具吸引力。

一般來說，我會建議儘可能使用Nifti格式。它的標準更加簡潔，不兼容性很少。尼巴​​貝爾也支持這一點。

Answer 3

這是我需要編寫的代碼的功能版本。它將從給定的二維像素陣列中寫入一個16位灰度DICOM圖像。根據DICOM標準，對於每個圖像和系列，UID應該是唯一的，這個代碼並不擔心，因爲我不知道UID實際上做了什麼。如果別人做，我會很樂意將它添加英寸

import dicom, dicom.UID 
from dicom.dataset import Dataset, FileDataset 
import numpy as np 
import datetime, time 

def write_dicom(pixel_array,filename): 
    """ 
    INPUTS: 
    pixel_array: 2D numpy ndarray. If pixel_array is larger than 2D, errors. 
    filename: string name for the output file. 
    """ 

    ## This code block was taken from the output of a MATLAB secondary 
    ## capture. I do not know what the long dotted UIDs mean, but 
    ## this code works. 
    file_meta = Dataset() 
    file_meta.MediaStorageSOPClassUID = 'Secondary Capture Image Storage' 
    file_meta.MediaStorageSOPInstanceUID = '1.3.6.1.4.1.9590.100.1.1.111165684411017669021768385720736873780' 
    file_meta.ImplementationClassUID = '1.3.6.1.4.1.9590.100.1.0.100.4.0' 
    ds = FileDataset(filename, {},file_meta = file_meta,preamble="\0"*128) 
    ds.Modality = 'WSD' 
    ds.ContentDate = str(datetime.date.today()).replace('-','') 
    ds.ContentTime = str(time.time()) #milliseconds since the epoch 
    ds.StudyInstanceUID = '1.3.6.1.4.1.9590.100.1.1.124313977412360175234271287472804872093' 
    ds.SeriesInstanceUID = '1.3.6.1.4.1.9590.100.1.1.369231118011061003403421859172643143649' 
    ds.SOPInstanceUID = '1.3.6.1.4.1.9590.100.1.1.111165684411017669021768385720736873780' 
    ds.SOPClassUID = 'Secondary Capture Image Storage' 
    ds.SecondaryCaptureDeviceManufctur = 'Python 2.7.3' 

    ## These are the necessary imaging components of the FileDataset object. 
    ds.SamplesPerPixel = 1 
    ds.PhotometricInterpretation = "MONOCHROME2" 
    ds.PixelRepresentation = 0 
    ds.HighBit = 15 
    ds.BitsStored = 16 
    ds.BitsAllocated = 16 
    ds.SmallestImagePixelValue = '\\x00\\x00' 
    ds.LargestImagePixelValue = '\\xff\\xff' 
    ds.Columns = pixel_array.shape[0] 
    ds.Rows = pixel_array.shape[1] 
    if pixel_array.dtype != np.uint16: 
     pixel_array = pixel_array.astype(np.uint16) 
    ds.PixelData = pixel_array.tostring() 

    ds.save_as(filename) 
    return 



if __name__ == "__main__": 
# pixel_array = np.arange(256*256).reshape(256,256) 
# pixel_array = np.tile(np.arange(256).reshape(16,16),(16,16)) 
    x = np.arange(16).reshape(16,1) 
    pixel_array = (x + x.T) * 32 
    pixel_array = np.tile(pixel_array,(16,16)) 
    write_dicom(pixel_array,'pretty.dcm')