method for upscaling data with no hardware constraints. to be used for benchmarks.

.gitignore

 out/*
 logs/*
 __pycache__/
+temp/*
\ No newline at end of file

bjobs/exciscope.sh

+#!/bin/sh
+#BSUB -J exciscope
+#BSUB -q gpuv100
+#BSUB -n 8
+#BSUB -R "span[hosts=1]"
+#BSUB -R "rusage[mem=8GB]"
+#BSUB -W 1:00
+#BSUB -B
+#BSUB -o logs/bjobs/exciscope_%J.out
+#BSUB -e logs/bjobs/exciscope_%J.err
+
+nvidia-smi
+module load cuda/11.6
+/appl/cuda/11.6.0/samples/bin/x86_64/linux/release/deviceQuery
+
+#### Parameters ####
+ANGLE_CHUNKS=32
+VOL_CHUNKS_PER_AXIS=3
+####################
+
+. env_setup.sh
+
+VOL_CHUNKS=$VOL_CHUNKS_PER_AXIS,$VOL_CHUNKS_PER_AXIS,$VOL_CHUNKS_PER_AXIS
+
+DATA_PATH="/work3/s214743/ct_data/exciscope_out"
+TIFF_DIR="/dtu/3d-imaging-center/projects/2024_QIM_platform/analysis/data/ses-240826-143058-woodzoom/02-tomo/Output/Binaries"
+
+/usr/bin/time -v python -u ct/chunking.py \
+    --tiff-dir "$TIFF_DIR" \
+    --output-dir "$DATA_PATH" \
+    --angle-chunks $ANGLE_CHUNKS \
+    --vol-chunks $VOL_CHUNKS \
+    --read-batch-size 2 \
+    --use-pickle-geometries
+
+du -lh $DATA_PATH
+# rm -rf $DATA_PATH
\ No newline at end of file

ct/upscale_data.py

+import os
+import cv2
+from tqdm import tqdm
+import numpy as np
+import sys
+from cil.framework import ImageGeometry, AcquisitionGeometry, ImageData
+from ct.chunking import create_reader
+import pickle
+
+
+DOUBLING_STEPS = int(sys.argv[1]) # number of projections will increase by the factor 2**DOUBLING_STEPS
+# DOUBLING_STEPS = 2
+ORIG_SIZE = 1000
+TESTING = False # True or False
+
+
+INPUT_FOLDER = '/dtu/3d-imaging-center/projects/2021_DANFIX_Casper/raw_data_3DIM/Casper_top_3_2 [2021-03-17 16.54.39]/'
+CONFIG_PATH = '/dtu/3d-imaging-center/projects/2021_DANFIX_Casper/raw_data_3DIM/Casper_top_3_2 [2021-03-17 16.54.39]/Casper_top_3_2_recon.xtekct'
+
+def get_input_tiff_pairs(input_folder):
+    input_tiff_pairs = []
+    idx = 0
+    for i, filename in enumerate(sorted(os.listdir(input_folder))):
+        if filename.lower().endswith(('.tif', '.tiff')):
+            idx += 1
+            input_tiff_pairs.append(([idx, os.path.join(input_folder, filename)]))
+        
+        # small subset to test on when testing
+        if SUBSET:
+            if input_folder == INPUT_FOLDER and idx == 9:
+                break
+    
+    return input_tiff_pairs
+
+def upscale_image_sizes(input_folder, output_folder):
+    input_tiff_pairs = get_input_tiff_pairs(input_folder)
+    for idx, input_path in tqdm(input_tiff_pairs, desc="Upscaling TIFFs"):
+        output_path = os.path.join(output_folder, f'projection_{idx:04d}.tif')
+        img = cv2.imread(input_path, cv2.IMREAD_UNCHANGED)
+        img_upscale = cv2.resize(img, (NEW_SIZE, NEW_SIZE), interpolation=cv2.INTER_LINEAR)
+        cv2.imwrite(output_path, img_upscale)
+
+def double_tiffs(output_folder, doubling_step):
+    input_tiff_pairs = get_input_tiff_pairs(output_folder)
+    for idx, input_path in tqdm(input_tiff_pairs, desc=f"Doubling step {doubling_step}"):
+        img = cv2.imread(input_path, cv2.IMREAD_UNCHANGED)
+
+        if idx == len(input_tiff_pairs):
+            # the last and first image are not guaranteed to be close to each other.
+            # safer to just repeat the last one
+            img_new = img
+        else:
+            idx_next = idx + 1
+            img_next = cv2.imread(os.path.join(output_folder, f'projection_{idx_next:04d}.tif'), cv2.IMREAD_UNCHANGED)
+            img_new = ((img.astype(np.float32) + img_next.astype(np.float32))/2).astype(np.uint16)
+
+        output_path = os.path.join(output_folder, f'projection_{idx:04d}_new.tif') # works since ord('.') > ord('_')
+        cv2.imwrite(output_path, img_new)
+    
+    input_tiff_pairs = get_input_tiff_pairs(output_folder)
+    for idx, input_path in input_tiff_pairs:
+        os.rename(input_path, os.path.join(output_folder, f'temp_{idx:04d}.tif'))
+    
+    for idx, _ in input_tiff_pairs:
+        os.rename(os.path.join(output_folder, f'temp_{idx:04d}.tif'), os.path.join(output_folder, f'projection_{idx:04d}.tif'))
+
+def get_new_angles(angles):
+    angles_shift1 = np.take(angles, 1+np.arange(len(angles)), mode='clip')
+    angles_new = angles + (angles_shift1 - angles)/2
+    angles_combined = np.sort(np.concatenate((angles, angles_new)))[::-1] # reversed to follow the structure of casper data
+    return angles_combined
+
+ag = create_reader(CONFIG_PATH).get_geometry()
+
+for doubling_step in range(DOUBLING_STEPS):
+    if TESTING and (doubling_step == 0):
+        SUBSET = True
+    else:
+        SUBSET = False
+    
+    NEW_SIZE = 2**(doubling_step+1) * ORIG_SIZE
+    OUTPUT_FOLDER = f'/dtu/3d-imaging-center/projects/2024_QIM_platform/analysis/data/casper_{NEW_SIZE}'
+    os.makedirs(OUTPUT_FOLDER, exist_ok=True)
+
+    upscale_image_sizes(INPUT_FOLDER, OUTPUT_FOLDER)
+    double_tiffs(OUTPUT_FOLDER, doubling_step)
+    
+    ag.set_panel(
+        num_pixels=2*ag.config.panel.num_pixels,
+        pixel_size=0.5*ag.config.panel.pixel_size,
+        origin=ag.config.panel.origin
+    )
+    ag.set_angles(
+        angles=get_new_angles(ag.angles),
+        initial_angle=ag.config.angles.initial_angle,
+        angle_unit=ag.config.angles.angle_unit
+    )
+    
+    ig = ag.get_ImageGeometry()
+    
+    with open(os.path.join(OUTPUT_FOLDER, 'ag.pkl'), 'wb') as f:
+        pickle.dump(ag, f)
+
+    with open(os.path.join(OUTPUT_FOLDER, 'ig.pkl'), 'wb') as f:
+        pickle.dump(ig, f)
+    
+    INPUT_FOLDER = OUTPUT_FOLDER
+
+
+
+
+
+
+
+### OLD
+# for i, filename in enumerate(tqdm(sorted(os.listdir(input_folder)), desc="Upscaling TIFFs")):
+#     if filename.lower().endswith(('.tif', '.tiff')):
+#         input_path = os.path.join(input_folder, filename)
+#         output_path = os.path.join(output_folder, filename)
+
+#         img = cv2.imread(input_path, cv2.IMREAD_UNCHANGED)  # reads 16-bit as-is
+#         if img is None:
+#             print(f"Failed to read {filename}")
+#             continue
+
+#         resized = cv2.resize(img, (new_size, new_size), interpolation=cv2.INTER_LINEAR)
+#         cv2.imwrite(output_path, resized)
+
+# for i, filename in enumerate(tqdm(sorted(os.listdir(output_folder)), desc="Upscaling number of projections")):
+#     if filename.lower().endswith(('.tif', '.tiff')):
+#         output_path = os.path.join(output_folder, filename)
+
+#         img = cv2.imread(input_path, cv2.IMREAD_UNCHANGED)  # reads 16-bit as-is
+#         if img is None:
+#             print(f"Failed to read {filename}")
+#             continue
+
+#         resized = cv2.resize(img, (new_size, new_size), interpolation=cv2.INTER_LINEAR)
+#         cv2.imwrite(output_path, resized)
+

profiling/lp_chunking.py

@@ -33,9 +33,11 @@ def main():
     processor = ChunkedCT(
         config_path='/dtu/3d-imaging-center/projects/2021_DANFIX_Casper/raw_data_3DIM/Casper_top_3_2 [2021-03-17 16.54.39]/Casper_top_3_2_recon.xtekct',
         tiff_dir=None,
-        output_dir='out/chunking/lp',
+        # output_dir='out/chunking/lp',
+        output_dir='/work3/s214743/ct_data/lp',
         angle_chunks=2,
         vol_chunks=(2,2,2),
+        verbose=0,
     )
     processor.write_projection_chunks()
     processor.reconstruct()