packages = ["pillow", "numpy", "pydicom", "pynrrd", "scipy"] Download DICOM Convert Software

Volume Slicer

Loading file(s)...
1.00x
0
0
0
from js import document from pyodide.ffi.wrappers import add_event_listener import numpy as np from PIL import Image from io import BytesIO import base64 import math from js import window import pydicom import nrrd from pyodide.http import pyfetch import tempfile import os # Global variable to store the loaded volume. volume = None volume_dims = None file_type = None # To track if we loaded NPZ, NPY, NRRD, or DICOM def extract_slice(): """ Extract a 2D slice from the 3D volume using the current slider values. The volume is assumed to be a 3D or 4D array with shape (nZ, nY, nX, channels) where the coordinates correspond to (x, y, z) with volume indices accessed as volume[z, y, x]. """ global volume, volume_dims if volume is None: return # Get slider values. z_pos = int(document.getElementById("z-slider").value) x_pos = int(document.getElementById("x-slider").value) y_pos = int(document.getElementById("y-slider").value) angle_x = float(document.getElementById("angle-x-slider").value) angle_y = float(document.getElementById("angle-y-slider").value) # Convert angles to radians. ax = math.radians(angle_x) ay = math.radians(angle_y) # Compute rotation matrices for x and y. Rx = np.array([ [1, 0, 0], [0, math.cos(ax), -math.sin(ax)], [0, math.sin(ax), math.cos(ax)], ]) Ry = np.array([ [math.cos(ay), 0, math.sin(ay)], [0, 1, 0], [-math.sin(ay), 0, math.cos(ay)], ]) # Compose rotations (first rotate about X, then Y). R = Ry.dot(Rx) # Determine the slicing plane. # The default slice has normal (0, 0, 1) with in-plane axes (1, 0, 0) and (0, 1, 0). normal = R.dot(np.array([0, 0, 1])) u_axis = R.dot(np.array([1, 0, 0])) v_axis = R.dot(np.array([0, 1, 0])) # Interpret volume shape. if volume.ndim == 3: nZ, nY, nX = volume.shape channels = 1 elif volume.ndim == 4: nZ, nY, nX = volume.shape[:3] channels = volume.shape[3] else: # For 2D DICOM files, we need to handle them differently nY, nX = volume.shape nZ = 1 channels = 1 volume = volume.reshape(1, nY, nX) # Define the center of the volume in (x, y, z) coordinates. center = np.array([(nX - 1) / 2.0, (nY - 1) / 2.0, (nZ - 1) / 2.0]) # Set output slice parameters. out_size = 256 # Resolution of the output slice (256 x 256) L = max(nX, nY, nZ) # Create a grid in plane coordinates (u, v) u = np.linspace(-L/2, L/2, out_size) v = np.linspace(-L/2, L/2, out_size) U, V = np.meshgrid(u, v, indexing="xy") # each of shape (out_size, out_size) # For each output pixel compute the corresponding (x, y, z) in volume: # p = center + (slice_offset * normal) + (u * u_axis) + (v * v_axis) P = ( np.array([x_pos - center[0], y_pos - center[1], z_pos - center[2]]).reshape(1, 1, 3) + center.reshape(1, 1, 3) + U[..., np.newaxis] * u_axis.reshape(1, 1, 3) + V[..., np.newaxis] * v_axis.reshape(1, 1, 3) ) # P has shape (out_size, out_size, 3) with (x, y, z) coordinates. X = P[..., 0] Y = P[..., 1] Z = P[..., 2] # Convert to volume indices; volume indices are accessed as volume[z, y, x]. x_idx = np.rint(X).astype(np.int32) y_idx = np.rint(Y).astype(np.int32) z_idx = np.rint(Z).astype(np.int32) # Create mask for valid indices valid_mask = ( (x_idx >= 0) & (x_idx < nX) & (y_idx >= 0) & (y_idx < nY) & (z_idx >= 0) & (z_idx < nZ) ) # Clip indices for safe array access x_idx_safe = np.clip(x_idx, 0, nX - 1) y_idx_safe = np.clip(y_idx, 0, nY - 1) z_idx_safe = np.clip(z_idx, 0, nZ - 1) # Sample the volume if channels == 1: slice_img = volume[z_idx_safe, y_idx_safe, x_idx_safe] slice_img = np.where(valid_mask, slice_img, 0) else: slice_img = volume[z_idx_safe, y_idx_safe, x_idx_safe, :] slice_img = np.where(valid_mask[..., np.newaxis], slice_img, 0) # Normalize to uint8 if needed. if slice_img.dtype != np.uint8: m = slice_img.min() M = slice_img.max() if M > m: slice_img = (255 * (slice_img - m) / (M - m)).astype(np.uint8) else: slice_img = slice_img.astype(np.uint8) # Convert the numpy slice to a Pillow image. if channels in [1, 3, 4]: pil_im = Image.fromarray(slice_img) else: pil_im = Image.fromarray(slice_img.astype(np.uint8)) # Convert the Pillow image to a PNG data URL. buffer = BytesIO() pil_im.save(buffer, format="PNG") encoded = base64.b64encode(buffer.getvalue()).decode("ascii") data_url = f"data:image/png;base64,{encoded}" # Update the output area with the new image. output_div = document.getElementById("output") output_div.innerHTML = "" img_elem = document.createElement("img") img_elem.src = data_url output_div.appendChild(img_elem) # Update slider display values. document.getElementById("z-value").innerText = f"{z_pos}" document.getElementById("x-value").innerText = f"{x_pos}" document.getElementById("y-value").innerText = f"{y_pos}" document.getElementById("angle-x-value").innerText = f"{angle_x:.1f}°" document.getElementById("angle-y-value").innerText = f"{angle_y:.1f}°" async def on_file_upload(e): """ When a file is uploaded, load it and update the volume. Supports NPZ, NPY, NRRD, and DICOM files. """ file_list = e.target.files if file_list.length == 0: return # Check what type of file we're dealing with first_file = file_list.item(0) file_ext = first_file.name.rsplit('.', 1)[-1].lower() if '.' in first_file.name else '' global volume, volume_dims, file_type try: # Show loading indicator document.getElementById("loading").style.display = "block" document.getElementById("error").style.display = "none" if file_ext == 'npz': # Handle NPZ files bytes_data = await get_bytes_from_file(first_file) npz_file = np.load(BytesIO(bytes_data)) keys = list(npz_file.files) if len(keys) == 0: raise ValueError("No arrays found in the NPZ file.") volume = npz_file[keys[0]] if volume.ndim != 4 and volume.ndim != 3: raise ValueError(f"Expected a 3D or 4D array. Got shape: {volume.shape}") if volume.ndim == 4: nZ, nY, nX = volume.shape[:3] else: nY, nX, nZ = volume.shape volume = volume.reshape(nY, nX, nZ) volume_dims = (nZ, nY, nX) file_type = "NPZ" elif file_ext == 'npy': # Handle NPY files bytes_data = await get_bytes_from_file(first_file) volume = np.load(BytesIO(bytes_data)) if volume.ndim != 4 and volume.ndim != 3: raise ValueError(f"Expected a 3D or 4D array. Got shape: {volume.shape}") if volume.ndim == 4: nZ, nY, nX = volume.shape[:3] else: nY, nX, nZ = volume.shape volume = volume.reshape(nY, nX, nZ) volume_dims = (nZ, nY, nX) file_type = "NPY" elif file_ext == 'nrrd': # Handle NRRD files bytes_data = await get_bytes_from_file(first_file) # Write bytes to a temporary file since nrrd.read() requires a file path temp_file = tempfile.NamedTemporaryFile( suffix='.nrrd', delete=False ) temp_path = temp_file.name temp_file.write(bytes_data) temp_file.close() try: data, header = nrrd.read(temp_path) # Rotate clockwise 90 degrees data = np.rot90(data, k=-1, axes=(0, 2)) data = np.flip(data, axis=2) # Stretch width (columns) to double from scipy import ndimage nZ, nY, nX = data.shape zoom_factors = (1, 1, 2) data = ndimage.zoom(data, zoom_factors, order=1) volume = data nZ, nY, nX = volume.shape volume_dims = (nZ, nY, nX) file_type = "NRRD" finally: # Clean up temp file os.unlink(temp_path) elif file_ext == 'dcm': # Handle DICOM files bytes_data = await get_bytes_from_file(first_file) dcm = pydicom.dcmread(BytesIO(bytes_data)) volume = dcm.pixel_array # Handle DICOM volume dimensions if volume.ndim == 4: nZ, nY, nX = volume.shape[:3] elif volume.ndim == 3: nZ, nY, nX = volume.shape else: # 2D DICOM nY, nX = volume.shape nZ = 1 volume = volume.reshape(1, nY, nX) volume_dims = (nZ, nY, nX) file_type = "DCM" else: raise ValueError(f"Unsupported file type: {file_ext}") # Save volume information to window window.volume_dims = volume_dims # Update the slice slider limits. z_slider = document.getElementById("z-slider") z_slider.min = 0 z_slider.max = volume_dims[0] - 1 z_slider.value = volume_dims[0] // 2 x_slider = document.getElementById("x-slider") x_slider.min = 0 x_slider.max = volume_dims[2] - 1 x_slider.value = volume_dims[2] // 2 y_slider = document.getElementById("y-slider") y_slider.min = 0 y_slider.max = volume_dims[1] - 1 y_slider.value = volume_dims[1] // 2 # Extract the initial slice extract_slice() # Hide loading indicator document.getElementById("loading").style.display = "none" except Exception as e: # Hide loading indicator and show error document.getElementById("loading").style.display = "none" error_div = document.getElementById("error") error_div.innerHTML = f"Error loading file: {str(e)}" error_div.style.display = "block" async def get_bytes_from_file(file): array_buf = await file.arrayBuffer() return array_buf.to_bytes() # Attach event listeners. add_event_listener( document.getElementById("file-upload"), "change", on_file_upload ) add_event_listener( document.getElementById("z-slider"), "input", lambda e: extract_slice() ) add_event_listener( document.getElementById("x-slider"), "input", lambda e: extract_slice() ) add_event_listener( document.getElementById("y-slider"), "input", lambda e: extract_slice() ) add_event_listener( document.getElementById("angle-x-slider"), "input", lambda e: extract_slice() ) add_event_listener( document.getElementById("angle-y-slider"), "input", lambda e: extract_slice() ) window.extract_slice = extract_slice