WIP: read volume (3d image) data

Limitations:
- does not support compressed/encapsulated pixel data
(see non-raw dicom transfer syntaxes, ex: JPEG2000)
- assumes slicethickness is spacing
- harcode pixel type and dimension in imageType

Author: Alexis Girault

examples/Dicom/src/index.js

@@ -17,12 +17,21 @@ const outputFileInformation = curry(async function outputFileInformation (output
   // Select DICOM serie
   outputTextArea.textContent = "Please select serie..."
   setupDicomForm(patients, async (serie) => {
+    console.time('customRead:')
+    const image1 = serie.getImageData()
+    console.log(image1)
+    console.warn(image1.data.length)
+    console.timeEnd('customRead:')
     outputTextArea.textContent = "Loading..."
 
     // Read DICOM serie
+    console.time('itkRead:')
     const files = Object.values(serie.images).map((image) => image.file)
     const { image, webWorker } = await readImageDICOMFileSeries(null, files)
     webWorker.terminate()
+    console.log(image)
+    console.warn(image.data.length)
+    console.timeEnd('itkRead:')
 
     // Display
     function replacer (key, value) {

examples/Dicom/src/parseDicomFiles.js

@@ -5,6 +5,18 @@ import "regenerator-runtime/runtime";
 
 import DICOM_TAG_DICT from './dicomTags'
 
+function concatenate(resultConstructor, arrays) {
+  const totalLength = arrays.reduce((total, arr) => {
+    return total + arr.length
+  }, 0);
+  const result = new resultConstructor(totalLength);
+  arrays.reduce((offset, arr) => {
+    result.set(arr, offset);
+    return offset + arr.length;
+  }, 0);
+  return result;
+}
+
 class DICOMEntity {
   constructor() {
     this.metaData = {}
@@ -120,6 +132,140 @@ class DICOMSeries extends DICOMEntity {
     }
     this.images[imageNumber] = new DICOMImage(metaData, file)
   }
+
+  getImageData() {
+    function numArrayFromString(str, separator = '\\') {
+      const strArray = str.split(separator)
+      return strArray.map(Number)
+    }
+
+    const slices = Object.values(this.images)
+    const meta = slices[0].metaData
+
+    // Origin
+    const origin = numArrayFromString(meta.ImagePositionPatient)
+
+    // Spacing
+    const spacing = numArrayFromString(meta.PixelSpacing)
+    spacing.push(Number(meta.SliceThickness)) // TODO: or SpacingBetweenSlices?
+
+    // Dimensions
+    const size = [
+      meta.Rows,
+      meta.Columns,
+      Object.keys(this.images).length,
+    ]
+
+    // Direction matrix (3x3)
+    const directionCosines = numArrayFromString(meta.ImageOrientationPatient)
+    const iDirCos = directionCosines.slice(0, 3)
+    const jDirCos = directionCosines.slice(3, 6)
+    const kDirCos = [
+      iDirCos[1] * jDirCos[2] - iDirCos[2] * jDirCos[1],
+      iDirCos[2] * jDirCos[0] - iDirCos[0] * jDirCos[2],
+      iDirCos[0] * jDirCos[1] - iDirCos[1] * jDirCos[0],
+    ]
+    const direction = {
+      data: [
+        iDirCos[0], jDirCos[0], kDirCos[0],
+        iDirCos[1], jDirCos[1], kDirCos[1],
+        iDirCos[2], jDirCos[2], kDirCos[2],
+      ],
+    }
+
+    // Pixel data type
+    const unsigned = (meta.PixelRepresentation === 0) && (meta.RescaleIntercept > 0)
+    const bits = meta.BitsAllocated
+    let ArrayType
+    let intType
+    switch (bits) {
+      case 8:
+        ArrayType = unsigned ? Uint8Array : Int8Array
+        intType = unsigned ? 'uint8_t' : 'int8_t'
+        break
+      case 16:
+        ArrayType = unsigned ? Uint16Array : Int16Array
+        intType = unsigned ? 'uint16_t' : 'int16_t'
+        break
+      case 32:
+        ArrayType = unsigned ? Uint32Array : Int32Array
+        intType = unsigned ? 'uint32_t' : 'int32_t'
+        break
+      default:
+        throw Error(`Unknown pixel bit type (${bits})`)
+    }
+
+    // Image info
+    const imageType = {
+      dimension: 3,
+      componentType: intType,
+      pixelType: 1, // TODO: based on meta.PhotometricInterpretation?
+      components: meta.SamplesPerPixel,
+    }
+
+    // Dataview on pixel data
+    const pixelDataArrays = slices.map((image) => {
+      const value = image.metaData.PixelData
+      if (value.buffer.constructor === ArrayType && value.offset === 0) {
+        return value.buffer
+      }
+      return new ArrayType(value.buffer, value.offset)
+    })
+
+    // Concatenate all pixel data
+    const data = pixelDataArrays.length === 1
+      ? pixelDataArrays[0]
+      : concatenate(ArrayType, pixelDataArrays)
+
+    // Masking bits
+    let maskFunction
+    if (meta.BitsStored !== bits) {
+      let mask = ''
+      for (let i = 0; i < bits; i += 1) {
+        if (i < meta.HighBit - meta.BitsStored || i > meta.HighBit) {
+          mask += '0'
+        } else {
+          mask += '1'
+        }
+      }
+      maskFunction = (value, index) => { data[index] = value & mask }
+    }
+
+    // Rescale
+    const b = Number(meta.RescaleIntercept)
+    const m = Number(meta.RescaleSlope)
+    const hasIntercept = !Number.isNaN(b) && b !== 0
+    const hasSlope = !Number.isNaN(m) && m !== 1
+    let rescaleFunction
+    if (hasIntercept && hasSlope) {
+      rescaleFunction = (value, index) => { data[index] = m * value + b }
+    } else if (hasIntercept) {
+      rescaleFunction = (value, index) => { data[index] = value + b }
+    } else if (hasSlope) {
+      rescaleFunction = (value, index) => { data[index] = m * value }
+    }
+
+    // Apply transformations if needed
+    if (maskFunction && rescaleFunction) {
+      data.forEach((_, index) => {
+        maskFunction(data[index], index)
+        rescaleFunction(data[index], index)
+      })
+    } else if (maskFunction) {
+      data.forEach(maskFunction)
+    } else if (rescaleFunction) {
+      data.forEach(rescaleFunction)
+    }
+
+    return {
+      imageType,
+      origin,
+      spacing,
+      direction,
+      size,
+      data,
+    }
+  }
 }
 
 class DICOMImage extends DICOMEntity {
@@ -147,6 +293,9 @@ class DICOMImage extends DICOMEntity {
       'BitsStored',
       'HighBit',
       'PixelRepresentation',
+      'PixelData',
+      'RescaleIntercept',
+      'RescaleSlope',
       ]
     }
 
@@ -223,60 +372,86 @@ async function parseDicomFiles(fileList, ignoreFailedFiles = false) {
           return
         }
 
-        if (element.fragments) {
-          console.warn(`${tagName} contains fragments which isn't supported`)
-          return
-        }
+        let value = undefined
+
+        if (tagName === 'PixelData') {
+          if (element.fragments) {
+            let bot = element.basicOffsetTable
+            // if basic offset table is empty, calculate it
+            if (bot.length === 0) {
+              bot = dicomParser.createJPEGBasicOffsetTable(dataSet, element)
+            }
 
-        let vr = element.vr
-        if (vr === undefined) {
-          if (tagInfo === undefined || tagInfo.vr === undefined) {
-            console.warn(`${tagName} vr is unknown, skipping`)
+            const imageFrames = []
+            for (let frameIndex = 0; frameIndex < bot.length; frameIndex += 1) {
+              imageFrames.push(dicomParser.readEncapsulatedImageFrame(dataSet, element, 0, bot))
+            }
+            const buffer = imageFrames.length === 1
+              ? imageFrames[0]
+              : concatenate(imageFrames[0].constructor, imageFrames)
+            value = {
+              buffer,
+              offset: 0,
+              length: buffer.length,
+              encapsulated: true,
+            }
+          } else {
+            value = {
+              buffer: dataSet.byteArray.buffer,
+              offset: element.dataOffset,
+              length: element.length,
+              encapsulated: false,
+            }
+          }
+        } else {
+          let vr = element.vr
+          if (vr === undefined) {
+            if (tagInfo === undefined || tagInfo.vr === undefined) {
+              console.warn(`${tagName} vr is unknown, skipping`)
+            }
+            vr = tagInfo.vr
           }
-          vr = tagInfo.vr
-        }
 
-        let value = undefined
-        switch (vr) {
-          case 'US':
-            value = dataSet.uint16(tag)
-            break
-          case 'SS':
-            value = dataSet.int16(tag)
-            break
-          case 'UL':
-            value = dataSet.uint32(tag)
-            break
-          case 'US':
-            value = dataSet.int32(tag)
-            break
-          case 'FD':
-            value = dataSet.double(tag)
-            break
-          case 'FL':
-            value = dataSet.float(tag)
-            break
-          case 'AT':
-            value = `(${dataSet.uint16(tag, 0)},${dataSet.uint16(tag, 1)})`
-            break
-          case 'OB':
-          case 'OW':
-          case 'UN':
-          case 'OF':
-          case 'UT':
-            // TODO: binary data? is this correct?
-            if (element.length === 2) {
+          switch (vr) {
+            case 'US':
               value = dataSet.uint16(tag)
-            } else if (element.length === 4) {
+              break
+            case 'SS':
+              value = dataSet.int16(tag)
+              break
+            case 'UL':
               value = dataSet.uint32(tag)
-            } else {
-              // don't store binary data, only meta data
-              return
-            }
-            break
-          default: //string
-            value = dataSet.string(tag)
-            break
+              break
+            case 'US':
+              value = dataSet.int32(tag)
+              break
+            case 'FD':
+              value = dataSet.double(tag)
+              break
+            case 'FL':
+              value = dataSet.float(tag)
+              break
+            case 'AT':
+              value = `(${dataSet.uint16(tag, 0)},${dataSet.uint16(tag, 1)})`
+              break
+            case 'OB':
+            case 'OW':
+            case 'UN':
+            case 'OF':
+            case 'UT':
+              // TODO: binary data? is this correct?
+              if (element.length === 2) {
+                value = dataSet.uint16(tag)
+              } else if (element.length === 4) {
+                value = dataSet.uint32(tag)
+              } else {
+                return
+              }
+              break
+            default: //string
+              value = dataSet.string(tag)
+              break
+          }
         }
 
         metaData[tagName] = value