Pixel Data - Part 3: Compression and decompression¶
In part 1 of this tutorial you learned how to access the pixel data as either the raw bytes
or a NumPy
ndarray
and in part 2 you learned how to create new pixel data and add it to a Dataset
.
In this final part you’ll learn how to compress and decompress datasets containing
Pixel Data.
Prerequisites
Installing using pip:
python -m pip install -U pydicom numpy pylibjpeg[all] pyjpegls
Installing on conda:
conda install numpy
conda install -c conda-forge pydicom
pip install pylibjpeg[all] pyjpegls
Compression of Pixel Data¶
pydicom can perform dataset compression for the the following transfer syntaxes:
JPEG-LS Lossless and JPEG-LS Near-lossless compression with pyjpegls.
JPEG 2000 Lossless and JPEG 2000 compression with pylibjpeg and pylibjpeg-openjpeg.
RLE Lossless, which doesn’t need any additional packages but can be sped up if pylibjpeg and pylibjpeg-rle are available.
For all other transfer syntaxes it’s entirely up to you to compress the Pixel Data in a manner conformant to the requirements of the DICOM Standard:
Each frame of pixel data must be compressed separately
All compressed frames must then be encapsulated.
The encapsulated byte stream is used to set the Pixel Data value
When the amount of compressed frame data is very large then it’s recommended (but not required) that an extended offset table also be included in the dataset
The VR for compressed Pixel Data is always OB
Compressing a dataset (with RLE Lossless)¶
Compression of an existing uncompressed dataset can be performed by passing the Transfer
Syntax UID of the compression method you’d like to use to Dataset.compress()
, or by using the compress()
function. We’ll be using RLE Lossless to start with, which is based on the
PackBits compression scheme:
>>> from pydicom import examples
>>> from pydicom.uid import RLELossless
>>> ds = examples.ct
>>> ds.file_meta.TransferSyntaxUID.is_compressed
False
>>> ds.compress(RLELossless)
If you’re creating a new dataset, or if you want to update the Pixel Data for an
existing dataset, you can pass an ndarray
along with the Transfer
Syntax UID:
import numpy as np
from pydicom import Dataset
from pydicom.uid import RLELossless
ds = Dataset()
ds.Rows = 320
ds.Columns = 480
ds.BitsAllocated = 8
ds.BitsStored = 8
ds.HighBit = ds.BitsStored - 1
ds.PixelRepresentation = 0
ds.SamplesPerPixel = 1
ds.PhotometricInterpretation = "MONOCHROME2"
arr = np.ones((ds.Rows, ds.Columns), dtype="uint8")
ds.compress(RLELossless, arr)
assert ds.file_meta.TransferSyntaxUID == RLELossless
assert isinstance(ds.PixelData, bytes)
In both cases this will compress the Dataset
in-place:
The Pixel Data will be set with the encapsulated RLE codestream
The Transfer Syntax UID will be set to RLE Lossless
A new SOP Instance UID value will be also be generated, but this can be disabled by passing
generate_instance_uid=False
.
When using an ndarray
the shape
,
dtype
and contents of arr must match the corresponding
Image Pixel module elements in the dataset,
such as Rows, Columns, Samples per Pixel, etc. If they don’t match you’ll get an exception:
>>> from pydicom import examples
>>> from pydicom.uid import RLELossless
>>> ds = examples.ct
>>> arr = np.zeros((ds.Rows, ds.Columns + 1), dtype='<i2')
>>> ds.compress(RLELossless, arr)
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
File ".../pydicom/src/pydicom/dataset.py", line 1957, in compress
encoded = [f for f in frame_iterator]
^^^^^^^^^^^^^^^^^^^^^^^^^^^
File ".../pydicom/pixels/encoders/base.py", line 678, in iter_encode
runner.validate()
File ".../pydicom/pixels/encoders/base.py", line 304, in validate
self._validate_array()
File ".../pydicom/pixels/encoders/base.py", line 333, in _validate_array
raise ValueError(
ValueError: Mismatch between the expected ndarray shape (128, 128) and the actual shape (128, 129)
When there are multiple plugins available for compressing the given transfer syntax a specific encoding plugin can be used by passing the plugin name via the encoding_plugin argument:
>>> ds.compress(RLELossless, encoding_plugin='pylibjpeg')
The RLE compression method is well supported by DICOM applications and can compress a wide range of images, however it’s usually less efficient than the JPEG family of compression schemes. More information on performing compression with RLE Lossless can be found in the RLE encoding guide.
Compressing with JPEG-LS¶
The JPEG-LS compression scheme is based on ISO/IEC 14495-1/ITU T.87. While it can compress 2- to 16-bit images and uses a lossy quality specification mechanism that’s easy to understand, it’s not well suited for lossy compression of signed integers and is generally not well supported by third-party applications, so keep that in mind if you decide to use it.
Lossless compression
Performing lossless compression is straightforward:
>>> from pydicom import examples
>>> from pydicom.uid import JPEGLSLossless
>>> ds = examples.ct
>>> ds.compress(JPEGLSLossless)
Lossy compression
Lossy compression is a bit more complicated, especially when the pixel data uses signed integers. First up though, we’ll use an example with unsigned pixel data.
Warning
pydicom makes no recommendations for specifying the image quality for lossy encoding methods. Any examples of lossy encoding are for illustration purposes only.
>>> from pydicom import examples
>>> from pydicom.uid import JPEGLSNearLossless
>>> ds = examples.rgb_color
>>> ds.PixelRepresentation
0
>>> ds.compress(JPEGLSNearLossless, jls_error=3)
The jls_error parameter is used to control the loss in image quality, and is
directly related to the JPEG-LS NEAR parameter, which is the absolute allowed error
in (unsigned) pixel data values. A jls_error of 3
therefore means that all
pixels will be within 3 intensity units of the original.
In our second lossy JPEG-LS example we’ll use a dataset with 16-bit signed integers, which is where the complication starts. The NEAR parameter is defined in terms of unsigned integers, so when used with signed values there can potentially be compression errors of up to the maximum bit-depth of the pixel data. To avoid this, the range of pixel values must be in the closed interval:
[-2**(ds.BitsStored - 1) + jls_error, 2**(ds.BitsStored - 1) - 1 - jls_error]
For example, with a Bits Stored of 8
and jls_error=3
the pixels must be in the
range [-125, 124].
>>> from pydicom import examples
>>> from pydicom.uid import JPEGLSNearLossless
>>> ds = examples.ct
>>> ds.PixelRepresentation
1
>>> ds.BitsStored
16
>>> arr = ds.pixel_array
>>> arr.min(), arr.max()
(128, 2191)
>>> ds.compress(JPEGLSNearLossless, jls_error=3)
In this example the pixel values are within the allowed range so we don’t need to do anything further. If that weren’t the case you’d have to rescale the values or use a different compression method such as JPEG 2000 (recommended).
More information on performing compression with JPEG-LS can be found in the JPEG-LS encoding guide.
Compressing with JPEG 2000¶
The JPEG 2000 compression scheme is based on ISO/IEC 15444-1/ITU T.800. The format is fairly well supported by third-party applications and it can compress images with a wide variety of properties, making it a good choice for compressing datasets.
Two transfer syntaxes are available that use JPEG 2000 compression; JPEG 2000 Lossless and JPEG 2000. While the DICOM Standard allows JPEG 2000 to be either lossy or lossless, when used for compression in pydicom it’s always treated as being lossy in order to simplify its usage.
Lossless compression
As with RLE and JPEG-LS, performing lossless compression is straightforward:
>>> from pydicom import examples
>>> from pydicom.uid import JPEG2000Lossless
>>> ds = examples.ct
>>> ds.compress(JPEG2000Lossless)
For RGB pixel data, JPEG 2000 can perform multiple component transformation (MCT) during the encoding process, which should improve the compression efficiency. This can be enabled or disabled by setting an appropriate Photometric Interpretation prior to compression:
"RGB"
to disable MCT"YBR_RCT"
to enable MCT for JPEG 2000 Lossless"YBR_ICT"
to enable MCT for JPEG 2000
>>> from pydicom import examples
>>> from pydicom.uid import JPEG2000Lossless
>>> ds = examples.rgb_color
>>> ds.PhotometricInterpretation
"RGB"
>>> ds.compress(JPEG2000Lossless) # No MCT applied
>>> len(ds.PixelData)
334412
>>> ds = examples.rgb_color
>>> ds.PhotometricInterpretation = "YBR_RCT"
>>> ds.compress(JPEG2000Lossless) # MCT applied
>>> len(ds.PixelData)
152342
Lossy compression
Lossy compression with JPEG 2000 is both more and less complicated then JPEG-LS; you don’t have to worry about the pixel values for signed integers, but specifying the image quality is less intuitive.
Warning
pydicom makes no recommendations for specifying the image quality for lossy encoding methods. Any examples of lossy encoding are for illustration purposes only.
>>> from pydicom import examples
>>> from pydicom.uid import JPEG2000
>>> ds = examples.ct
>>> ds.compress(JPEG2000, j2k_cr=[5, 2]) # 2 quality layers
With JPEG 2000 image quality is specified with either the j2k_cr or j2k_psnr parameters:
j2k_cr is a
list[float]
of compression ratios to use for each quality layer and is directly related to OpenJPEG’s -r compression ratio option. There must be at least one layer and the minimum allowable compression ratio is1
. When using multiple layers they should be ordered in decreasing value from left to right.j2k_psnr is a
list[float]
of the peak signal-to-noise ratios (in dB) to use for each quality layer and is directly related to OpenJPEG’s -q quality option. There must be at least one layer and when using multiple layers they should be ordered in increasing value from left to right.
Choosing appropriate quality settings for JPEG 2000 is far beyond the scope of this tutorial, but whatever you end up selecting should be thoroughly tested with a representative sample of expected pixel data.
More information on performing compression with JPEG 2000 can be found in the JPEG 2000 encoding guide.
Encapsulating data compressed by third-party packages¶
You can also use pydicom with third-party compression packages to encapsulate
the compressed Pixel Data, provided they meet the requirements of the
corresponding transfer syntax. The encapsulate()
or
encapsulate_extended()
functions are used to encapsulate the
compressed data.
from pydicom import examples
from pydicom.encaps import encapsulate, encapsulate_extended
from pydicom.uid import JPEGBaseline8Bit
# Fetch an example dataset
ds = examples.ct
# Use third-party package to compress
# Let's assume it compresses to JPEG Baseline
frames: list[bytes] = third_party_compression_func(...)
# Set the *Transfer Syntax UID* appropriately
ds.file_meta.TransferSyntaxUID = JPEGBaseline8Bit
# For *Samples per Pixel* 1 the *Photometric Interpretation* is unchanged
# Basic encapsulation
ds.PixelData = encapsulate(frames)
ds["PixelData"].VR = "OB" # always for encapsulated pixel data
ds.save_as("ct_compressed_basic.dcm")
# Extended encapsulation
result: tuple[bytes, bytes, bytes] = encapsulate_extended(frames)
ds.PixelData = result[0]
ds.ExtendedOffsetTable = result[1]
ds.ExtendedOffsetTableLength = result[2]
ds.save_as("ct_compressed_ext.dcm")
Decompression of Pixel Data¶
Datasets with a compressed Transfer Syntax UID can be decompressed with
Dataset.decompress()
or the
decompress()
function.
>>> from pydicom import examples
>>> ds = examples.jpeg2k
>>> ds.decompress()
This will decompress the Dataset
in-place:
The Pixel Data will be set using the uncompressed pixel data.
The Transfer Syntax UID will be changed to Explicit VR Little Endian.
The Image Pixel module elements will be updated as required to match the uncompressed pixel data.
A new SOP Instance UID value will be also be generated, but this can be disabled by passing
generate_instance_uid=False
.
Dataset decompression uses the same backend as accessing compressed Pixel Data, so the same customization options of the decoding process apply. For example, to use a specific plugin you can pass its name via the decoding_plugin argument:
>>> from pydicom import examples
>>> ds = examples.jpeg2k
>>> ds.decompress(decoding_plugin="pylibjpeg")
If the dataset’s Pixel Data is in the YCbCr color space it will also be converted
to RGB by default. This can be disabled by passing as_rgb=False
:
import numpy as np
from pydicom import examples
from pydicom.pixels import convert_color_space, pixel_array
from pydicom.uid import JPEG2000Lossless
# Original dataset in RGB
ds = examples.rgb_color
assert ds.PhotometricInterpretation == "RGB"
# Convert to YCbCr and compress
ybr = convert_color_space(ds.pixel_array, "RGB", "YBR_FULL")
ds.PhotometricInterpretation = "YBR_FULL"
ds.compress(JPEG2000Lossless, ybr)
assert ds.PhotometricInterpretation == "YBR_FULL"
# RGB reference - needed because converting RGB -> YBR -> RGB is lossy
rgb = convert_color_space(ybr, "YBR_FULL", "RGB")
# Decompress with conversion to RGB
ds.decompress()
assert ds.PhotometricInterpretation == "RGB"
assert np.array_equal(rgb, pixel_array(ds, raw=True))
# Decompress without conversion to RGB
ds.PhotometricInterpretation = "YBR_FULL"
ds.compress(JPEG2000Lossless, ybr)
ds.decompress(as_rgb=False)
assert ds.PhotometricInterpretation == "YBR_FULL"
assert np.array_equal(ybr, pixel_array(ds, raw=True))
Conclusion¶
In part 3 of this tutorial you’ve learned how to use pydicom to compress and decompress datasets and how to encapsulate pixel data that has been compressed by third-party packages. Having made it to the end of the pixel data tutorial you should now be comfortable using pydicom to perform pixel data related tasks.