diff --git a/image_analysis/azint/AzIntEngineGPU.cu b/image_analysis/azint/AzIntEngineGPU.cu index 52879eee..dce1d43c 100644 --- a/image_analysis/azint/AzIntEngineGPU.cu +++ b/image_analysis/azint/AzIntEngineGPU.cu @@ -91,8 +91,6 @@ void gpu_azim( AzIntEngineGPU::AzIntEngineGPU(const AzimuthalIntegrationMapping &integration, std::shared_ptr stream) : AzIntEngine(integration), stream(stream), - gpu_azint_correction(npixel), - gpu_pixel_to_bin(npixel), gpu_sum(azint_bins), gpu_sum2(azint_bins), gpu_count(azint_bins), @@ -100,22 +98,22 @@ AzIntEngineGPU::AzIntEngineGPU(const AzimuthalIntegrationMapping &integration, s cpu_sum2_reg(azint_sum2), cpu_count_reg(azint_count) { + int device = 0; + cuda_err(cudaGetDevice(&device)); // this worker's GPU, not necessarily 0 cudaDeviceProp prop{}; - cuda_err(cudaGetDeviceProperties(&prop, 0)); + cuda_err(cudaGetDeviceProperties(&prop, device)); threads = 128; blocks = 4 * prop.multiProcessorCount; shared_size = prop.sharedMemPerBlock; shared_needed = azint_bins * (2 * sizeof(float) + sizeof(uint32_t)); - // On this engine's stream, like every other operation it issues: the streams are non-blocking, so a - // NULL-stream copy is no longer ordered against the kernels that read what it uploads. The one-time - // synchronise leaves the constructor with the uploads settled rather than in flight. - cuda_err(cudaMemcpyAsync(gpu_azint_correction, integration.Corrections().data(), sizeof(float) * npixel, - cudaMemcpyHostToDevice, *stream)); - cuda_err(cudaMemcpyAsync(gpu_pixel_to_bin, integration.GetPixelToBin().data(), sizeof(uint16_t) * npixel, - cudaMemcpyHostToDevice, *stream)); - cuda_err(cudaStreamSynchronize(*stream)); + // Geometry-only, so shared per GPU: the first engine on this device uploads them, the rest reuse + // them. Keyed by the mapping's own vectors, which outlive every engine built from it. + gpu_azint_correction = SharedDeviceTable(integration.Corrections().data(), npixel, + integration.Corrections().data(), *stream); + gpu_pixel_to_bin = SharedDeviceTable(integration.GetPixelToBin().data(), npixel, + integration.GetPixelToBin().data(), *stream); } void AzIntEngineGPU::Run(const ImagePreprocessorBuffer &image, AzimuthalIntegrationProfile &profile) { @@ -127,12 +125,12 @@ void AzIntEngineGPU::Run(const ImagePreprocessorBuffer &image, AzimuthalIntegrat if (shared_needed < shared_size) { gpu_azim_shared<<>>( - gpu_pixel_to_bin,gpu_azint_correction,image.getGPUBuffer(), gpu_sum, gpu_sum2, + gpu_pixel_to_bin->get(),gpu_azint_correction->get(),image.getGPUBuffer(), gpu_sum, gpu_sum2, gpu_count, npixel, azint_bins ); } else { gpu_azim<<>>( - gpu_pixel_to_bin,gpu_azint_correction,image.getGPUBuffer(), gpu_sum, gpu_sum2, + gpu_pixel_to_bin->get(),gpu_azint_correction->get(),image.getGPUBuffer(), gpu_sum, gpu_sum2, gpu_count, npixel, azint_bins ); } diff --git a/image_analysis/azint/AzIntEngineGPU.h b/image_analysis/azint/AzIntEngineGPU.h index 08394373..33932a13 100644 --- a/image_analysis/azint/AzIntEngineGPU.h +++ b/image_analysis/azint/AzIntEngineGPU.h @@ -5,6 +5,7 @@ #include "AzIntEngine.h" #include "../indexing/CUDAMemHelpers.h" +#include "../indexing/CudaSharedTables.h" class AzIntEngineGPU : public AzIntEngine { std::shared_ptr stream; @@ -13,8 +14,10 @@ class AzIntEngineGPU : public AzIntEngine { size_t shared_needed; size_t shared_size; - CudaDevicePtr gpu_azint_correction; - CudaDevicePtr gpu_pixel_to_bin; + // Geometry-only tables: one copy per GPU, shared with every other engine on it (see + // CudaSharedTables.h) instead of one copy per worker thread. + std::shared_ptr> gpu_azint_correction; + std::shared_ptr> gpu_pixel_to_bin; CudaDevicePtr gpu_sum; CudaDevicePtr gpu_sum2; diff --git a/image_analysis/image_preprocessing/ImagePreprocessorGPU.cu b/image_analysis/image_preprocessing/ImagePreprocessorGPU.cu index b8180cbd..8844134d 100644 --- a/image_analysis/image_preprocessing/ImagePreprocessorGPU.cu +++ b/image_analysis/image_preprocessing/ImagePreprocessorGPU.cu @@ -93,25 +93,23 @@ ImagePreprocessorGPU::ImagePreprocessorGPU(const DiffractionExperiment &experime std::shared_ptr stream) : ImagePreprocessor(experiment), stream(stream), - gpu_mask(npixels), gpu_decompressed_image(npixels * sizeof(uint32_t)), // Overshoot - if input image is 1- or 2-byte, then it is still fine, while memory loss is minimal gpu_stats(1), cpu_stats(1), cpu_stats_reg(cpu_stats) { - // Setup mask + // Setup mask. The same for every worker, so it is uploaded once per GPU and shared; keyed on the + // PixelMask's own vector, which the derived table is a pure function of. std::vector mask_vec(npixels); for (int i = 0; i < npixels; i++) mask_vec[i] = (mask.GetMask().at(i) != 0); + gpu_mask = SharedDeviceTable(mask.GetMask().data(), npixels, mask_vec.data(), *stream); - // On this engine's stream, like every other operation it issues: the streams are non-blocking, so a - // NULL-stream copy is no longer ordered against the kernels that read the mask. Synchronise before - // leaving the constructor - mask_vec is a local and the copy must not outlive it. - cudaMemcpyAsync(gpu_mask, mask_vec.data(), npixels, cudaMemcpyHostToDevice, *stream); - cudaStreamSynchronize(*stream); - - // Setup GPU settings + // Setup GPU settings. The current device, not device 0: workers are pinned round-robin across GPUs, + // so device 0's SM count can belong to a different card than the one these kernels launch on. + int device = 0; + cudaGetDevice(&device); cudaDeviceProp prop{}; - cudaGetDeviceProperties(&prop, 0); + cudaGetDeviceProperties(&prop, device); threads = 128; blocks = 4 * prop.multiProcessorCount; @@ -154,7 +152,7 @@ ImageStatistics ImagePreprocessorGPU::Analyze(ImagePreprocessorBuffer &processed cudaMemcpyAsync(gpu_stats, cpu_stats.data(), sizeof(ImageStatistics), cudaMemcpyHostToDevice, *stream); preprocess_kernel <<< blocks, threads, 0, *stream >>>( reinterpret_cast(gpu_decompressed_image.get()), - gpu_mask, + gpu_mask->get(), processed_image.getGPUBuffer(), gpu_stats, sat_value, diff --git a/image_analysis/image_preprocessing/ImagePreprocessorGPU.h b/image_analysis/image_preprocessing/ImagePreprocessorGPU.h index 652c429d..3d4d3b92 100644 --- a/image_analysis/image_preprocessing/ImagePreprocessorGPU.h +++ b/image_analysis/image_preprocessing/ImagePreprocessorGPU.h @@ -5,12 +5,14 @@ #include "ImagePreprocessor.h" #include "../indexing/CUDAMemHelpers.h" +#include "../indexing/CudaSharedTables.h" class ImagePreprocessorGPU : public ImagePreprocessor { std::shared_ptr stream; int threads; int blocks; - CudaDevicePtr gpu_mask; + // Geometry-only, so one copy per GPU shared with every other engine on it (CudaSharedTables.h). + std::shared_ptr> gpu_mask; CudaDevicePtr gpu_decompressed_image; CudaDevicePtr gpu_stats; diff --git a/image_analysis/indexing/CudaSharedTables.h b/image_analysis/indexing/CudaSharedTables.h new file mode 100644 index 00000000..57c8d755 --- /dev/null +++ b/image_analysis/indexing/CudaSharedTables.h @@ -0,0 +1,77 @@ +// SPDX-FileCopyrightText: 2026 Filip Leonarski, Paul Scherrer Institute +// SPDX-License-Identifier: GPL-3.0-only + +#pragma once + +#include +#include +#include +#include + +#include "CUDAMemHelpers.h" + +// Read-only lookup tables that depend only on the detector geometry (pixel -> azimuthal bin, the +// per-pixel correction factors, the pixel mask). One analysis engine is built per worker thread, so +// each of those used to upload its own copy: on a 18 Mpx detector that is ~220 MB per thread, and +// 32 threads spent ~7 GB of device memory on identical data. +// +// Upload once per GPU instead and hand every engine on that GPU a shared pointer to the same table. +// The cache is keyed by (device, key) because a worker thread is pinned round-robin to a device +// (pin_gpu()), so on a multi-GPU node each device keeps its own copy - a kernel may only read memory +// resident on the device it runs on. `key` identifies the table's source data; use the address of the +// host vector that produced it, which lives in the experiment / integration mapping and therefore +// outlives every engine. +// +// Entries are held weakly, so the tables are released once the last engine using them is gone. + +namespace jfjoch_cuda_shared_tables { + struct Registry { + std::mutex m; + std::map, std::weak_ptr> tables; + }; + + inline Registry ®istry() { + static Registry r; + return r; + } + + // Not called cuda_err: the .cu files that include this header define their own such helper in an + // anonymous namespace, and a second one at global scope would make every call ambiguous. + inline void check(cudaError_t val) { + if (val != cudaSuccess) + throw JFJochException(JFJochExceptionCategory::GPUCUDAError, cudaGetErrorString(val)); + } +} + +// Return the device-resident copy of `host` (`count` elements) for the calling thread's GPU, +// uploading it on `stream` the first time it is asked for. +template +std::shared_ptr> SharedDeviceTable(const void *key, size_t count, const T *host, + cudaStream_t stream) { + int device = 0; + jfjoch_cuda_shared_tables::check(cudaGetDevice(&device)); + + auto ® = jfjoch_cuda_shared_tables::registry(); + // The upload happens while the lock is held: another worker must not obtain the pointer before + // its content is on the device. + std::lock_guard lock(reg.m); + auto &slot = reg.tables[{device, key}]; + if (auto cached = slot.lock()) + return std::static_pointer_cast>(cached); + + // Free on the device that allocated it - the last engine to drop the table may well be a worker + // pinned to a different GPU. + std::shared_ptr> table(new CudaDevicePtr(count), [device](CudaDevicePtr *p) { + int current = 0; + cudaGetDevice(¤t); + cudaSetDevice(device); + delete p; + cudaSetDevice(current); + }); + jfjoch_cuda_shared_tables::check( + cudaMemcpyAsync(table->get(), host, count * sizeof(T), cudaMemcpyHostToDevice, stream)); + jfjoch_cuda_shared_tables::check(cudaStreamSynchronize(stream)); + + slot = std::shared_ptr(table); + return table; +} diff --git a/image_analysis/spot_finding/AdaptiveSpotFinderGPU.cu b/image_analysis/spot_finding/AdaptiveSpotFinderGPU.cu index 262e7446..78379d2b 100644 --- a/image_analysis/spot_finding/AdaptiveSpotFinderGPU.cu +++ b/image_analysis/spot_finding/AdaptiveSpotFinderGPU.cu @@ -170,8 +170,6 @@ AdaptiveSpotFinderGPU::AdaptiveSpotFinderGPU(const AzimuthalIntegrationMapping & stream(std::move(in_stream)), nbins(in_mapping.GetBinNumber()), npix(in_mapping.GetPixelToBin().size()), - gpu_pixel_to_bin(npix), - gpu_corrections(npix), gpu_sum(nbins), gpu_sum2(nbins), gpu_count(nbins), @@ -203,14 +201,12 @@ AdaptiveSpotFinderGPU::AdaptiveSpotFinderGPU(const AzimuthalIntegrationMapping & shared_clip = static_cast(nbins) * (2 * sizeof(float) + sizeof(uint32_t)); use_shared = (shared_plain < prop.sharedMemPerBlock); - // On this engine's stream, like every other operation it issues: the streams are non-blocking, so a - // NULL-stream copy is no longer ordered against the kernels that read what it uploads. The one-time - // synchronise leaves the constructor with the uploads settled rather than in flight. - cuda_err(cudaMemcpyAsync(gpu_pixel_to_bin, mapping.GetPixelToBin().data(), sizeof(uint16_t) * npix, - cudaMemcpyHostToDevice, *stream)); - cuda_err(cudaMemcpyAsync(gpu_corrections, mapping.Corrections().data(), sizeof(float) * npix, - cudaMemcpyHostToDevice, *stream)); - cuda_err(cudaStreamSynchronize(*stream)); + // Both tables are functions of the detector geometry alone, so they are uploaded once per GPU and + // shared: the azimuthal-integration engine in the same worker reads the very same two arrays. + gpu_pixel_to_bin = SharedDeviceTable(mapping.GetPixelToBin().data(), npix, + mapping.GetPixelToBin().data(), *stream); + gpu_corrections = SharedDeviceTable(mapping.Corrections().data(), npix, + mapping.Corrections().data(), *stream); } void AdaptiveSpotFinderGPU::ReducePass(const ImagePreprocessorBuffer &image, float clip_k, @@ -218,12 +214,12 @@ void AdaptiveSpotFinderGPU::ReducePass(const ImagePreprocessorBuffer &image, flo if (use_shared) { const size_t shared = accumulate_corrected ? shared_plain : shared_clip; reduce_rings_shared<<>>( - gpu_pixel_to_bin, gpu_corrections, image.getGPUBuffer(), gpu_mean, gpu_sigma, + gpu_pixel_to_bin->get(), gpu_corrections->get(), image.getGPUBuffer(), gpu_mean, gpu_sigma, clip_k, accumulate_corrected, gpu_sum, gpu_sum2, gpu_count, gpu_sum_corr, gpu_sum2_corr, npix, nbins); } else { reduce_rings_global<<>>( - gpu_pixel_to_bin, gpu_corrections, image.getGPUBuffer(), gpu_mean, gpu_sigma, + gpu_pixel_to_bin->get(), gpu_corrections->get(), image.getGPUBuffer(), gpu_mean, gpu_sigma, clip_k, accumulate_corrected, gpu_sum, gpu_sum2, gpu_count, gpu_sum_corr, gpu_sum2_corr, npix, nbins); } @@ -327,7 +323,7 @@ void AdaptiveSpotFinderGPU::Detect(const ImagePreprocessorBuffer &image, cuda_err(cudaMemcpyAsync(gpu_thr, host_thr.data(), sizeof(float) * nbins, cudaMemcpyHostToDevice, *stream)); cuda_err(cudaMemsetAsync(gpu_strong, 0, OutputByteSize(), *stream)); flag_strong<<>>( - image.getGPUBuffer(), gpu_pixel_to_bin, gpu_thr, gpu_strong, npix, nbins); + image.getGPUBuffer(), gpu_pixel_to_bin->get(), gpu_thr, gpu_strong, npix, nbins); cuda_err(cudaMemcpyAsync(output_buffer.data(), gpu_strong, OutputByteSize(), cudaMemcpyDeviceToHost, *stream)); cuda_err(cudaStreamSynchronize(*stream)); } diff --git a/image_analysis/spot_finding/AdaptiveSpotFinderGPU.h b/image_analysis/spot_finding/AdaptiveSpotFinderGPU.h index 9481de42..999327ba 100644 --- a/image_analysis/spot_finding/AdaptiveSpotFinderGPU.h +++ b/image_analysis/spot_finding/AdaptiveSpotFinderGPU.h @@ -30,6 +30,7 @@ #include "../../common/AzimuthalIntegrationProfile.h" #include "../../common/AzimuthalIntegrationMapping.h" #include "../indexing/CUDAMemHelpers.h" +#include "../indexing/CudaSharedTables.h" class AdaptiveSpotFinderGPU : public ImageSpotFinder { const AzimuthalIntegrationMapping &mapping; @@ -46,9 +47,10 @@ class AdaptiveSpotFinderGPU : public ImageSpotFinder { size_t shared_clip = 0; // per-block shared bytes for a sigma-clip pass (raw rings only) bool use_shared = true; // false -> nbins too large for shared memory, use the global-atomics kernel - // Static mapping inputs (uploaded once). - CudaDevicePtr gpu_pixel_to_bin; - CudaDevicePtr gpu_corrections; + // Static mapping inputs: geometry-only, so one copy per GPU shared with every other engine on it + // (see CudaSharedTables.h) rather than one copy per worker thread. + std::shared_ptr> gpu_pixel_to_bin; + std::shared_ptr> gpu_corrections; // Raw per-ring accumulators (re-zeroed each pass) + derived stats used to clip and threshold. // double, like the CPU engine's ring accumulators: the ring sigma is the cancelling difference