OSCR

Group-Patch Joint Compression: Compressing Dynamic B<sub>0</sub> and Static RF Spatial Modulations Across k-Space Subregion Groups for Highly Accelerated MRI.

Code ↔ Paper

10 matches between paragraphs of the paper and lines of its authors' code, computed by the harvester (lexical-v1). Click a colored paragraph or line to see its counterpart.

The 10 matches
  1. [1] § Theory › Multidimensional Compressibility and Periodic Redundancy ↔ Recon/3D_CS_Wave_9T_RAM1T/Run_B0SENSE_recon.m, lines 107–140 · score 0.78 · variable density Poisson, phase encoded steps, receiver channels, disc, aliased, FFT
  2. [2] § Methods › Pulse Sequence and Data Acquisition ↔ Pulseq_sequence/writeWaveCAIPI.m, lines 47–93 · score 0.77 · dual ACS, retrospectively undersampled, Wave CAIPI, readout oversampling, strengths, phase encoding
  3. [3] § Methods › Two Data Compression Approaches ↔ Recon/2D_CS_LocalB0Coils_9T/Run_B0SENSE_recon.m, lines 44–70 · score 0.77 · local B0 coils, oversampled readout, modulation period, receiver channels, rapid B0 modulation, phase encoded
  4. [4] § Methods › Pulse Sequence and Data Acquisition ↔ Pulseq_sequence/writeWaveCAIPI.m, lines 47–93 · score 0.76 · slice oversampling, wave gradients, slice thickness, Wave CAIPI, 230, TE
  5. [5] § Methods › Pulse Sequence and Data Acquisition ↔ Recon/2D_CS_LocalB0Coils_9T/Run_B0SENSE_recon.m, lines 44–70 · score 0.76 · local B0 coil, B0 field modulations, readout oversampling, rapid B0 modulations, strengths, FLASH
  6. [6] § Methods › Two Data Compression Approaches ↔ Recon/3D_L2_WaveCAIPI_3T/Run_B0SENSE_recon.m, lines 206–256 · score 0.73 · RF sensitivity maps, composite encoding matrix, Wave CAIPI, modulation period, compression matrices, B0 modulation
  7. [7] § Methods › Two Image Reconstruction Algorithms ↔ Recon/3D_CS_Wave_9T_RAM1T/Run_B0SENSE_recon.m, lines 50–80 · score 0.66 · wavelet regularization, undersampling factor, compressed sensing, global, B0 modulation, thresholding
  8. [8] § Results › Visualization of Compressed Spatial Encoding Functions ↔ Recon/3D_L2_WaveCAIPI_3T/Run_B0SENSE_recon.m, lines 206–256 · score 0.66 · RF sensitivity maps, composite encoding, B0 modulation period, B0 encoding, joint compression, matrix
  9. [9] § Theory › From RF Array to Joint B0 ‐RF Compression ↔ Recon/2D_CS_LocalB0Coils_9T/Run_B0SENSE_recon.m, lines 175–222 · score 0.61 · composite spatial encoding, receiver channel, compression matrices, virtual, segmented, joint compression
  10. [10] § Results › RF‐Only Versus B0 ‐RF Compression for In Vivo 3D Scans at 9.4 T ↔ Recon/3D_L2_WaveCAIPI_3T/Run_B0SENSE_recon.m, lines 109–141 · score 0.61 · L2 Wave CAIPI, rapid B0 modulations, recon, compression factors, FLASH, joint compression

Paper

Loaded from Europe PMC by your browser, not stored by OSCR: doi.org · Europe PMC

The paper is loaded when this pane is shown.

The authors' code

