Skip to content
Merged
Show file tree
Hide file tree
Changes from 1 commit
Commits
File filter

Filter by extension

Filter by extension


Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
8 changes: 4 additions & 4 deletions .devcontainer/Dockerfile
Original file line number Diff line number Diff line change
Expand Up @@ -43,7 +43,7 @@ ENTRYPOINT [ "/usr/local/share/docker-init.sh" ]
CMD [ "sleep", "infinity" ]

ARG MINIFORGE_NAME=Miniforge3
ARG MINIFORGE_VERSION=24.3.0-0
ARG MINIFORGE_VERSION=25.11.0-1
ARG CONDA_DIR=/opt/conda

# Based on https://github.com/conda-forge/miniforge-images/blob/master/ubuntu/Dockerfile
Expand All @@ -63,8 +63,8 @@ RUN wget --no-hsts --quiet https://github.com/conda-forge/miniforge/releases/dow
# Create a conda environment from the environment file in the repo root.
COPY --from=file-normalizer --chown=$USER_UID:conda /data/environment.yml /tmp/build/
RUN umask 0002 \
&& ${CONDA_DIR}/bin/mamba env create -f /tmp/build/environment.yml \
&& ${CONDA_DIR}/bin/mamba clean -fy \
&& ${CONDA_DIR}/bin/conda env create -f /tmp/build/environment.yml \
&& ${CONDA_DIR}/bin/conda clean -fy \
&& sudo chown -R :conda ${CONDA_DIR}/envs

# Add a file that is to be sourced from .bashrc and from the devops pipeline stages
Expand All @@ -90,7 +90,7 @@ RUN . /opt/conda/etc/profile.d/conda.sh \


# Install the yardl tool
ARG YARDL_VERSION=0.6.4
ARG YARDL_VERSION=0.6.6
RUN wget --quiet "https://github.com/microsoft/yardl/releases/download/v${YARDL_VERSION}/yardl_${YARDL_VERSION}_linux_x86_64.tar.gz" \
&& tar -xzf "yardl_${YARDL_VERSION}_linux_x86_64.tar.gz" \
&& mv yardl "/opt/conda/envs/${CONDA_ENVIRONMENT_NAME}/bin/" \
Expand Down
5 changes: 2 additions & 3 deletions .github/actions/configure-mrd-build-environment/action.yml
Original file line number Diff line number Diff line change
Expand Up @@ -17,13 +17,12 @@ runs:
tar -xzf "yardl_${YARDL_VERSION}_linux_x86_64.tar.gz"
pwd >> $GITHUB_PATH

- name: Setup Mambaforge
- name: Setup Miniforge
uses: conda-incubator/setup-miniconda@v3
with:
miniforge-version: latest
# Do not specify environment file - see Cache step below
activate-environment: mrd
use-mamba: true

- name: Get Date
id: get-date
Expand All @@ -43,7 +42,7 @@ runs:

- name: Update Environment
shell: bash -el {0}
run: mamba env update -n mrd -f environment.yml
run: conda env update -n mrd -f environment.yml
if: steps.cache-conda.outputs.cache-hit != 'true'

- name: Set up MATLAB
Expand Down
2 changes: 1 addition & 1 deletion .github/workflows/mrd_build.yml
Original file line number Diff line number Diff line change
Expand Up @@ -11,7 +11,7 @@ on:
env:
# Increase this to manually reset the conda environment cache
CONDA_CACHE_NUMBER: 0
YARDL_VERSION: 0.6.4
YARDL_VERSION: 0.6.6

defaults:
run:
Expand Down
2 changes: 1 addition & 1 deletion VERSION
Original file line number Diff line number Diff line change
@@ -1 +1 @@
2.1.2-preview
2.2.0
2 changes: 1 addition & 1 deletion cpp/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -28,7 +28,7 @@ endif ()
message(STATUS "Build type: ${CMAKE_BUILD_TYPE}")

