Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Solutions for the high density mesh for MEG with DUNEuro: integration points and memory #365

Merged
Show file tree
Hide file tree
Changes from all commits
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
23 changes: 22 additions & 1 deletion toolbox/forward/bst_duneuro.m
Original file line number Diff line number Diff line change
Expand Up @@ -76,6 +76,26 @@
MegChannels = [MegChannels; iChan, sChan.Loc(:,iInteg)', sChan.Orient(:,iInteg)', sChan.Weight(iInteg)];
end
end
% In the case where the MEG integration points are used
if cfg.UseIntegrationPoint == 0
% loop over the integration Points
% chan_loc = figure_3d('GetChannelPositions', cfg, cfg.iMeg); % <= this function is not sufficient, we need also the weights.
MegChannelsTemp = [];
for iChan = 1 : length(cfg.iMeg)
group = MegChannels(MegChannels(:,1) == iChan,:);
groupPositive = group(group(:,end)>0,:);
groupNegative = group(group(:,end)<0,:);
if ~isempty(groupPositive)
equivalentPositionPostive = sum(repmat(abs(groupPositive(:,end)),[1 3]) .* groupPositive(:,2:4));
MegChannelsTemp = [MegChannelsTemp; iChan equivalentPositionPostive groupPositive(1,5:7) sum(groupPositive(:,end))];
end
if ~isempty(groupNegative)
equivalentPositionNegative = sum(repmat(abs(groupNegative(:,end)),[1 3]) .* groupNegative(:,2:4));
MegChannelsTemp = [MegChannelsTemp; iChan equivalentPositionNegative groupNegative(1,5:7) sum(groupNegative(:,end))];
end
end
MegChannels = MegChannelsTemp;
end
end

%% ===== HEAD MODEL =====
Expand Down Expand Up @@ -392,6 +412,7 @@
fprintf(fid, '[meg]\n');
fprintf(fid, 'intorderadd = %d\n', cfg.MegIntorderadd);
fprintf(fid, 'type = %s\n', cfg.MegType);
fprintf(fid, 'cache.enable = %s\n',bool2str(cfg.EnableCacheMemory) );
% [coils]
fprintf(fid, '[coils]\n');
fprintf(fid, 'filename = %s\n', CoilFile);
Expand Down Expand Up @@ -518,7 +539,7 @@
nbChannel = length(channelIndex);
weighted_B = zeros(nbChannel,size(Bfull,2));
for iCh = 1 : nbChannel
communChannel = find(iCh==MegChannels);
communChannel = find(iCh==MegChannels(:,1));
BcommunChannel = Bfull(communChannel(:),:);
WcommunChannel = MegChannels(communChannel(:), 8: end);
weighted_B(iCh,:) = sum (BcommunChannel.*WcommunChannel,1);
Expand Down
4 changes: 4 additions & 0 deletions toolbox/forward/duneuro_defaults.m
Original file line number Diff line number Diff line change
Expand Up @@ -76,6 +76,10 @@
cfgDef.BstEegLfFile = 'eeg_lf.dat';
cfgDef.BstMegLfFile = 'meg_lf.dat';

% [MEG computation Options]
cfgDef.UseIntegrationPoint = 1;
cfgDef.EnableCacheMemory = 0;
cfgDef.MegPerBlockOfSensor = 0; % ToDo
% Use default values if not set
if (nargin == 0) || isempty(cfg)
cfg = cfgDef;
Expand Down
52 changes: 49 additions & 3 deletions toolbox/forward/panel_duneuro.m
Original file line number Diff line number Diff line change
Expand Up @@ -41,7 +41,8 @@
OPTIONS = sProcess;
% Check if there is only MEG, for simplified model by default
isMegOnly = ~ismember('duneuro', {OPTIONS.EEGMethod, OPTIONS.ECOGMethod, OPTIONS.SEEGMethod});
% PROCESS CALL: panel_duneuro('CreatePanel', sProcess, sFiles)
isMeg = isequal(OPTIONS.MEGMethod, 'duneuro');
% PROCESS CALL: panel_duneuro('CreatePanel', sProcess, sFiles)
else
OPTIONS = sProcess.options.duneuro.Value;
% List of sensors
Expand All @@ -52,6 +53,7 @@
Modalities = setdiff(Modalities, 'MEG');
end
end
isMeg = any(ismember({'MEG', 'MEG MAG', 'MEG GRAD'}, Modalities));
isMegOnly = all(ismember(Modalities, {'MEG', 'MEG MAG', 'MEG GRAD'}));
% Get FEM files
sSubject = bst_get('Subject', sFiles(1).SubjectFile);
Expand Down Expand Up @@ -221,7 +223,29 @@
jCheckSrcForceInGM = [];
end
c.gridy = 3;
jPanelRight.add(jPanelInput, c);
jPanelRight.add(jPanelInput, c);

