Coordinates of voxel with maximum intensity - sLORETA

AdamKalina · June 5, 2020, 7:15pm

Hello, I would like to get coordinates of maximal source voxel in sLORETA results indicated in ImageGridAmp.

I have found similar question here: Maximal source coordinates and tried to adopt the solution, but the results are somewhat different (but close) from what I see in exported nifti files or in MriViewer after I press "M" in given time (basically I would like to know how the shortcut works). I dont think so it is problem in coordinates conversion. I even tried to compute norms for sources in ImageGridAmp, but that got me even further from the maximum showed in MriViewer.

I am using MRI volume based unconstrained head models (so I am unable to use scouts and I am not interested in vertex coordinates).

Thanks you for your help,
Adam

Francois · June 8, 2020, 8:11am

I would like to get coordinates of maximal source voxel in sLORETA results indicated in ImageGridAmp. [...] but the results are somewhat different (but close) from what I see in exported nifti files or in MriViewer

The ImageGridAmp values are calculated over a volume grid that you defined in the options of the forward model. It is much more sparse than the volume definition you have in the MRI (up to 256x256x256 voxels). To display it on top of the T1 MRI, we have to re-interpolate these values on this much finer grid of voxels. Therefore what you see in the MRI Viewer (and in the exported .nii files) are not the "real" values computed by the minimum norm method, but their interpolation over the finer grid of MRI voxels. The values are not expected to match exactly, but should be in the same range.

in MriViewer after I press "M" in given time (basically I would like to know how the shortcut works)

github.com

brainstorm-tools/brainstorm3/blob/master/toolbox/gui/figure_mri.m#L611


      
              MouseWheel_Callback(hFig, [], [], 1);

          case 'downarrow'

              MouseWheel_Callback(hFig, [], [], -1);

          % === DATABASE NAVIGATOR ===

          case {'f1', 'f2', 'f3', 'f4', 'f6'}

              bst_figures('NavigatorKeyPress', hFig, keyEvent);

          % CTRL+D : Dock figure

          case 'd'

              if ismember('control', keyEvent.Modifier)

                  isDocked = strcmpi(get(hFig, 'WindowStyle'), 'docked');

                  bst_figures('DockFigure', hFig, ~isDocked);

              end

          % CTRL+I : Save as image

          case 'i'

              if ismember('control', keyEvent.Modifier)

                  out_figure_image(hFig);

              end

          % CTRL+J : Open as image

          case 'j'

              if ismember('control', keyEvent.Modifier)

                  out_figure_image(hFig, 'Viewer');

github.com

brainstorm-tools/brainstorm3/blob/master/toolbox/gui/figure_mri.m#L1174


      
          end

          

          % === LOAD LANDMARKS ===

          % Load landmarks/fiducials

          [sMri, Handles] = LoadLandmarks(sMri, Handles);

          % Save mri and Handles

          GlobalData.Mri(iMri) = sMri;

          bst_figures('SetFigureHandles', hFig, Handles);

          

          % === CONFIGURE SLIDERS ===

          jSliders = [Handles.jSliderSagittal, Handles.jSliderCoronal, Handles.jSliderAxial];

          % Reset sliders callbacks

          java_setcb(jSliders, 'StateChangedCallback', []);

          % Configure each slider

          for i = 1:3

              % Set min and max bounds

              jSliders(i).setMinimum(1);

              jSliders(i).setMaximum(cubeSize(i));

          end

          % Set default location to middle of the volume

          SetLocation('voxel', sMri, Handles, TessInfo(iTess).CutsPosition);

AdamKalina · June 8, 2020, 1:34pm

Oh yeah, the interpolation, how could I forgot. Thank you!

While using breakpoints in Brainstorm code I reconstructed script doing what I need without displaying the MriViewer. I know it is very limited, but can you see any principal mistake?
% sMri = T1
% sResults = full kernel results
% t = desired point of time

isVolumeGrid = 1;
t = 0.0356;

[~,iTime]=min(abs(sResults.Time-t)); %in case desired time is not in TimeVector
tTime=sResults.Time(iTime);

ResultsValues = double(sResults.ImageGridAmp(:, iTime));
NormResultsValues = bst_source_orient([], sResults.nComponents, sResults.GridAtlas, ResultsValues, 'rms'); %counts norm

%MriInterp = grid_interp_mri_seeg(sResults.GridLoc, sMri); %SEEG
MriInterp = grid_interp_mri(sResults.GridLoc, sMri, sResults.SurfaceFile, 1, [], [], 1); % HDEEG

OverlayCube = tess_interp_mri_data(MriInterp, size(sMri.Cube(:,:,:,1)), NormResultsValues, isVolumeGrid);
[valMax, iMax] = max(OverlayCube(:));

%OverlayCube = 256256256
x = mod(iMax,256);
n = round(iMax/256);
z = round(n/256);
y = mod(n,256);

pcs = cs_convert(sMri,'voxel', 'mri', [x,y,z]);

Francois · June 8, 2020, 2:23pm

I'm sorry, this is a bit complicated to review - it is not exactly the purpose of this forum to help you proofread your personal Matlab scripts.
If there is some Brainstorm functional call that doesn't work as expected and you can't understand why, please let us know. With some precise examples and errors we can reproduce, we would be able to help you.

AdamKalina · June 8, 2020, 3:19pm

Oky, sorry and thank you for your help.