if(NOT CMAKE_CXX_STANDARD)
set(CMAKE_CXX_STANDARD 17)
set(CMAKE_CXX_STANDARD 20)
endif()
message(STATUS "C++ Standard: ${CMAKE_CXX_STANDARD}")

Expand Down
2 changes: 1 addition & 1 deletion cpp/conda/build.sh
Original file line number Diff line number Diff line change
Expand Up @@ -9,7 +9,7 @@ echo 'Building mrd conda package using CMake...'

cmake -GNinja \
-DCMAKE_BUILD_TYPE=Release \
-DCMAKE_CXX_STANDARD=17 \
-DCMAKE_CXX_STANDARD=20 \
-DCMAKE_INSTALL_PREFIX=${PREFIX} \
-DCMAKE_OSX_SYSROOT=${CONDA_BUILD_SYSROOT} \
../
Expand Down
21 changes: 9 additions & 12 deletions cpp/conda/meta.yaml
Original file line number Diff line number Diff line change
Expand Up @@ -17,21 +17,18 @@ requirements:
- cmake
- ninja
host:
- fftw=3.3.9
- fmt=12.0.0
- fftw=3.3.*
- howardhinnant_date=3.0.3
- hdf5=1.10.6
- imagemagick=7.0.11_13
- ismrmrd>=1.14.2
- hdf5=1.14.*
- imagemagick>=7.1.2_8
- ismrmrd>=1.15.0
- nlohmann_json=3.12.0
- xtensor=0.25.0
- xtensor-fftw=0.2.5
- xtensor=0.27.1
run:
- fmt=12.0.0
- hdf5=1.10.6
- imagemagick=7.0.11_13
- ismrmrd>=1.14.2
- fftw=3.3.9
- fftw=3.3.*
- hdf5=1.14.*
- imagemagick=7.1.2_8
- ismrmrd>=1.15.0

about:
home: https://ismrmrd.github.io/mrd
Expand Down
13 changes: 7 additions & 6 deletions cpp/mrd-tools/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -10,7 +10,9 @@ add_executable(mrd_stream_to_hdf5 mrd_stream_to_hdf5.cc)
target_link_libraries(mrd_stream_to_hdf5 mrd_generated)

add_executable(mrd_stream_recon mrd_stream_recon.cc)
target_link_libraries(mrd_stream_recon fftw3f mrd_generated)
target_link_libraries(mrd_stream_recon mrd_generated fftw3f)
# Suppress false positive warnings from GCC 15+ with xtensor complex operations
target_compile_options(mrd_stream_recon PRIVATE "-Wno-array-bounds" "-Wno-stringop-overflow")