% ==== PANEL RIGHT: MEG COMPUTATIONS OPTIONS ====
jPanelMegComputationOption = gui_river([1,1], [0,6,6,6], 'MEG computation options');
if isMeg
% Use integration Points, recommended for high mesh density
jCheckUseIntegrationPoint = gui_component('checkbox', jPanelMegComputationOption, 'br', 'Use MEG integration points', [], '', [], []);
% Enable MEG cache memory for high mesh density if users do not
% high memory, or want to use the integration points
jCheckEnableCacheMemory = gui_component('checkbox', jPanelMegComputationOption, 'br', 'Enable cache memory', [], '', [], []);
% Enable the MEG Computation per block of sensors
... jCheckMegPerBlockOfSensor = gui_component('checkbox', jPanelMegComputationOption, 'br', 'Compute per block of sensors [Todo]', [], '', [], []);
% Set jCheckUseIntegrationPoint to 1 as default option
if (OPTIONS.UseIntegrationPoint)
jCheckUseIntegrationPoint.setSelected(1);
end
c.gridy = 4;
jPanelRight.add(jPanelMegComputationOption, c);
else
jCheckUseIntegrationPoint = [];
jCheckEnableCacheMemory = [];
jCheckMegPerBlockOfSensor = [];
end

% ==== PANEL RIGHT: OUTPUT OPTIONS ====
jPanelOutput = gui_river([1,1], [0,6,6,6], 'Output options');
Expand All @@ -230,7 +254,7 @@
if (OPTIONS.BstSaveTransfer)
jCheckSaveTransfer.setSelected(1);
end
c.gridy = 4;
c.gridy = 5;
jPanelRight.add(jPanelOutput, c);

% ===== GLUE =====
Expand Down Expand Up @@ -287,6 +311,9 @@
'jTextSrcShrink', jTextSrcShrink, ...
'jCheckSrcForceInGM', jCheckSrcForceInGM, ...
'jCheckSaveTransfer', jCheckSaveTransfer, ...
'jCheckUseIntegrationPoint', jCheckUseIntegrationPoint,...
'jCheckEnableCacheMemory', jCheckEnableCacheMemory,...
...'jCheckMegPerBlockOfSensor', jCheckMegPerBlockOfSensor,...
'UseTensor', OPTIONS.UseTensor);
ctrl.FemNames = OPTIONS.FemNames;
% Create the BstPanel object that is returned by the function
Expand Down Expand Up @@ -339,6 +366,7 @@ function UpdatePanel(isForced)
jPanelOptSub.setVisible(ExpertMode && jRadioSrcModelSub.isSelected());
jPanelInput.setVisible(ExpertMode);
jPanelOutput.setVisible(ExpertMode);
jPanelMegComputationOption.setVisible(ExpertMode);
% Update expert button
if ExpertMode
jButtonExpert.setText('Hide details');
Expand Down Expand Up @@ -428,6 +456,24 @@ function UpdatePanel(isForced)
end
% Output options
s.BstSaveTransfer = ctrl.jCheckSaveTransfer.isSelected();
% MEG Computation options
if ~isempty(ctrl.jCheckUseIntegrationPoint)
s.UseIntegrationPoint = ctrl.jCheckUseIntegrationPoint.isSelected();
else
s.UseIntegrationPoint = 1;
end

if ~isempty(ctrl.jCheckEnableCacheMemory)
s.EnableCacheMemory = ctrl.jCheckEnableCacheMemory.isSelected();
else
s.EnableCacheMemory = 0;
end

% if ~isempty(ctrl.jCheckMegPerBlockOfSensor)
% s.MegPerBlockOfSensor = ctrl.jCheckMegPerBlockOfSensor.isSelected();
% else
% s.MegPerBlockOfSensor = 0;
% end
end


Expand Down