MATLAB · 529 lines · 19 KB · BSD-3-Clause · 3 matches

  1. %% Descriptions
  3. % This script is a basic demo of the joint compression technique
  4. % for 2D GRE data (2025-06-18, in-vivo, 9.4T human scanner with 8-channel local B0 coils)
  5. % with rapid B0 field modulation (readout oversampling),
  6. % and 32 RF receivers,
  7. % using compressed sensing reconstruction (FISTA).
  9. % For comparison, data with 4x RF array compression is also provided.
  10. % Please manually change the loaded data name to compare the speed.
  12. % In this implementation, the 2D readout encoding matrix is explicitly
  13. % written. This can be optimized by e.g., only writing readout encoding matrix
  14. % across a B0 modulation period, and reusing it for different periods, to
  15. % save memory. Nevertheless, this issue is not crucial for 2D.
  17. % Please cite the paper:
  18. % R.Tian, K.Scheffler, "Group-Patch Joint Compression: Compressing Dynamic B0 and Static RF Spatial Modulations Across k-Space Subregion Groups for Highly Accelerated MRI," Magnetic Resonance in Medicine, (2026), https://doi.org/10.1002/mrm.70439.
  19. % See Supporting Information Figure S3 for this dataset.
  21. % For questions, please contact Rui Tian
  22. % ([email hidden]/[email hidden])
  24. %% Environment
  26. % Add ToolBox to the path, for FFT, Wavelet((c) Michael Lustig 2007), and TV((c) Michael Lustig 2007) functions
  27. addpath(genpath('./'))
  29. %% Load data and parameters
  31. load('data_param_2D_JCS.mat')
  32. % load('data_param_2D_RFcompression.mat')
  34. fprintf('\n');
  35. disp('Data Loaded.')
  37. % Fully-sampled k-space data: dat_exp[oversampled readout,phase-encoding,1,1,RF receivers]
  38. % Auto-calibrated B0 modulation maps: maps_B0Evo[oversampled readout,read pixel,phase pixel]
  39. % maps_SENSE[read pixel,phase pixel,1,1,RF receivers]
  40. % idx_z2D: slice number for multi-slice experiments
  41. % img_fg: mask
  42. % para: parameters
  44. %% Specify some local parameters
  46. % In this script, you may not need to tune these below
  47. nRO = para.nty_unOS; % readout k-space point number without oversampling
  48. nx = para.ntx; % phase encoding k-space point number
  49. nC = para.nC; % RF receiver channel number
  50. nOS = para.nyOS_flash; % readout oversampling factor
  51. idx_us_x = para.idx_usx; % undersampling factor along phase encoding
  52. nRecoPix_yx = para.nRecoPix_yx; % in-plane pixel number in recon grid
  53. nNgCh = para.nNgCh; % local B0 coil channel
  54. nMR = nRO*nOS; % readout k-space point number with oversampling
  55. nGridMR = nRecoPix_yx*nOS; % For recon grid, readout k-space point with oversampling
  56. bSENSE = para.bSENSE; % 1: using multiple RF receivers for joint sampling acceleration with SENSE
  57. bNgAcc = para.bNgAcc; % 1: rapid B0 field modulation ON
  58. bJointComp = para.bCompRF==2; % 1: Joint compression ON, 0: Joint compression OFF
  59. bNLReco = para.bNLReco; % 1: Compressed-sensing with L1 regularization is ON, as this script
  60. bWaThres = para.bWaThres; % 1: global wavelet threshold, 2: level-dependent wavelet threshold
  61. para.kperiod; % This is the sinusoidal period (45 oversampled readout points) of rapid B0 modulation for this dataset.
  63. % You might tune the following to change regularization strengths
  64. lambda_TV = para.lambda_TV; % coefficient for TV regularization term. higher->stronger
  65. lambda_walet = para.lambda_walet; % coefficient for wavelet regularization term. higher->stronger
  66. if bJointComp==1 % if joint compression is ON
  67. % This dataset has B0 modulation with period of 45 readout k-space points (including 8x oversampling)
  68. ntSeg_cs = para.ntSeg_cs; % the readout k-space point number for one compression group, in a B0 modulation period.
  69. nJointCs = para.nJointCs; % the joint compression factor
  70. end
  72. %% RF array compression
  74. if para.bCompRF ==1
  76. disp('RF-only compression starts......')
  78. maps_SENSE_cal = squeeze(permute(maps_SENSE,[1,2,3,5,4]));%[read,phase,(3Dslice),RF,2Dslice]
  80. tstRFCs = tic;
  81. [matA] = ComputeMat_CSRF(maps_SENSE_cal,img_fg,para.nCompRF);
  82. para.tRFcs_transformation = toc(tstRFCs);
  83. fprintf('\n');
  84. disp(['Time for calculating RF compression matrix: ',num2str(para.tRFcs_transformation),' s'])
  86. dat_exp = permute(dat_exp,[1,2,5,3,4]);
  87. maps_SENSE = permute(maps_SENSE,[1,2,5,3,4]);
  89. tstRFCs = tic;
  90. [dat_exp,maps_SENSE] = CompressRF(dat_exp,maps_SENSE,matA);
  92. para.tRFcs_apply = toc(tstRFCs);
  93. fprintf('\n');
  94. disp(['Time for applying RF compression matrix: ',num2str(para.tRFcs_apply),' s'])
  95. fprintf('\n')
  97. dat_exp = permute(dat_exp,[1,2,4,5,3]);
  98. maps_SENSE = permute(maps_SENSE,[1,2,4,5,3]);
  100. end
  102. %% Swap input data dimensions for further processing
  104. % dat_exp[ny_OS,nx,nz3D,nz2D->1,nC] -> [ny_OS,nPEx_intlp,nPEz_intlp,nC,nz2D->1]
  105. dat_exp = permute(dat_exp,[1,2,3,5,4]);
  107. %% Retrospective undersampling mask
  109. % Initialize the mask (interpolated size)
  110. mask_kxz = zeros(nRecoPix_yx,1);
  112. disp(['The undersampling factor: ', num2str(idx_us_x)])
  114. % Define the sampling point in the mask
  115. mask_kxz(1:idx_us_x:end) = 1;
  116. mask_kxz = single(mask_kxz);
  117. img_fg = single(img_fg);
  119. % mask index in image and k-space
  120. [ind_kx,ind_kz] = find(mask_kxz==1);
  121. ind_k = find(mask_kxz==1);
  122. [ind_y,ind_x,ind_z] = ind2sub(size(img_fg),find(img_fg));
  123. ind_fg = find(img_fg==1);
  124. sz_kmask = size(mask_kxz);
  125. sz_fg = size(img_fg);
  127. % Apply undersampling mask to fully-sampled k-space data
  128. dat_exp = dat_exp.*reshape(mask_kxz,1,nRecoPix_yx,1);
  130. % Keep the undersampling data in image domain for optional check.
  131. img_exp_us = SymFft(dat_exp,[1 2 3]);
  132. img_exp_us_sos = sqrt(sum(img_exp_us.*conj(img_exp_us),4));
  134. % Reduce memory by discarding omitted phase-encoded steps
  135. dat_exp = reshape(dat_exp,[size(dat_exp,1),size(dat_exp,2)*size(dat_exp,3),size(dat_exp,4)]);
  136. dat_exp = dat_exp(:,ind_k,:);
  138. % Reorganize the data dimension for multiple receivers
  139. if bSENSE == 1
  140. dat_exp = permute(dat_exp,[1 3 2]);
  141. end
  143. %% Prepare readout encoding matrix
  145. ERead = single(1);
  147. % If rapid B0 modulation ON, take the auto-calibrated B0 modulation matrix
  148. if bNgAcc==1
  149. ERead = ERead.*maps_B0Evo; clear maps_B0Evo;
  150. end
  152. % Encoding matrix for
  153. % Cartesian sampling by linear gradient
  154. % - > oversampled frequency-encoding
  155. EFe = dftmtx(nGridMR);
  156. EFe = circshift(conj(EFe),[floor(nGridMR/2) floor(nGridMR/2)]);
  157. cen_PE = floor(nGridMR/2)+1;
  158. min_PE = cen_PE - floor(nRecoPix_yx/2);
  159. max_PE = cen_PE + ceil(nRecoPix_yx/2)-1;
  160. EFe = EFe(:,min_PE:max_PE); % [nty_RO,nRecoPix_yx]
  161. ERead = ERead.*EFe;
  162. ERead = reshape(ERead,nMR,[]);
  163. ERead = single(ERead(:,ind_fg));
  165. % Encoding matrix for RF receive sensitivity encoding
  166. sense_maps_reshape=0;
  167. para.recnC = size(maps_SENSE,5);
  168. if bJointComp == 0 % if no joint compression, readout amtrix combined with sense matrix
  169. sense_maps_reshape = reshape(maps_SENSE,[1,nRecoPix_yx*nRecoPix_yx,para.recnC]);
  170. ERead = ERead.*sense_maps_reshape(:,ind_fg,:); clear sense_maps_reshape;
  171. ERead = permute(ERead,[1 3 2]);
  172. ERead = reshape(ERead,[],length(ind_fg));
  173. end
  175. %% Joint compression to reduce the reeadout model size
  177. if bJointComp == 1 % Joint compression ON
  179. ntnC_cs = ntSeg_cs*nC/nJointCs; % the number of virtual k-space points compressed from (a reaodut segment + all receiver channels)
  181. if mod(nGridMR,ntSeg_cs)~=0 % Here, the total readout point number = integer x readout point number in a segment, for convenience processing
  182. error('Compression in readout time domain is not divided into interger number of regions')
  183. end
  185. % Reorganize sensitivity matrix
  186. maps_SENSE = reshape(permute(maps_SENSE,[5,1,2,3,4]),[1,nC,nRecoPix_yx*nRecoPix_yx]);
  187. maps_SENSE = single(maps_SENSE(:,:,ind_fg));
  189. % Initiate compression matrix
  190. matA = single(zeros(ntnC_cs,ntSeg_cs*nC,para.kperiod/ntSeg_cs));
  192. disp('Joint compression starts......')
  195. % Calculate compression matrices from composite spatial encoding
  196. % maps.
  197. tstJointCs = tic;
  198. for idt_seg = 1:(para.kperiod/ntSeg_cs) % loop over different compression groups in a B0 modulation period
  200. % index for a readout segment, for a compression group
  201. t_st = 1+(idt_seg-1)*ntSeg_cs;
  202. t_end = t_st +(ntSeg_cs-1);
  204. % select and reshape parts of ERead, for a compression group
  205. E_seg = reshape(ERead(t_st:t_end,:),ntSeg_cs,1,[]);
  207. % construct the composite encoding matrix, for a compression
  208. % group
  209. E_joint_pix = bsxfun(@times,maps_SENSE,E_seg);
  210. E_joint_pix = reshape(E_joint_pix,ntSeg_cs*nC,length(ind_fg));
  211. E_joint_H_pix = E_joint_pix';
  213. % scaling coefficient
  214. coeff_scale = reshape(sum(conj(E_joint_pix) .* E_joint_pix, 1),1,length(ind_fg)); % [1 x p]
  216. matP = (E_joint_pix./coeff_scale)*E_joint_H_pix;
  218. [U,S,V] = svd(matP,'econ');
  220. VH = V';
  222. matA(:,:,idt_seg)= VH(1:ntnC_cs,:);
  224. % There there a few residue data points not belonging to
  225. % sets of every ntSeg_cs readout points
  226. if mod(nGridMR,ntSeg_cs)~=0
  228. % index for readout segment
  229. t_end = nGridMR;
  230. t_st = t_end - (ntSeg_cs+mod(nGridMR,ntSeg_cs)) +1;
  232. % select and reshape parts of ERead
  233. E_seg = reshape(ERead(t_st:t_end,:),ntSeg_cs,1,[]);
  235. % initialize the matrix to be decomposed
  236. matP = zeros(nC*ntSeg_cs,nC*ntSeg_cs);
  238. % update the k-space compression matrix over all pixels
  239. E_joint_pix = bsxfun(@times,maps_SENSE,E_seg);
  241. E_joint_pix = reshape(E_joint_pix,ntSeg_cs*nC,length(ind_fg));
  242. E_joint_H_pix = E_joint_pix';
  244. % scaling vector, pixel dependent
  245. coeff_scale = reshape(sum(conj(E_joint_pix) .* E_joint_pix, 1),1,length(ind_fg)); % [1 x p]
  247. % square matrix P for a compression group
  248. matP = (E_joint_pix./coeff_scale)*E_joint_H_pix;
  250. % svd decompositing the square matrix P. the
  251. % compressibility can be checked vis S
  252. [U,S,V] = svd(matP,'econ');
  254. VH = V';
  256. % select the subsets of singular vectors
  257. ntnC_cs_last = round(length(t_st:t_end)*nC/nJointCs);
  258. matA_res= VH(1:ntnC_cs_last,:);
  259. end
  260. end
  262. fprintf('\n')
  263. X = (['The (joint) compression matrices are calculated, for compression factor: ',num2str(nC*ntSeg_cs/ntnC_cs)]);
  264. disp(X)
  265. fprintf('\n')
  267. para.tcs_transformation = toc(tstJointCs);
  268. disp(['Time for calculating compression matrices: ',num2str(para.tcs_transformation),' s'])
  269. fprintf('\n')
  272. % Apply the compression matrices
  273. tstJointCs_apply = tic;
  275. % If total readout points before compression ~= integer x readout
  276. % segment points in a compression group
  277. if mod(nGridMR,ntSeg_cs)~=0 % if remained readout points do not form a compression group
  278. tres_ed = nGridMR;
  279. tres_st = tres_ed - (ntSeg_cs+mod(nGridMR,ntSeg_cs)) +1;
  280. ERead_res = ERead(tres_st:tres_ed,:);
  281. dat_exp_res = dat_exp(tres_st:tres_ed,:,:);
  283. tsegs_st = 1;
  284. tsegs_ed = tres_st-1;
  285. else % simpler case
  286. tsegs_st = 1;
  287. tsegs_ed = nGridMR;
  288. end
  290. if round(length(tsegs_st:tsegs_ed)/para.kperiod)~=(length(tsegs_st:tsegs_ed)/para.kperiod)
  291. if mod(para.kperiod-mod(length(tsegs_st:tsegs_ed),para.kperiod),ntSeg_cs)~=0
  292. error('Error! Requires reprogramming of compressed groups. The remained readout pts cannot be compressed by given compression matrices.')
  293. end
  294. dat_exp = padarray(dat_exp,[para.kperiod-mod(length(tsegs_st:tsegs_ed),para.kperiod),0,0],0,'post');
  295. ERead = padarray(ERead,[para.kperiod-mod(length(tsegs_st:tsegs_ed),para.kperiod),0],0,'post');
  296. tsegs_ed = tsegs_ed+para.kperiod-mod(length(tsegs_st:tsegs_ed),para.kperiod);
  297. end
  299. % Reorganize data to be compressed
  300. dat_exp = dat_exp(tsegs_st:tsegs_ed,:,:);
  301. dat_exp = reshape(dat_exp,...
  302. ntSeg_cs,...
  303. para.kperiod/ntSeg_cs,...
  304. length(tsegs_st:tsegs_ed)/para.kperiod,...
  305. nC,...
  306. length(ind_k));
  307. dat_exp = permute(dat_exp,[1,4,2,3,5]);
  308. dat_exp = reshape(dat_exp,...
  309. ntSeg_cs*nC,...
  310. 1,...
  311. para.kperiod/ntSeg_cs,...
  312. length(tsegs_st:tsegs_ed)/para.kperiod,...
  313. length(ind_k)); % [nToBeComps,1,sets of comp. matrix,different mod. periods,phase encoding]
  315. % data is compressed
  316. dat_exp_cs = pagemtimes(matA,dat_exp);
  318. % Reorganize readout encoding matrix
  319. ERead = reshape(ERead,ntSeg_cs,para.kperiod/ntSeg_cs,length(tsegs_st:tsegs_ed)/para.kperiod,length(ind_fg));
  320. ERead = permute(ERead,[1,3,4,2]); % ->[ntSeg_cs,length(tsegs_st:tsegs_ed)/para.kperiod,length(ind_fg),para.kperiod/ntSeg_cs]
  321. ERead_cs = single(zeros(ntnC_cs,1,length(tsegs_st:tsegs_ed)/para.kperiod,length(ind_fg),para.kperiod/ntSeg_cs)); %[compressed#,comp.seg.,yxz_mask]
  323. % compress readout encoding matrix
  324. for idx_mod = 1:(para.kperiod/ntSeg_cs)
  326. E_seg = ERead(:,:,:,idx_mod);
  327. E_seg = reshape(E_seg,...
  328. ntSeg_cs,...
  329. 1,...
  330. length(tsegs_st:tsegs_ed)/para.kperiod,...
  331. length(ind_fg));
  332. E_seg = E_seg.*reshape(maps_SENSE,[1,nC,1,length(ind_fg)]);
  333. E_seg = reshape(E_seg,ntSeg_cs*nC,...
  334. 1,...
  335. length(tsegs_st:tsegs_ed)/para.kperiod,...
  336. length(ind_fg)); % [nToBeComps,1,sets of comp. matrix,different mod. periods,phase encoding]
  338. ERead_cs(:,:,:,:,idx_mod) = pagemtimes(matA(:,:,idx_mod),E_seg);
  339. end
  340. ERead_cs = permute(ERead_cs,[1,5,3,4,2]);
  341. ERead_cs = reshape(ERead_cs,ntnC_cs/ntSeg_cs*length(tsegs_st:tsegs_ed),length(ind_fg));
  343. % manage if total readout points before compression ~= integer x readout
  344. % segment points in a compression group
  345. if nGridMR~=length(tsegs_st:tsegs_ed)
  346. nRes_temp = (length(tsegs_st:tsegs_ed)-nGridMR)*nC/nJointCs;
  347. ERead_cs = ERead_cs(1:(end-nRes_temp),:);
  348. dat_exp_cs = reshape(dat_exp_cs,[],size(dat_exp_cs,5));
  349. dat_exp_cs = dat_exp_cs(1:(end-nRes_temp),:);
  350. end
  352. para.tcs_apply = toc(tstJointCs_apply);
  353. disp(['Time for applying compression matrices: ',num2str(para.tcs_apply),' s'])
  354. fprintf('\n')
  356. X = (['The (joint) compression matrices are applied, for compression factor: ',num2str(nC*ntSeg_cs/ntnC_cs)]);
  357. disp(X)
  358. fprintf('\n')
  360. clear ERead ; % save memory
  361. clear dat_exp;
  362. clear E_seg;
  364. nMR = size(ERead_cs,1);nC=1;% para.recnC=1;
  366. else
  367. % keep using ..._cs to be compatible with the joint compressed
  368. % variables
  369. ERead_cs = ERead; clear ERead;
  370. dat_exp_cs = dat_exp; clear dat_exp;
  371. end
  373. %% Reorganize data dimensions for compressed sensing recon. The implementation may be further optimized. For 2D CS, this is fine so far.
  375. disp('Reorganizes data dimensions for CS recon......')
  377. % Map compressed encoding matrix with only sampled phase encoded steps back
  378. % to the matrix that include the zero-zilled phase encoded steps
  379. ERead_full = single(zeros(size(ERead_cs,1),size(img_fg,1)*size(img_fg,2)*size(img_fg,3)));
  380. ERead_full(:,ind_fg) = ERead_cs;
  381. ERead_full = reshape(ERead_full,size(ERead_cs,1),size(img_fg,1),size(img_fg,2),size(img_fg,3));
  382. ERead_cs = ERead_full;clear ERead_full;
  384. fprintf('\n');
  385. mem_GB = size(ERead_cs,1)*size(img_fg,1)*size(img_fg,2)*size(img_fg,3)*8*1e-9;
  386. X = ['Composite readout matrix used for recon: ',num2str(mem_GB),'GB'];
  387. disp(X)
  389. % Prepare the adjoint encoding matrix for CS recon, to avoid regenerating it during iteration.
  390. ERead_cs_H = reshape(ERead_cs,size(ERead_cs,1),[]);
  391. ERead_cs_H = ERead_cs_H';
  392. ERead_cs_H = reshape(ERead_cs_H,size(img_fg,1),size(img_fg,2),size(img_fg,3),size(ERead_cs_H,2));
  393. ERead_cs_H = permute(ERead_cs_H,[1,4,2,3]);
  395. %% FISTA recon of joint encoding of rapid B0 modulation & SENSE
  397. % reshape data into 1D
  398. dat_exp_cs = dat_exp_cs(:);
  400. % Extract the max. eigenvalue for normalization. You may need to calculate
  401. % or empirically get one value for your dataset.
  402. norm_A = para.sReco;
  404. % normalization of readout encoding matrix
  405. L=1.1;
  406. ERead_cs = ERead_cs./norm_A;
  407. ERead_cs_H = ERead_cs_H./norm_A;
  409. % initialize variable
  410. u0 = single(zeros(length(ind_y),1));
  411. z0 = single(zeros(length(ind_y),1));
  412. t0 = single(1);
  414. % iteration step 40. empirical. you can tune it.
  415. nIter = 40;
  417. % Recon starts
  418. fprintf('\n');
  419. disp('Iterative FISTA reconstruction starts......')
  420. if bNLReco==1 % nonlinear regularization is applied
  421. ct_tv = 0;
  422. t_nl = tic;
  423. for idx_tv = lambda_TV
  424. ct_tv = ct_tv +1;
  425. ct_wave = 0;
  426. for idx_wave = lambda_walet % 2e-4 for level dependent
  427. ct_wave = ct_wave+1;
  428. for idx_Iter=1:nIter
  430. % FISTA iteration
  432. if idx_Iter == 1 % Take initialized values
  433. u_this = u0;
  434. z_this = z0;
  435. t_this = t0;
  436. else % Take updated values from the last iteration
  437. u_this = u_next;
  438. z_this = z_next;
  439. t_this = t_next;
  440. end
  442. img0 = Afun_NGSENSE_fft(dat_exp_cs,'transp',ERead_cs_H,sz_kmask,ind_k,ind_kx,ind_kz,sz_fg,ind_fg,ind_y,ind_x,ind_z);
  443. data_est = Afun_NGSENSE_fft(u_this,'notransp',ERead_cs,sz_kmask,ind_k,ind_kx,ind_kz,sz_fg,ind_fg,ind_y,ind_x,ind_z);
  444. img_est = Afun_NGSENSE_fft(data_est,'transp',ERead_cs_H,sz_kmask,ind_k,ind_kx,ind_kz,sz_fg,ind_fg,ind_y,ind_x,ind_z);
  446. [deri_TV,cost_TV] = apply_TV3D(u_this,sz_fg,ind_fg,idx_tv);
  448. grad = (1/L).*(img0-img_est);
  449. z_next = u_this+grad;
  450. z_next = z_next-(1/L).*deri_TV(img_fg==1);
  452. z_next = apply_Wavelet(z_next,size(img_fg),img_fg,bWaThres,idx_wave);
  453. z_next = z_next(img_fg==1);
  455. % Update the inversion time
  456. t_next = 0.5*(1+sqrt(1+4*t_this^2));
  458. % Update the input to the next FISTA iteration
  459. u_next = z_next + (t_this-1)/t_next*(z_next-z_this);
  461. end
  462. end
  463. end
  465. tdur_nl = toc(t_nl);
  466. fprintf('\n');
  467. disp(['Nonlinear recon takes total time: ',num2str(tdur_nl), 's']);
  469. recon_temp = single(zeros([nRecoPix_yx,nRecoPix_yx]));
  470. recon_temp(ind_fg) = u_next;
  471. recon = recon_temp;
  473. else % no nonlinear regularization is applied
  475. t_nl = tic;
  476. for idx_Iter=1:nIter
  478. % FISTA iteration
  479. if idx_Iter == 1 % Take initialized values
  480. u_this = u0;
  481. z_this = z0;
  482. t_this = t0;
  483. else % Take updated values from the last iteration
  484. u_this = u_next;
  485. z_this = z_next;
  486. t_this = t_next;
  487. end
  489. img0 = Afun_NGSENSE_fft(dat_exp_cs,'transp',ERead_cs_H,sz_kmask,ind_k,ind_kx,ind_kz,sz_fg,ind_fg,ind_y,ind_x,ind_z);
  490. data_est = Afun_NGSENSE_fft(u_this,'notransp',ERead_cs,sz_kmask,ind_k,ind_kx,ind_kz,sz_fg,ind_fg,ind_y,ind_x,ind_z);
  491. img_est = Afun_NGSENSE_fft(data_est,'transp',ERead_cs_H,sz_kmask,ind_k,ind_kx,ind_kz,sz_fg,ind_fg,ind_y,ind_x,ind_z);
  493. grad = (1/L).*(img0-img_est);
  494. z_next = u_this+grad;
  496. % Update the inversion time
  497. t_next = 0.5*(1+sqrt(1+4*t_this^2));
  499. %Update the input to the next FISTA iteration
  500. u_next = z_next + (t_this-1)/t_next*(z_next-z_this);
  502. end
  504. tdur_nl = toc(t_nl);
  505. fprintf('\n');
  506. disp(['L2 recon takes total time: ',num2str(tdur_nl), 's']);
  508. recon_temp = single(zeros([nRecoPix_yx,nRecoPix_yx]));
  509. recon_temp(ind_fg) = u_next;
  510. recon = recon_temp;
  511. end
  514. % Display the recon.
  515. figure
  516. imagesc(abs(rot90(recon,1)))
  517. colormap(gray)
  518. if bJointComp==1
  519. title(['2D CS recon, joint compress ',num2str(round(para.tcs_transformation+para.tcs_apply,1)),' s, recon ',num2str(round(tdur_nl,1)), ' s'])
  520. else
  521. title(['2D CS recon, RF compress ',num2str(round(para.tRFcs_transformation+para.tRFcs_apply,1)),' s, recon ',num2str(round(tdur_nl,1)), ' s'])
  522. end
  524. fprintf('\n');
  525. disp('This script only provides you a basic demo about the benefits of the joint compression technique.')
  526. disp('The implementation can be further optimized.')

Run_B0SENSE_recon.m at commit 9a7ecb0, under BSD-3-Clause · at the source

Overview

Authors: Rui Tian1, Klaus Scheffler1,2
ORCID iDs: Rui Tian, Klaus Scheffler
  1. High‐Field MR Center, Max Planck Institute for Biological Cybernetics Tübingen Germany
  2. Department for Biomedical Magnetic Resonance University of Tübingen Tübingen Germany
Institutions: Max Planck Institute for Biological Cybernetics (Germany); University of Tübingen (Germany)
Journal: Magnetic resonance in medicine, volume 96, issue 5, pages 2106-2126
Dates: received 13 October 2025; accepted 6 May 2026; published online 23 July 2026; in print November 2026
Type: Research article · Language: English
License: CC BY
Identifiers: DOI 10.1002/mrm.70439 · PMID 42490765 · PMCID PMC13527271 · OpenAlex W7170148670
Open access: hybrid, a free copy (OpenAlex)
Preprint: arxiv:2601.05418
Status: code verified
Categories: structural MRI / diffusion (modality), human (organism), methods / tools (subfield)
Methods: Spectral & time-frequency, fMRI & imaging
Keywords: compressed sensing reconstruction, compression, FRONSAC, local B0 coil array, parallel imaging, Wave‐CAIPI
MeSH: Brain*, Data Compression*, Image Processing, Computer-Assisted*, Magnetic Resonance Imaging*, Algorithms, Compression Algorithms, Fourier Analysis, Humans, Radio Waves (* major topic)
Journal subjects: Imaging Methodology
Topic: Advanced MRI Techniques and Applications (Radiology, Nuclear Medicine and Imaging, Medicine), according to OpenAlex
Funding: European Research Council (834940)
Citations: cited by 1 paper (Europe PMC); 92 references in the paper

Abstract

Purpose: To accelerate MRI further, rapid B0 field modulations can be applied during oversampled readout to capture additional physical information, as in Wave‐CAIPI/FRONSAC/local B0 coils modulation techniques. These methods, however, turn the Fourier readout into a non‐Fourier‐encoded dimension that cannot be reconstructed by FFT, posing significant reconstruction challenges especially in compressed‐sensing or neural‐network frameworks.

Theory and Methods: Because the rapid B0 modulations still vary slowly relative to the oversampled ADC dwell time, we exploit this encoding redundancy by compressing k‐space patch‐by‐patch across subregions, each of which is jointly encoded by a distinct subset of B0 and RF (receive) spatial encoding functions. For each subset, a compression matrix is computed once and reused to compress all patches encoded by the same B0‐RF spatial modulations. This can be implemented by feeding subsets of B0 and RF spatial encoding maps into an adapted conventional RF array compression algorithm, mimicking an expanded set of virtual receiver channels. This approach was evaluated on human brain scans at 9.4 T/3 T.

Results: The proposed group‐patch joint compression achieves substantially higher compression factors than conventional RF‐only compression, while minimally compromising encoding efficiency. Typically, joint compression factors of 11×–20× led to negligible encoding loss, dramatically reducing reconstruction time and peak memory usage. For example, compressed‐sensing reconstruction took 1.4–5.1 s/2D slice, 177 s–10.1 min/3D volume on a high‐memory CPU node.

Conclusion: Given joint encoding of dynamic B0 and static RF fields, compressing multidimensional k‐space patches in separate groups outperforms compressing RF receivers alone. This substantially mitigates a fundamental computational bottleneck when combining rapid B0 and RF‐receiver modulations.

Reproduced under the paper's license (CC BY), from the paper cited above.

Repositories

Its files are read in the Code ↔ Paper reader above, with 10 matches between paragraphs and lines of code.

Zenodo 20533698

License: BSD-3-Clause
State: the link answers, verified on 27 September 2026
Evidence: files inventoried
Size: 1 file
Software Heritage: not checked
Found in: “Data Availability Statement”
Not found: README, license file, CITATION.cff, environment file, tests, continuous integration, documentation
Availability: 1 check, the latest on 27 September 2026: the link answers (HTTP 200)
  • 27 September 2026: the link answers (HTTP 200)
At the source: zenodo.org/records/20533698

ruitianwater/GroupPatchJCS_v1.0

License: BSD-3-Clause
State: the link answers, verified on 27 September 2026
Evidence: files inventoried
Commit: 9a7ecb08031b5f533870af47aee6527a51fedb18, 5 June 2026
Languages: MATLAB (209), Python (2), C++ (2), C/C++ (1)
Size: 274 files, 214 scripts
Software Heritage: not archived
Found in: “Data Availability Statement”
Holds: README, license file, tests, documentation
Not found: CITATION.cff, environment file, continuous integration
Tools: Wavelet Toolbox (2 files), Image Processing Toolbox (1 file)
Availability: 1 check, the latest on 27 September 2026: the link answers
  • 27 September 2026: the link answers
216 files
At the source: github.com/ruitianwater/GroupPatchJCS_v1.0

harmonizedmri.github.io/projects

License: none: the authors keep all their rights
State: the link answers, verified on 27 September 2026
Evidence: the link answers
Software Heritage: not checked
Found in: “Data Availability Statement”
Not found: README, license file, CITATION.cff, environment file, tests, continuous integration, documentation
Availability: 1 check, the latest on 27 September 2026: the link answers (HTTP 200)
  • 27 September 2026: the link answers (HTTP 200)
At the source: harmonizedmri.github.io/projects/

The paper's code and data availability statement is in the Data section.

Tracing map

Proposed by the machine: these links were found in the paper and verified at the source, without human review. The map will receive a Zenodo DOI once one of the paper's authors has validated it with their ORCID.

What the map holds:

  • 3 repositories of the authors' code, each at its verified commit, with its license and how the link was found in the paper;
  • 214 scripts, each with its path and the digest of its content;
  • 10 matches between paragraphs of the paper and lines of the code (method lexical-v1);
  • neither the text of the paper nor the code itself.

Its JSON (tracing-map.json) is deposited on Zenodo with its DOI once the map is validated.

Data

No dataset and no data link were found in the paper.

Data Availability Statement

The example code or data that support the findings of this study is openly available in Zenodo at https://zenodo.org/records/20533698, in Github at https://github.com/ruitianwater/GroupPatchJCS_v1.0, and collected at the HarmonizedMRI website at https://harmonizedmri.github.io/projects/.

Reproduced under the paper's license (CC BY), from the paper cited above.

Versions

The history of this record: each version stored by the harvester or made by a correction of its authors or of the maintainers of its code, and what changed in its facts. The texts of the paper (its abstract, its availability statements) are not part of it; versions that changed only those are not listed.

Version 2, 28 September 2026

  • Publisher: n/a → Wiley

Version 1, 27 September 2026: the first record

Recorded: type, language, journal, volume, issue, pages, dates, 2 authors, 6 keywords, 9 MeSH terms, 1 funder, 80 references.

Cite

This paper

Tian, R., & Scheffler, K. (2026). Group-Patch Joint Compression: Compressing Dynamic B<sub>0</sub> and Static RF Spatial Modulations Across k-Space Subregion Groups for Highly Accelerated MRI. Magnetic resonance in medicine, 96(5), 2106-2126. https://doi.org/10.1002/mrm.70439

BibTeX

@article{tian2026group,
author = {Tian, Rui and Scheffler, Klaus},
title = {{Group-Patch Joint Compression: Compressing Dynamic B\<sub\>0\</sub\> and Static RF Spatial Modulations Across k-Space Subregion Groups for Highly Accelerated MRI}},
journal = {Magnetic resonance in medicine},
year = {2026},
month = jul,
volume = {96},
number = {5},
pages = {2106--2126},
publisher = {Wiley},
issn = {0740-3194},
doi = {10.1002/mrm.70439},
url = {https://doi.org/10.1002/mrm.70439},
pmid = {42490765},
pmcid = {PMC13527271}
}

RIS

TY - JOUR
AU - Tian, Rui
AU - Scheffler, Klaus
TI - Group-Patch Joint Compression: Compressing Dynamic B<sub>0</sub> and Static RF Spatial Modulations Across k-Space Subregion Groups for Highly Accelerated MRI
T2 - Magnetic resonance in medicine
J2 - Magn Reson Med
PY - 2026
DA - 2026/07/23
VL - 96
IS - 5
SP - 2106
EP - 2126
SN - 0740-3194
PB - Wiley
DO - 10.1002/mrm.70439
UR - https://doi.org/10.1002/mrm.70439
LA - en
ER -

CSL-JSON

{
"id": "10.1002/mrm.70439",
"type": "article-journal",
"title": "Group-Patch Joint Compression: Compressing Dynamic B<sub>0</sub> and Static RF Spatial Modulations Across k-Space Subregion Groups for Highly Accelerated MRI",
"container-title": "Magnetic resonance in medicine",
"author": [
{
"family": "Tian",
"given": "Rui"
},
{
"family": "Scheffler",
"given": "Klaus"
}
],
"container-title-short": "Magn Reson Med",
"volume": "96",
"issue": "5",
"page": "2106-2126",
"DOI": "10.1002/mrm.70439",
"PMID": "42490765",
"PMCID": "PMC13527271",
"ISSN": "0740-3194",
"publisher": "Wiley",
"URL": "https://doi.org/10.1002/mrm.70439",
"language": "en",
"issued": {
"date-parts": [
[
2026,
7,
23
]
]
}
}

The tracing map gets a citation of its own once an author has validated it and it has a DOI.

Similar papers

The papers with a page that share the most with this one: the tools found in their code, their categories, datasets, cited references and authors, the rarest counting most.

[1] doi:10.1002/mrm.70469 [code]
Laterally Oscillating Trajectory for Undersampling Slices: LOTUS.
Journal: Magnetic resonance in medicine
In common: Image Processing Toolbox, structural MRI / diffusion, 8 references
[2] doi:10.1002/mrm.70433 [code]
Single-Shot 2D Radial Echo Planar Imaging for Functional MRI.
Journal: Magnetic resonance in medicine
In common: 5 references
[3] doi:10.1002/mrm.70407 [code]
Time-Conditioned Zero-Shot Self-Supervised Reconstruction for Accelerated 3D Ultra-Low-Field MRI.
Journal: Magnetic resonance in medicine
In common: structural MRI / diffusion, 4 references
[4] doi:10.1002/hbm.70553 [code]
Axon Diameter Mapping in the Living Human Brain with Ultra-High-Gradient Diffusion MRI at 500 mT/m Gradient Strength.
Journal: Human brain mapping
In common: Image Processing Toolbox, structural MRI / diffusion, 3 references
[5] doi:10.1109/tmi.2026.3664328 [code]
Axon Diameter Mapping From Myelin Water Diffusion MRI.
Journal: IEEE transactions on medical imaging
In common: Image Processing Toolbox, structural MRI / diffusion, 3 references
[6] doi:10.1002/mrm.70496 [code]
A Deep Nonlinear Subspace Modeling and Reconstruction for Diffusion-Weighted Imaging Using Denoising Auto-Encoder.
Journal: Magnetic resonance in medicine
In common: methods / tools, structural MRI / diffusion, 3 references
[7] doi:10.1002/mrm.70490 [code]
Dependence of the Extra-Cellular Diffusion Coefficient on the Fractions of Neurites and Cell Bodies in Gray Matter.
Journal: Magnetic resonance in medicine
In common: structural MRI / diffusion, 3 references
[8] doi:10.1002/nbm.70305 [code]
Reducing Noise Induced by Cardiac Pulsatility in Brain Maps of R<sub>2</sub>* and Magnetic Susceptibility Using Tailored k-space Sampling.
Journal: NMR in biomedicine
In common: Image Processing Toolbox, structural MRI / diffusion, 2 references
[9] doi:10.7554/elife.92805 [code]
Brain-wide mapping of layer-specific functional connectivity in the human cortex at 3T using draining-vein-suppressed fMRI.
Journal: eLife
In common: Image Processing Toolbox, 2 references
[10] doi:10.1002/mrm.70416 [code]
Fast and Robust Diffusion Posterior Sampling for MR Image Reconstruction Using the Preconditioned Unadjusted Langevin Algorithm.
Journal: Magnetic resonance in medicine
In common: methods / tools, structural MRI / diffusion, 2 references

Contribute

The authors of this paper can claim it, correct its record and validate its tracing map, and the maintainers of its code (its owner, or a public member of its organization) correct what it says of their repository; anyone signed in can ask for its removal. Every request goes to OSCR's own machine, which answers it; your account page follows them.

Sign in with ORCID to claim this paper as one of its authors, correct its record or validate its tracing map: when the paper's metadata lists your ORCID iD, you are recognized at once. Maintainers of its code: sign in with GitHub, then claim the repository on your account page.

Request its removal

To ask OSCR to remove this record, the copies of its authors' scripts or its tracing map, use the removal request page: signed in, you say who you are, what to remove and why, then review and confirm the request. Published rules decide every request (how).

Discussion, reproductions, activity

Discussion: questions and error reports about this paper and its code, from signed-in readers and its authors. It opens with sign-in.

Reproductions: reports from readers who ran the authors' code: what they reproduced, with which environment, commit and data. It opens with sign-in.

Activity: what happens around this paper: new versions of its record, its map's validation, discussions and reproductions. It opens with sign-in.