install(TARGETS
mrd_hdf5_to_stream
Expand Down Expand Up @@ -38,17 +40,16 @@ endif()


find_package(ImageMagick COMPONENTS Magick++)
find_package(fmt)

if (ImageMagick_FOUND AND fmt_FOUND)
message(STATUS "ImageMagick and fmt found.")
if (ImageMagick_FOUND)
message(STATUS "ImageMagick found.")
include_directories(${ImageMagick_INCLUDE_DIRS})

add_executable(mrd_image_stream_to_png mrd_image_stream_to_png.cc)
target_compile_options(mrd_image_stream_to_png PRIVATE "-DMAGICKCORE_QUANTUM_DEPTH=8" "-DMAGICKCORE_HDRI_ENABLE=0")
target_link_libraries(mrd_image_stream_to_png mrd_generated ${ImageMagick_LIBRARIES} fmt::fmt)
target_link_libraries(mrd_image_stream_to_png mrd_generated ${ImageMagick_LIBRARIES})

install(TARGETS mrd_image_stream_to_png DESTINATION bin)
else()
message(STATUS "ImageMagick or fmt not found. Skipping mrd_image_stream_to_png.")
message(STATUS "ImageMagick not found. Skipping mrd_image_stream_to_png.")
endif()
2 changes: 1 addition & 1 deletion cpp/mrd-tools/converters.cc
Original file line number Diff line number Diff line change
Expand Up @@ -4,7 +4,7 @@
#include <ismrmrd/version.h>

#include <date/date.h>
#include <xtensor/xview.hpp>
#include <xtensor/views/xview.hpp>

namespace mrd::converters {

Expand Down
114 changes: 114 additions & 0 deletions cpp/mrd-tools/fftw_wrappers.h
Original file line number Diff line number Diff line change
@@ -0,0 +1,114 @@
#pragma once

#include <complex>
#include <cstring>
#include <fftw3.h>
#include <xtensor.hpp>
#include <xtensor/misc/xmanipulation.hpp>

namespace fftw_wrappers {

namespace detail {

// Helper function for 1D FFTW operations
template <class E>
inline auto fftw_1d_impl(E&& e, int direction, bool normalize) -> xt::xarray<std::complex<float>> {
auto size = e.size();

std::vector<size_t> shape = {size};
xt::xarray<std::complex<float>> result = xt::zeros<std::complex<float>>(shape);
xt::xarray<std::complex<float>> input = e;

fftwf_complex* in = reinterpret_cast<fftwf_complex*>(input.data());
fftwf_complex* out = reinterpret_cast<fftwf_complex*>(result.data());

fftwf_plan plan = fftwf_plan_dft_1d(size, in, out, direction, FFTW_ESTIMATE);
fftwf_execute(plan);
fftwf_destroy_plan(plan);

if (normalize) {
result /= static_cast<float>(size);
}

return result;
}

// Helper function for 2D FFTW operations
template <class E>
inline auto fftw_2d_impl(E&& e, int direction, bool normalize) -> xt::xarray<std::complex<float>> {
auto shape = e.shape();

if (shape.size() != 2) {
throw std::runtime_error(direction == FFTW_FORWARD ? "fft2 requires 2D input" : "ifft2 requires 2D input");
}

auto n0 = shape[0];
auto n1 = shape[1];

std::vector<size_t> result_shape = {n0, n1};
xt::xarray<std::complex<float>> result = xt::zeros<std::complex<float>>(result_shape);
xt::xarray<std::complex<float>> input = e;

fftwf_complex* in = reinterpret_cast<fftwf_complex*>(input.data());
fftwf_complex* out = reinterpret_cast<fftwf_complex*>(result.data());

fftwf_plan plan = fftwf_plan_dft_2d(n0, n1, in, out, direction, FFTW_ESTIMATE);
fftwf_execute(plan);
fftwf_destroy_plan(plan);

if (normalize) {
result /= static_cast<float>(n0 * n1);
}

return result;
}
} // namespace detail

// FFT shift functions
template <class E>
inline auto fftshift(E&& e) -> xt::xarray<typename std::decay_t<E>::value_type> {
xt::xarray<typename std::decay_t<E>::value_type> result = e;

for (std::size_t ax = 0; ax < result.dimension(); ++ax) {
auto size = result.shape()[ax];
auto shift = size / 2;
result = xt::roll(result, shift, ax);
}
return result;
}

template <class E>
inline auto ifftshift(E&& e) -> xt::xarray<typename std::decay_t<E>::value_type> {
xt::xarray<typename std::decay_t<E>::value_type> result = e;

for (std::size_t ax = 0; ax < result.dimension(); ++ax) {
auto size = result.shape()[ax];
auto shift = (size + 1) / 2; // Different from fftshift for odd sizes
result = xt::roll(result, shift, ax);
}
return result;
}

// 1D FFT/IFFT wrappers - now using shared implementation
template <class E>
inline auto fft(E&& e) -> xt::xarray<std::complex<float>> {
return detail::fftw_1d_impl(std::forward<E>(e), FFTW_FORWARD, false);
}

template <class E>
inline auto ifft(E&& e) -> xt::xarray<std::complex<float>> {
return detail::fftw_1d_impl(std::forward<E>(e), FFTW_BACKWARD, true);
}

// 2D FFT/IFFT wrappers - now using shared implementation
template <class E>
inline auto fft2(E&& e) -> xt::xarray<std::complex<float>> {
return detail::fftw_2d_impl(std::forward<E>(e), FFTW_FORWARD, false);
}

template <class E>
inline auto ifft2(E&& e) -> xt::xarray<std::complex<float>> {
return detail::fftw_2d_impl(std::forward<E>(e), FFTW_BACKWARD, true);
}

} // namespace fftw_wrappers
2 changes: 0 additions & 2 deletions cpp/mrd-tools/ismrmrd_to_mrd.cc
Original file line number Diff line number Diff line change
Expand Up @@ -6,11 +6,9 @@
#include <ismrmrd/meta.h>
#include <ismrmrd/serialization_iostream.h>
#include <ismrmrd/version.h>
#include <xtensor/xview.hpp>

#include <date/date.h>


void print_usage(std::string program_name) {
std::cerr << "Usage: " << program_name << std::endl;
std::cerr << " -i|--input <input ISMRMRD stream> (default: stdin)" << std::endl;
Expand Down
13 changes: 3 additions & 10 deletions cpp/mrd-tools/mrd_image_stream_to_png.cc
Original file line number Diff line number Diff line change
@@ -1,14 +1,7 @@
// NOTE:
// GCC 13+ appears to emit false positive warnings for array-bounds
// when using xt::views (specifically the calls to `xt::amax` below).
// The statements are correct and there's no overflow, so we can suppress the warnings.
#pragma GCC diagnostic push
#pragma GCC diagnostic ignored "-Warray-bounds"
#include "mrd/binary/protocols.h"
#pragma GCC diagnostic pop

#include <Magick++.h>
#include <fmt/format.h>
#include <format>
#include <iostream>

void print_usage(std::string program_name) {
Expand Down Expand Up @@ -97,8 +90,8 @@ int main(int argc, char** argv) {
image.depth(8);
image.type(Magick::GrayscaleType);

// Use fmt to generate filename from prefix and increment with 6 digits (e.g. prefix_000001.png)
std::string filename = fmt::format("{}{:05d}.png", prefix, image_count++);
// Generate filename from prefix and increment with 5 digits (e.g. prefix_00001.png)
std::string filename = std::format("{}{:05d}.png", prefix, image_count++);
image.write(filename);
if (verbose) {
std::cerr << "Generated image " << filename << std::endl;
Expand Down
11 changes: 5 additions & 6 deletions cpp/mrd-tools/mrd_phantom.cc
Original file line number Diff line number Diff line change
@@ -1,13 +1,12 @@
#include "mrd/binary/protocols.h"
#include "shepp_logan_phantom.h"
#include "fftw_wrappers.h"

#include <iostream>
#include <random>
#include <xtensor-fftw/basic.hpp>
#include <xtensor-fftw/helper.hpp>
#include <xtensor/xio.hpp>
#include <xtensor/xrandom.hpp>
#include <xtensor/xview.hpp>
#include <xtensor/io/xio.hpp>
#include <xtensor/generators/xrandom.hpp>
#include <xtensor/views/xview.hpp>

using namespace mrd;

Expand Down Expand Up @@ -329,7 +328,7 @@ int main(int argc, char** argv) {
// Apply FFT
coil_images = fftshift(coil_images);
for (unsigned int c = 0; c < ncoils; c++) {
auto tmp1 = xt::fftw::fft2(xt::xarray<std::complex<float>>(xt::view(coil_images, c, 0, xt::all(), xt::all())));
auto tmp1 = fftw_wrappers::fft2(xt::xarray<std::complex<float>>(xt::view(coil_images, c, 0, xt::all(), xt::all())));
tmp1 /= std::sqrt(1.0f * tmp1.size());
xt::view(coil_images, c, 0, xt::all(), xt::all()) = tmp1;
}
Expand Down
29 changes: 9 additions & 20 deletions cpp/mrd-tools/mrd_stream_recon.cc
Original file line number Diff line number Diff line change
@@ -1,24 +1,9 @@

#include "fftw_wrappers.h"
#include "mrd/binary/protocols.h"
#include "mrd/protocols.h"
#include "mrd/types.h"

#include <xtensor-fftw/basic.hpp>
#include <xtensor-fftw/helper.hpp>

// NOTE:
// GCC 13+ appears to emit false positive warnings for array-bounds and stringop-overflow
// when using xt::view with xt::range (specifically in `remove_oversampling` below).
// The assignment is correct and there's no overflow, so we can suppress the warnings.
#pragma GCC diagnostic push
#pragma GCC diagnostic ignored "-Warray-bounds"
#pragma GCC diagnostic ignored "-Wstringop-overflow"
#include <xtensor/xview.hpp>
#pragma GCC diagnostic pop

xt::xtensor<std::complex<float>, 4> fftshift(xt::xtensor<std::complex<float>, 4> x) {
return xt::roll(xt::roll(x, x.shape(3) / 2, 3), x.shape(2) / 2, 2);
}
#include <xtensor/views/xview.hpp>

void remove_oversampling(mrd::Acquisition& acq, const mrd::EncodingType& enc) {
auto new_shape = acq.data.shape();
Expand All @@ -30,16 +15,20 @@ void remove_oversampling(mrd::Acquisition& acq, const mrd::EncodingType& enc) {
xt::xtensor<std::complex<float>, 2> new_data = xt::zeros<std::complex<float>>(new_shape);
for (size_t c = 0; c < acq.Coils(); c++) {
auto ft_line = xt::xarray<std::complex<float>>(xt::view(acq.data, c, xt::all()));
ft_line = xt::fftw::fftshift(xt::fftw::ifft(xt::fftw::ifftshift(ft_line)));
ft_line = fftw_wrappers::fftshift(fftw_wrappers::ifft(fftw_wrappers::ifftshift(ft_line)));
ft_line *= std::sqrt(1.0f * ft_line.size());
ft_line = xt::view(ft_line, xt::range(x_pad, rNx + x_pad));
ft_line = xt::fftw::fftshift(xt::fftw::fft(xt::fftw::ifftshift(ft_line)));
ft_line = fftw_wrappers::fftshift(fftw_wrappers::fft(fftw_wrappers::ifftshift(ft_line)));
ft_line /= std::sqrt(1.0f * ft_line.size());
xt::view(new_data, c, xt::all()) = ft_line;
}
acq.data = new_data;
}

xt::xtensor<std::complex<float>, 4> fftshift(xt::xtensor<std::complex<float>, 4> x) {
return xt::roll(xt::roll(x, x.shape(3) / 2, 3), x.shape(2) / 2, 2);
}

void print_usage(std::string program_name) {
std::cerr << "Usage: " << program_name << std::endl;
std::cerr << " -i|--input <input MRD stream> (default: stdin)" << std::endl;
Expand Down Expand Up @@ -190,7 +179,7 @@ int main(int argc, char** argv) {
auto slice = xt::view(buffer, c, s, xt::all(), xt::all(), xt::all(), xt::all());
slice = fftshift(slice);
for (unsigned int coil = 0; coil < slice.shape(0); coil++) {
auto tmp1 = xt::fftw::ifft2(xt::xarray<std::complex<float>>(xt::view(slice, coil, 0, xt::all(), xt::all())));
auto tmp1 = fftw_wrappers::ifft2(xt::xarray<std::complex<float>>(xt::view(slice, coil, 0, xt::all(), xt::all())));
tmp1 *= std::sqrt(1.0f * tmp1.size());
xt::view(slice, coil, 0, xt::all(), xt::all()) = tmp1;
}
Expand Down
Loading
Loading