1# -*- coding: utf-8 -*-
2# imageio is distributed under the terms of the (new) BSD License.
3
4"""Read DICOM files.
5
6Backend Library: internal
7
8A format for reading DICOM images: a common format used to store
9medical image data, such as X-ray, CT and MRI.
10
11This format borrows some code (and ideas) from the pydicom project. However,
12only a predefined subset of tags are extracted from the file. This allows
13for great simplifications allowing us to make a stand-alone reader, and
14also results in a much faster read time.
15
16By default, only uncompressed and deflated transfer syntaxes are supported.
17If gdcm or dcmtk is installed, these will be used to automatically convert
18the data. See https://github.com/malaterre/GDCM/releases for installing GDCM.
19
20This format provides functionality to group images of the same
21series together, thus extracting volumes (and multiple volumes).
22Using volread will attempt to yield a volume. If multiple volumes
23are present, the first one is given. Using mimread will simply yield
24all images in the given directory (not taking series into account).
25
26Parameters
27----------
28progress : {True, False, BaseProgressIndicator}
29 Whether to show progress when reading from multiple files.
30 Default True. By passing an object that inherits from
31 BaseProgressIndicator, the way in which progress is reported
32 can be costumized.
33
34"""
35
36# todo: Use pydicom:
37# * Note: is not py3k ready yet
38# * Allow reading the full meta info
39# I think we can more or less replace the SimpleDicomReader with a
40# pydicom.Dataset For series, only ned to read the full info from one
41# file: speed still high
42# * Perhaps allow writing?
43
44import os
45import sys
46import logging
47import subprocess
48
49from ..core import Format, BaseProgressIndicator, StdoutProgressIndicator
50from ..core import read_n_bytes
51
52_dicom = None # lazily loaded in load_lib()
53
54logger = logging.getLogger(__name__)
55
56
57def load_lib():
58 global _dicom
59 from . import _dicom
60
61 return _dicom
62
63
64# Determine endianity of system
65sys_is_little_endian = sys.byteorder == "little"
66
67
68def get_dcmdjpeg_exe():
69 fname = "dcmdjpeg" + ".exe" * sys.platform.startswith("win")
70 for dir in (
71 "c:\\dcmtk",
72 "c:\\Program Files",
73 "c:\\Program Files\\dcmtk",
74 "c:\\Program Files (x86)\\dcmtk",
75 ):
76 filename = os.path.join(dir, fname)
77 if os.path.isfile(filename):
78 return [filename]
79
80 try:
81 subprocess.check_call([fname, "--version"])
82 return [fname]
83 except Exception:
84 return None
85
86
87def get_gdcmconv_exe():
88 fname = "gdcmconv" + ".exe" * sys.platform.startswith("win")
89 # Maybe it's on the path
90 try:
91 subprocess.check_call([fname, "--version"])
92 return [fname, "--raw"]
93 except Exception:
94 pass
95 # Select directories where it could be
96 candidates = []
97 base_dirs = [r"c:\Program Files"]
98 for base_dir in base_dirs:
99 if os.path.isdir(base_dir):
100 for dname in os.listdir(base_dir):
101 if dname.lower().startswith("gdcm"):
102 suffix = dname[4:].strip()
103 candidates.append((suffix, os.path.join(base_dir, dname)))
104 # Sort, so higher versions are tried earlier
105 candidates.sort(reverse=True)
106 # Select executable
107 filename = None
108 for _, dirname in candidates:
109 exe1 = os.path.join(dirname, "gdcmconv.exe")
110 exe2 = os.path.join(dirname, "bin", "gdcmconv.exe")
111 if os.path.isfile(exe1):
112 filename = exe1
113 break
114 if os.path.isfile(exe2):
115 filename = exe2
116 break
117 else:
118 return None
119 return [filename, "--raw"]
120
121
122class DicomFormat(Format):
123 """See :mod:`imageio.plugins.dicom`"""
124
125 def _can_read(self, request):
126 # If user URI was a directory, we check whether it has a DICOM file
127 if os.path.isdir(request.filename):
128 files = os.listdir(request.filename)
129 for fname in sorted(files): # Sorting make it consistent
130 filename = os.path.join(request.filename, fname)
131 if os.path.isfile(filename) and "DICOMDIR" not in fname:
132 with open(filename, "rb") as f:
133 first_bytes = read_n_bytes(f, 140)
134 return first_bytes[128:132] == b"DICM"
135 else:
136 return False
137 # Check
138 return request.firstbytes[128:132] == b"DICM"
139
140 def _can_write(self, request):
141 # We cannot save yet. May be possible if we will used pydicom as
142 # a backend.
143 return False
144
145 # --
146
147 class Reader(Format.Reader):
148 _compressed_warning_dirs = set()
149
150 def _open(self, progress=True):
151 if not _dicom:
152 load_lib()
153 if os.path.isdir(self.request.filename):
154 # A dir can be given if the user used the format explicitly
155 self._info = {}
156 self._data = None
157 else:
158 # Read the given dataset now ...
159 try:
160 dcm = _dicom.SimpleDicomReader(self.request.get_file())
161 except _dicom.CompressedDicom as err:
162 # We cannot do this on our own. Perhaps with some help ...
163 cmd = get_gdcmconv_exe()
164 if not cmd and "JPEG" in str(err):
165 cmd = get_dcmdjpeg_exe()
166 if not cmd:
167 msg = err.args[0].replace("using", "installing")
168 msg = msg.replace("convert", "auto-convert")
169 err.args = (msg,)
170 raise
171 else:
172 fname1 = self.request.get_local_filename()
173 fname2 = fname1 + ".raw"
174 try:
175 subprocess.check_call(cmd + [fname1, fname2])
176 except Exception:
177 raise err
178 d = os.path.dirname(fname1)
179 if d not in self._compressed_warning_dirs:
180 self._compressed_warning_dirs.add(d)
181 logger.warning(
182 "DICOM file contained compressed data. "
183 + "Autoconverting with "
184 + cmd[0]
185 + " (this warning is shown once for each directory)"
186 )
187 dcm = _dicom.SimpleDicomReader(fname2)
188
189 self._info = dcm._info
190 self._data = dcm.get_numpy_array()
191
192 # Initialize series, list of DicomSeries objects
193 self._series = None # only created if needed
194
195 # Set progress indicator
196 if isinstance(progress, BaseProgressIndicator):
197 self._progressIndicator = progress
198 elif progress is True:
199 p = StdoutProgressIndicator("Reading DICOM")
200 self._progressIndicator = p
201 elif progress in (None, False):
202 self._progressIndicator = BaseProgressIndicator("Dummy")
203 else:
204 raise ValueError("Invalid value for progress.")
205
206 def _close(self):
207 # Clean up
208 self._info = None
209 self._data = None
210 self._series = None
211
212 @property
213 def series(self):
214 if self._series is None:
215 pi = self._progressIndicator
216 self._series = _dicom.process_directory(self.request, pi)
217 return self._series
218
219 def _get_length(self):
220 if self._data is None:
221 dcm = self.series[0][0]
222 self._info = dcm._info
223 self._data = dcm.get_numpy_array()
224
225 nslices = self._data.shape[0] if (self._data.ndim == 3) else 1
226
227 if self.request.mode[1] == "i":
228 # User expects one, but lets be honest about this file
229 return nslices
230 elif self.request.mode[1] == "I":
231 # User expects multiple, if this file has multiple slices, ok.
232 # Otherwise we have to check the series.
233 if nslices > 1:
234 return nslices
235 else:
236 return sum([len(serie) for serie in self.series])
237 elif self.request.mode[1] == "v":
238 # User expects a volume, if this file has one, ok.
239 # Otherwise we have to check the series
240 if nslices > 1:
241 return 1
242 else:
243 return len(self.series) # We assume one volume per series
244 elif self.request.mode[1] == "V":
245 # User expects multiple volumes. We have to check the series
246 return len(self.series) # We assume one volume per series
247 else:
248 raise RuntimeError("DICOM plugin should know what to expect.")
249
250 def _get_slice_data(self, index):
251 nslices = self._data.shape[0] if (self._data.ndim == 3) else 1
252
253 # Allow index >1 only if this file contains >1
254 if nslices > 1:
255 return self._data[index], self._info
256 elif index == 0:
257 return self._data, self._info
258 else:
259 raise IndexError("Dicom file contains only one slice.")
260
261 def _get_data(self, index):
262 if self._data is None:
263 dcm = self.series[0][0]
264 self._info = dcm._info
265 self._data = dcm.get_numpy_array()
266
267 nslices = self._data.shape[0] if (self._data.ndim == 3) else 1
268
269 if self.request.mode[1] == "i":
270 return self._get_slice_data(index)
271 elif self.request.mode[1] == "I":
272 # Return slice from volume, or return item from series
273 if index == 0 and nslices > 1:
274 return self._data[index], self._info
275 else:
276 L = []
277 for serie in self.series:
278 L.extend([dcm_ for dcm_ in serie])
279 return L[index].get_numpy_array(), L[index].info
280 elif self.request.mode[1] in "vV":
281 # Return volume or series
282 if index == 0 and nslices > 1:
283 return self._data, self._info
284 else:
285 return (
286 self.series[index].get_numpy_array(),
287 self.series[index].info,
288 )
289 # mode is `?` (typically because we are using V3). If there is a
290 # series (multiple files), index referrs to the element of the
291 # series and we read volumes. If there is no series, index
292 # referrs to the slice in the volume we read "flat" images.
293 elif len(self.series) > 1:
294 # mode is `?` and there are multiple series. Each series is a ndimage.
295 return (
296 self.series[index].get_numpy_array(),
297 self.series[index].info,
298 )
299 else:
300 # mode is `?` and there is only one series. Each slice is an ndimage.
301 return self._get_slice_data(index)
302
303 def _get_meta_data(self, index):
304 if self._data is None:
305 dcm = self.series[0][0]
306 self._info = dcm._info
307 self._data = dcm.get_numpy_array()
308
309 nslices = self._data.shape[0] if (self._data.ndim == 3) else 1
310
311 # Default is the meta data of the given file, or the "first" file.
312 if index is None:
313 return self._info
314
315 if self.request.mode[1] == "i":
316 return self._info
317 elif self.request.mode[1] == "I":
318 # Return slice from volume, or return item from series
319 if index == 0 and nslices > 1:
320 return self._info
321 else:
322 L = []
323 for serie in self.series:
324 L.extend([dcm_ for dcm_ in serie])
325 return L[index].info
326 elif self.request.mode[1] in "vV":
327 # Return volume or series
328 if index == 0 and nslices > 1:
329 return self._info
330 else:
331 return self.series[index].info
332 else: # pragma: no cover
333 raise ValueError("DICOM plugin should know what to expect.")