diff --git a/AtlasViewerGUI.m b/AtlasViewerGUI.m index 9c68835..8548ae1 100644 --- a/AtlasViewerGUI.m +++ b/AtlasViewerGUI.m @@ -57,8 +57,10 @@ function ParseArgs(args) atlasViewer.dirnameAtlas = getAtlasDir(args); end if length(args)>3 - atlasViewer.handles.dataTree = args{4}; - atlasViewer.dataTree = get(atlasViewer.handles.dataTree, 'userdata'); + if sum(isgraphics(args{4})) + atlasViewer.handles.dataTree = args{4}; + atlasViewer.dataTree = get(atlasViewer.handles.dataTree, 'userdata'); + end end % Change current folder to dirnameSubj and load data @@ -651,29 +653,12 @@ function pushbuttonRegisterProbeToSurface_Callback(hObject, eventdata, handles) refpts = set_eeg_active_pts(refpts, [], false); % Finish registration -if isPreRegisteredProbe(probe, refpts) +if 0 %isPreRegisteredProbe(probe, refpts) - % check if probe is already on the surface - probe_on_surface = true; - if ~isempty(atlasViewer.probe.registration.al) - for u = 1:size(atlasViewer.probe.registration.al,1) - ii = ismember(atlasViewer.refpts.labels,atlasViewer.probe.registration.al{u,2}); - refpt_pos = atlasViewer.refpts.pos(ii,:); - opt_pos = atlasViewer.probe.optpos_reg(atlasViewer.probe.registration.al{u,1},:); - dist_refpt_to_opt = sqrt(sum(refpt_pos-opt_pos).^2); - if dist_refpt_to_opt > 1 - probe_on_surface = false; - break; - end - end - end - % Register probe by simply pulling (or pushing) optodes toward surface % toward (or away from) center of head. - if ~ probe_on_surface probe = pullProbeToHeadsurf(probe, headobj); probe.hOptodesIdx = 1; - end else @@ -1281,15 +1266,13 @@ function menuItemGenerateMCInput_Callback(hObject, eventdata, handles) imgrecon = atlasViewer.imgrecon; hbconc = atlasViewer.hbconc; probe = atlasViewer.probe; -headvol = atlasViewer.headvol; pialsurf = atlasViewer.pialsurf; -headsurf = atlasViewer.headsurf; dirnameSubj = atlasViewer.dirnameSubj; axesv = atlasViewer.axesv; try if isempty(eventdata) || strcmp(eventdata.EventName,'Action') - fwmodel = genSensitivityProfile(fwmodel,probe,headvol,pialsurf,headsurf,dirnameSubj); + fwmodel = genSensitivityProfile(fwmodel, probe, dirnameSubj); if isempty(fwmodel.Adot) return; end @@ -1512,7 +1495,7 @@ function menuItemFs2Viewer_Callback(hObject, eventdata, handles) % -------------------------------------------------------------------- -function menuItemShowRefpts_Callback(hObject, eventdata, handles) +function menuItemShowRefpts_Callback(hObject, ~, ~) global atlasViewer switch(get(hObject, 'tag')) @@ -1592,7 +1575,6 @@ function pushbuttonCopyFigure_Callback(~, ~, ~) axis equal axis vis3d set(gca, 'unit','normalized'); -p = get(gca, 'position'); set(gca, 'unit','normalized', 'position', [.20, .30, .40, .40]); % colormap is a propery of figure not axes. Since we don't want to @@ -1666,10 +1648,10 @@ function menuItemLighting_Callback(hObject, eventdata, handles) axesv = atlasViewer.axesv; -if strcmp(get(hObject,'checked'), 'on'); +if strcmp(get(hObject,'checked'), 'on') set(hObject,'checked', 'off'); val=0; -elseif strcmp(get(hObject,'checked'), 'off'); +elseif strcmp(get(hObject,'checked'), 'off') set(hObject,'checked', 'on'); val=1; end @@ -1686,12 +1668,12 @@ function menuItemLighting_Callback(hObject, eventdata, handles) % -------------------------------------------------------------------- -function menuItemFindRefpts_Callback(hObject, eventdata, handles) +function menuItemFindRefpts_Callback(~, ~, ~) FindRefptsGUI(); % -------------------------------------------------------------------- -function menuProbePlacementVariation_Callback(hObject, eventdata, handles) +function menuProbePlacementVariation_Callback(~, ~, ~) plotProbePlacementVariation(); @@ -1958,7 +1940,7 @@ function menuCalcOptProp_Callback(~, ~, ~) set(hLabelsSurf,'FaceVertexCData',faceVertexCData); faceVertexAlphaData(iFace(iFaceMin)) = ones(length(iFace(iFaceMin)),1); set(hLabelsSurf,'FaceVertexAlphaData',faceVertexAlphaData); - iFaces = [iFaces iFace(iFaceMin)]; + iFaces = [iFaces, iFace(iFaceMin)]; end if all(ishandles(hProjectionRays)) set(hProjectionRays,'color','k'); @@ -2245,13 +2227,8 @@ function menuItemGetSensitivityatMNICoordinates_Callback(~, ~, handles) % -------------------------------------------------------------------- -function checkbox_Display_MNI_Projection_Callback(hObject, eventdata, handles) +function checkbox_Display_MNI_Projection_Callback(hObject, ~, ~) global atlasViewer - -fwmodel = atlasViewer.fwmodel; -headvol = atlasViewer.headvol; -labelssurf = atlasViewer.labelssurf; -headvol = atlasViewer.headvol; headsurf = atlasViewer.headsurf; if get(hObject,'value')==1 % if checkbox is checked @@ -2268,12 +2245,12 @@ function checkbox_Display_MNI_Projection_Callback(hObject, eventdata, handles) foo = num2str(coordinate_mni); fooi = sprintf(foo(1,:)); if no_mni == 1 - defaultanswer = {[fooi]}; + defaultanswer = {fooi}; else for i = 2:no_mni foon = ['; ' sprintf(foo(i,:))]; defaultanswer = {[fooi foon]}; - fooi = [fooi foon]; + fooi = [fooi, foon]; end %defaultanswer = {[sprintf(foo(1,:)) ';' sprintf(foo(2,:))]} end else @@ -2300,7 +2277,6 @@ function checkbox_Display_MNI_Projection_Callback(hObject, eventdata, handles) h2 = findobj('Marker','o'); set(h2,'Visible','off'); - headvol = atlasViewer.headvol; vertices = atlasViewer.labelssurf.mesh.vertices; % Project MNI in MC space to head surface and pial surface @@ -2379,7 +2355,6 @@ function menuItemSaveAnatomy_Callback(~, ~, ~) pialsurf = atlasViewer.pialsurf; labelssurf = atlasViewer.labelssurf; refpts = atlasViewer.refpts; -dirnameSubj = atlasViewer.dirnameSubj; saveHeadvol(headvol); saveHeadsurf(headsurf, headvol.T_2mc); @@ -2390,12 +2365,11 @@ function menuItemSaveAnatomy_Callback(~, ~, ~) % -------------------------------------------------------------------- -function menuItemLoadPrecalculatedProfile_Callback(hObject, eventdata, handles) +function menuItemLoadPrecalculatedProfile_Callback(~, ~, handles) global atlasViewer fwmodel = atlasViewer.fwmodel; headvol = atlasViewer.headvol; -pialsurf = atlasViewer.pialsurf; probe = atlasViewer.probe; T_vol2mc = headvol.T_2mc; @@ -2482,7 +2456,7 @@ function popupmenuImageDisplay_Callback(hObject, ~, handles) % -------------------------------------------------------------------- -function popupmenuImageDisplay_CreateFcn(hObject, eventdata, handles) +function popupmenuImageDisplay_CreateFcn(hObject, ~, ~) global popupmenuorder; popupmenuorder = struct(... diff --git a/ForwardModel/genSensitivityProfile.m b/ForwardModel/genSensitivityProfile.m index 16f4a61..64ebba4 100644 --- a/ForwardModel/genSensitivityProfile.m +++ b/ForwardModel/genSensitivityProfile.m @@ -1,4 +1,4 @@ -function fwmodel = genSensitivityProfile(fwmodel,probe,headvol,pialsurf,headsurf,dirnameSubj) +function fwmodel = genSensitivityProfile(fwmodel, probe, dirnameSubj) % Compute volume sensitivity matrix from Monte Carlo simulations and then % projects volume sensitivity onto pial mesh to obtain surface sensitivity @@ -12,6 +12,9 @@ % usage: % call function Adot in subject directory or registered atlas % directory. +global logger + +logger = InitLogger(logger); fwmodel = resetSensitivity(fwmodel,probe,dirnameSubj); if dirnameSubj(end)~='/' && dirnameSubj(end)~='\' @@ -76,13 +79,6 @@ Adot = single(zeros(nMeas,nNode,nWav)); % Adot_scalp = single(zeros(nMeas,nNode_scalp,nWav)); - - % We don't want to reserve the whole matrix in memory (too large). - % rather append each channel's sensitivity to a file. - A = single(zeros(nx,ny,nz)); - As = single(zeros(nx,ny,nz)); - Ad = single(zeros(nx,ny,nz)); - if fwmodel.AdotVolFlag fprintf('Warning: option to generate Adot 3 pt file enabled - May run out of memory\n'); fid1 = fopen([dirnameOut 'AdotVol.3pt'],'wb'); @@ -108,11 +104,11 @@ % load 2pt for given measurement from mc2 mcextreme output iS = probe.ml(iM,1); fileS = sprintf('%sfw%d.s%d.%s', dirnameOut, iW, iS, mc_output_ext); - As = loadMCFuncPtr(fileS, [nx ny nz 1]); + As = loadMCFuncPtr(fileS, [nx, ny, nz, 1]); iD = probe.ml(iM,2); fileD = sprintf('%sfw%d.d%d.%s', dirnameOut, iW, iD, mc_output_ext); - Ad = loadMCFuncPtr(fileD, [nx ny nz 1]); + Ad = loadMCFuncPtr(fileD, [nx, ny, nz, 1]); if fwmodel.normalizeFluence @@ -134,7 +130,7 @@ if sum_p~=0 As(idx_p) = As(idx_p) * (1 - abs(sum_n)) / sum_p; else - disp(sprintf('No photons launched into tissue from Src %d',iS)) + logger.Write('No photons launched into tissue from Src %d\n',iS) As(idx_p) = 0; end @@ -156,7 +152,7 @@ if sum_p~=0 Ad(idx_p) = Ad(idx_p) * (1 - abs(sum_n)) / sum_p; else - fprintf('No photons launched into tissue form Det %d\n',iD) + logger.Write('No photons launched into tissue form Det %d\n',iD) Ad(idx_p) = 0; end @@ -186,12 +182,18 @@ normfactor_d = As(floor(detpos{1}),floor(detpos{2}),floor(detpos{3})); end + [~, fs] = fileparts(fileS); + [~, fd] = fileparts(fileD); + logger.Write('Ch %d. Loading fluence files [%s, %s]\n', iM, fs, fd); + logger.Write('Generating sensitivity for ch %d = [%d, %d, %d], (s%d pos = [%0.1f, %0.1f, %0.1f], d%d pos = [%0.1f, %0.1f, %0.1f])\n', ... + iM, probe.ml(iM,1), probe.ml(iM,2), iW, iS, srcpos{1}, srcpos{2}, srcpos{3}, iD, detpos{1}, detpos{2}, detpos{3}); + % Get Adot normfactor = (normfactor_s + normfactor_d)/2; if normfactor~=0 A = (As.*Ad)/normfactor; else - fprintf('No photons detected between Src %d and Det %d',iS,iD); + logger.Write('No photons detected between Src %d and Det %d\n',iS,iD); A = zeros(size(As)); end if fwmodel.AdotVolFlag @@ -203,6 +205,9 @@ % memory function only available on windows memory; end + + + logger.Write('\n'); % Extract sensitivity values for the nodes of the low res mesh % (Adot) from the sensitivity matrix volume (A) diff --git a/ForwardModel/getFwmodel.m b/ForwardModel/getFwmodel.m index 893eaf7..c9517c5 100644 --- a/ForwardModel/getFwmodel.m +++ b/ForwardModel/getFwmodel.m @@ -123,33 +123,7 @@ % Look at MC input files for forward model params only % if not using precalculated fluence files if isempty(fwmodel.fluenceProfFnames) - - % Wavelength 1 - if exist([dirnameOut, 'fw1.s1.inp'],'file') - config = read_tMCimg_inp([dirnameOut, 'fw1.s1.inp']); - fwmodel.nphotons = config.phot_num; - for ii=1:size(config.tiss_prop,1) - fwmodel.headvol.tiss_prop(ii).scattering(:,1) = config.tiss_prop(ii,1); - fwmodel.headvol.tiss_prop(ii).anisotropy(:,1) = config.tiss_prop(ii,2); - fwmodel.headvol.tiss_prop(ii).absorption(:,1) = config.tiss_prop(ii,3); - fwmodel.headvol.tiss_prop(ii).refraction(:,1) = config.tiss_prop(ii,4); - end - fwmodel.nWavelengths = 1; - end - - % Wavelength 2 - if exist([dirnameOut, 'fw2.s1.inp'],'file') - config = read_tMCimg_inp([dirnameOut, 'fw2.s1.inp']); - fwmodel.nphotons = config.phot_num; - for ii=1:size(config.tiss_prop,1) - fwmodel.headvol.tiss_prop(ii).scattering(:,2) = config.tiss_prop(ii,1); - fwmodel.headvol.tiss_prop(ii).anisotropy(:,2) = config.tiss_prop(ii,2); - fwmodel.headvol.tiss_prop(ii).absorption(:,2) = config.tiss_prop(ii,3); - fwmodel.headvol.tiss_prop(ii).refraction(:,2) = config.tiss_prop(ii,4); - end - fwmodel.nWavelengths = 2; - end - + fwmodel.nWavelengths = length(fwmodel.headvol.tiss_prop(1).scattering); end end diff --git a/ForwardModel/get_tiss_prop.m b/ForwardModel/get_tiss_prop.m index 244176c..fffea06 100644 --- a/ForwardModel/get_tiss_prop.m +++ b/ForwardModel/get_tiss_prop.m @@ -46,38 +46,48 @@ tiss_prop = struct([]); -% Default tissue property values -SCATTERING_SKIN_DEF_VAL = 0.6600; -SCATTERING_SKULL_DEF_VAL = 0.8600; -SCATTERING_DM_DEF_VAL = 0.6600; -SCATTERING_CSF_DEF_VAL = 0.0100; -SCATTERING_GM_DEF_VAL = 1.1000; -SCATTERING_WM_DEF_VAL = 1.1000; -SCATTERING_OTHER_DEF_VAL = 0.8600; +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% Default tissue properties for wavelength 760nm and 860nm +% Reference: Sensors. 13 February 2023 +% Impact of Anatomical Variability on Sensitivity Pro?le in fNIRS�MRI Integration +% Augusto Bonilauri 1, * , Francesca Sangiuliano Intra 2, Francesca Baglio 3and Giuseppe Baselli 1 +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -ABSORPTION_SKIN_DEF_VAL = 0.0191; -ABSORPTION_SKULL_DEF_VAL = 0.0136; -ABSORPTION_DM_DEF_VAL = 0.0191; -ABSORPTION_CSF_DEF_VAL = 0.0026; -ABSORPTION_GM_DEF_VAL = 0.0186; -ABSORPTION_WM_DEF_VAL = 0.0186; -ABSORPTION_OTHER_DEF_VAL = 0.0191; +% Scattering +SCATTERING_SKIN_DEF_VAL = [0.6727, 0.5818]; +SCATTERING_SKULL_DEF_VAL = [0.8545, 0.7636]; +SCATTERING_DM_DEF_VAL = [0.6600, 0.2727]; +SCATTERING_CSF_DEF_VAL = [0.2727, 0.2727]; +SCATTERING_GM_DEF_VAL = [0.7600, 0.6165]; +SCATTERING_WM_DEF_VAL = [1.0825, 0.9188]; +SCATTERING_OTHER_DEF_VAL = [0.8600, 0.6727]; -ANISOTROPY_SKIN_DEF_VAL = 0.0010; -ANISOTROPY_SKULL_DEF_VAL = 0.0010; -ANISOTROPY_DM_DEF_VAL = 0.0010; -ANISOTROPY_CSF_DEF_VAL = 0.0010; -ANISOTROPY_GM_DEF_VAL = 0.0010; -ANISOTROPY_WM_DEF_VAL = 0.0010; -ANISOTROPY_OTHER_DEF_VAL = 0.0010; +% Absorption +ABSORPTION_SKIN_DEF_VAL = [0.0170, 0.0019]; +ABSORPTION_SKULL_DEF_VAL = [0.0116, 0.0139]; +ABSORPTION_DM_DEF_VAL = [0.0191, 0.6727]; +ABSORPTION_CSF_DEF_VAL = [0.0040, 0.0040]; +ABSORPTION_GM_DEF_VAL = [0.0180, 0.0192]; +ABSORPTION_WM_DEF_VAL = [0.0167, 0.0208]; +ABSORPTION_OTHER_DEF_VAL = [0.0180, 0.6727]; -REFRACTION_SKIN_DEF_VAL = 1.0000; -REFRACTION_SKULL_DEF_VAL = 1.0000; -REFRACTION_DM_DEF_VAL = 1.0000; -REFRACTION_CSF_DEF_VAL = 1.0000; -REFRACTION_GM_DEF_VAL = 1.0000; -REFRACTION_WM_DEF_VAL = 1.0000; -REFRACTION_OTHER_DEF_VAL = 1.0000; +% Anisotropy +ANISOTROPY_SKIN_DEF_VAL = [0.0010, 0.0010]; +ANISOTROPY_SKULL_DEF_VAL = [0.0010, 0.0010]; +ANISOTROPY_DM_DEF_VAL = [0.0010, 0.0010]; +ANISOTROPY_CSF_DEF_VAL = [0.0010, 0.0010]; +ANISOTROPY_GM_DEF_VAL = [0.0010, 0.0010]; +ANISOTROPY_WM_DEF_VAL = [0.0010, 0.0010]; +ANISOTROPY_OTHER_DEF_VAL = [0.0010, 0.0010]; + +% Refraction +REFRACTION_SKIN_DEF_VAL = [1.0000, 1.0000]; +REFRACTION_SKULL_DEF_VAL = [1.0000, 1.0000]; +REFRACTION_DM_DEF_VAL = [1.0000, 1.0000]; +REFRACTION_CSF_DEF_VAL = [1.0000, 1.0000]; +REFRACTION_GM_DEF_VAL = [1.0000, 1.0000]; +REFRACTION_WM_DEF_VAL = [1.0000, 1.0000]; +REFRACTION_OTHER_DEF_VAL = [1.0000, 1.0000]; %%% Extract args @@ -143,7 +153,7 @@ %%% Parse tissue names and tissue property names and find their values -propval = []; +propval = {}; m=length(propnames); n=length(names); for i=1:m @@ -152,84 +162,84 @@ for j=1:n switch lower(names{j}) case {'skin', 'scalp'} - propval(j,i) = ANISOTROPY_SKIN_DEF_VAL; + propval{j,i} = ANISOTROPY_SKIN_DEF_VAL; case {'skull', 'bone'} - propval(j,i) = ANISOTROPY_SKULL_DEF_VAL; + propval{j,i} = ANISOTROPY_SKULL_DEF_VAL; case {'dura', 'dura mater', 'dm'} - propval(j,i) = ANISOTROPY_DM_DEF_VAL; + propval{j,i} = ANISOTROPY_DM_DEF_VAL; case {'csf', 'cerebral spinal fluid'} - propval(j,i) = ANISOTROPY_CSF_DEF_VAL; + propval{j,i} = ANISOTROPY_CSF_DEF_VAL; case {'gm', 'gray matter', 'brain'} - propval(j,i) = ANISOTROPY_GM_DEF_VAL; + propval{j,i} = ANISOTROPY_GM_DEF_VAL; case {'wm', 'white matter'} - propval(j,i) = ANISOTROPY_WM_DEF_VAL; + propval{j,i} = ANISOTROPY_WM_DEF_VAL; case 'other' - propval(j,i) = ANISOTROPY_OTHER_DEF_VAL; + propval{j,i} = ANISOTROPY_OTHER_DEF_VAL; otherwise - propval(j,i) = ANISOTROPY_OTHER_DEF_VAL; + propval{j,i} = ANISOTROPY_OTHER_DEF_VAL; end end case 'scattering' for j=1:n switch lower(names{j}) case {'skin', 'scalp'} - propval(j,i) = SCATTERING_SKIN_DEF_VAL; + propval{j,i} = SCATTERING_SKIN_DEF_VAL; case {'skull', 'bone'} - propval(j,i) = SCATTERING_SKULL_DEF_VAL; + propval{j,i} = SCATTERING_SKULL_DEF_VAL; case {'dura', 'dura mater', 'dm'} - propval(j,i) = SCATTERING_DM_DEF_VAL; + propval{j,i} = SCATTERING_DM_DEF_VAL; case {'csf', 'cerebral spinal fluid'} - propval(j,i) = SCATTERING_CSF_DEF_VAL; + propval{j,i} = SCATTERING_CSF_DEF_VAL; case {'gm', 'gray matter', 'brain'} - propval(j,i) = SCATTERING_GM_DEF_VAL; + propval{j,i} = SCATTERING_GM_DEF_VAL; case {'wm', 'white matter'} - propval(j,i) = SCATTERING_WM_DEF_VAL; + propval{j,i} = SCATTERING_WM_DEF_VAL; case 'other' - propval(j,i) = SCATTERING_OTHER_DEF_VAL; + propval{j,i} = SCATTERING_OTHER_DEF_VAL; otherwise - propval(j,i) = SCATTERING_OTHER_DEF_VAL; + propval{j,i} = SCATTERING_OTHER_DEF_VAL; end end case 'absorption' for j=1:n switch lower(names{j}) case {'skin', 'scalp'} - propval(j,i) = ABSORPTION_SKIN_DEF_VAL; + propval{j,i} = ABSORPTION_SKIN_DEF_VAL; case {'skull', 'bone'} - propval(j,i) = ABSORPTION_SKULL_DEF_VAL; + propval{j,i} = ABSORPTION_SKULL_DEF_VAL; case {'dura', 'dura mater', 'dm'} - propval(j,i) = ABSORPTION_DM_DEF_VAL; + propval{j,i} = ABSORPTION_DM_DEF_VAL; case {'csf', 'cerebral spinal fluid'} - propval(j,i) = ABSORPTION_CSF_DEF_VAL; + propval{j,i} = ABSORPTION_CSF_DEF_VAL; case {'gm', 'gray matter', 'brain'} - propval(j,i) = ABSORPTION_GM_DEF_VAL; + propval{j,i} = ABSORPTION_GM_DEF_VAL; case {'wm', 'white matter'} - propval(j,i) = ABSORPTION_WM_DEF_VAL; + propval{j,i} = ABSORPTION_WM_DEF_VAL; case 'other' - propval(j,i) = ABSORPTION_OTHER_DEF_VAL; + propval{j,i} = ABSORPTION_OTHER_DEF_VAL; otherwise - propval(j,i) = ABSORPTION_OTHER_DEF_VAL; + propval{j,i} = ABSORPTION_OTHER_DEF_VAL; end end case 'refraction' for j=1:n switch lower(names{j}) case {'skin', 'scalp'} - propval(j,i) = REFRACTION_SKIN_DEF_VAL; + propval{j,i} = REFRACTION_SKIN_DEF_VAL; case {'skull', 'bone'} - propval(j,i) = REFRACTION_SKULL_DEF_VAL; + propval{j,i} = REFRACTION_SKULL_DEF_VAL; case {'dura', 'dura mater', 'dm'} - propval(j,i) = REFRACTION_DM_DEF_VAL; + propval{j,i} = REFRACTION_DM_DEF_VAL; case {'csf', 'cerebral spinal fluid'} - propval(j,i) = REFRACTION_CSF_DEF_VAL; + propval{j,i} = REFRACTION_CSF_DEF_VAL; case {'gm', 'gray matter', 'brain'} - propval(j,i) = REFRACTION_GM_DEF_VAL; + propval{j,i} = REFRACTION_GM_DEF_VAL; case {'wm', 'white matter'} - propval(j,i) = REFRACTION_WM_DEF_VAL; + propval{j,i} = REFRACTION_WM_DEF_VAL; case 'other' - propval(j,i) = REFRACTION_OTHER_DEF_VAL; + propval{j,i} = REFRACTION_OTHER_DEF_VAL; otherwise - propval(j,i) = REFRACTION_OTHER_DEF_VAL; + propval{j,i} = REFRACTION_OTHER_DEF_VAL; end end end @@ -240,6 +250,6 @@ for j=1:length(names) tiss_prop(j).name = names{j}; for i=1:length(propnames) - eval(sprintf('tiss_prop(j).%s = propval(j,i);',propnames{i})); + eval(sprintf('tiss_prop(j).%s = propval{j,i};',propnames{i})); end end diff --git a/Group/DataTree/AcquiredData/DataFiles/DataFilesClass.m b/Group/DataTree/AcquiredData/DataFiles/DataFilesClass.m index a118540..fab1709 100644 --- a/Group/DataTree/AcquiredData/DataFiles/DataFilesClass.m +++ b/Group/DataTree/AcquiredData/DataFiles/DataFilesClass.m @@ -16,6 +16,8 @@ properties (Access = private) lookupTable excludedFolders + savedFilename + changed end methods @@ -30,6 +32,8 @@ obj.err = -1; obj.errmsg = {}; obj.dirFormats = struct('type',0, 'choices',{{}}); + obj.lookupTable = []; + obj.changed = false; global logger global cfg @@ -38,7 +42,6 @@ cfg = InitConfig(cfg); obj.logger = logger; - obj.lookupTable = []; obj.excludedFolders = {}; skipconfigfile = false; @@ -76,7 +79,7 @@ obj.rootdir = filesepStandard(obj.rootdir,'full'); % Configuration parameters - obj.config = struct('RegressionTestActive','', 'AskToFixNameConflicts',1, 'DerivedDataDir',''); + obj.config = struct('RegressionTestActive','', 'AskToFixNameConflicts',1, 'DerivedDataDir','', 'DerivedDataRootDir',''); if skipconfigfile==false str = cfg.GetValue('Regression Test Active'); if strcmp(str,'true') @@ -110,108 +113,295 @@ [obj.rootdir, 'imagerecon']; [obj.rootdir, 'anatomical']; }; + obj.config.DerivedDataRootDir = [obj.rootdir, obj.config.DerivedDataDir, '/', f, '/']; + obj.savedFilename = [obj.config.DerivedDataRootDir, 'DatasetFiles', '_', obj.filetype, '.mat']; if nargin==0 return; end obj.err = 0; - - obj.InitFolderFromats(); + obj.InitFolderFormats(); obj.GetDataSet(); + obj.SaveDataSet(); + end + + + + % ----------------------------------------------------------------------------------- + function Compare(obj) + if obj.IsEmpty() + return + end + obj2 = DataFilesClass(); + obj2.rootdir = obj.rootdir; + obj2.filetype = obj.filetype; + obj2.dirFormats = obj.dirFormats; + obj2.GetDataSet('skipsaveddata') + + %%%% Create list of all the saved names for sorting later + filesSaved = cell(length(obj.files),1); + filesSavedErr = cell(length(obj.filesErr),1); + for ii = 1:length(obj.files) + filesSaved{ii} = obj.files(ii).name; + end + for ii = 1:length(obj.filesErr) + filesSavedErr{ii} = obj.filesErr(ii).name; + end + + %%%% Move any files from saved error list that have been updated to nonerror list + N = length(obj.filesErr); + for ii = N:-1:1 + for jj = 1:length(obj2.files) + if strcmp(obj.filesErr(ii).name, obj2.files(jj).name) + if obj.filesErr(ii).datenum ~= obj2.files(jj).datenum + obj.filesErr(ii) = []; + filesSavedErr{ii} = []; + obj.files(end+1) = obj2.files(jj); + if ~obj2.files(jj).isdir + obj.nfiles = obj.nfiles+1; + end + filesSaved{end+1} = obj2.files(jj).name; %#ok + obj.logger.Write('%d. Mismatch in dates for saved data in %s. Will re-validate this file.\n', jj, obj.filesErr(ii).name); + obj.changed = true; + end + break + end + end + end + + + %%%% Now check saved object nonerror list + N = length(obj.files); + for ii = N:-1:1 + found = false; + for jj = 1:length(obj2.files) + if strcmp(obj.files(ii).name, obj2.files(jj).name) + found = true; + if obj.files(ii).datenum ~= obj2.files(jj).datenum + obj.files(ii) = obj2.files(jj).copy(); + obj.logger.Write('%d. Mismatch in dates for saved data in %s. Will re-validate this file.\n', jj, obj.files(ii).name); + obj.changed = true; + end + break + end + end + if ~found + obj.logger.Write('%d. Deleting file %s from saved data.\n', ii, obj.files(ii).name); + if ~obj.files(ii).isdir + obj.nfiles = obj.nfiles-1; + end + obj.files(ii) = []; + filesSaved(ii) = []; + obj.changed = true; + end + end + + % We don't need to check if files were added if there's no + % change and the number of current files found is same as saved + % number of files. Exist early. + if obj.changed == false + if (length(obj.files) + length(obj.filesErr)) == length(obj2.files) + return + end + end + + + %%%% Find which files were added to current dataset and add + %%%% them to saved object + for ii = 1:length(obj2.files) + found = false; + for jj = 1:length(obj.files) + if strcmp(obj2.files(ii).name, obj.files(jj).name) + found = true; + break + end + end + + % If the current file filesCurr(ii) is not found in saved dataset then it's a new file. + % Add new file from current dataset to saved dataset + if ~found + obj.files(end+1) = obj2.files(ii); + if ~obj2.files(ii).isdir + obj.nfiles = obj.nfiles+1; + end + filesSaved{end+1} = obj2.files(ii).name; %#ok + obj.logger.Write('%d. Adding new file %s to saved data.\n', jj, obj2.files(ii).name); + obj.changed = true; + end + end + + + % Sort non-error files + if obj.changed + [~, order] = sort(filesSaved); + obj.files = obj.files(order); + obj.logger.Write('\n'); + end + + end + + + + % ----------------------------------------------------------------------------------- + function SaveDataSet(obj) + if ~obj.changed + return + end + obj.changed = false; + obj.logger.Write('Saving DataFilesClass object in %s\n', obj.savedFilename); + if ~exist(obj.config.DerivedDataRootDir, 'dir') + mkdir(obj.config.DerivedDataRootDir) + end + save(obj.savedFilename, '-mat', 'obj'); + end + + + + % ----------------------------------------------------------------------------------- + function b = LoadDataSet(obj) + b = false; + if ~exist(obj.savedFilename, 'file') + return + end + dataset = load(obj.savedFilename); + dataset.obj.Compare(); + if isempty(dataset.obj.files) && isempty(dataset.obj.filesErr) + return + end + obj.logger.Write('Loaded saved DataFilesClass object: %s\n', obj.savedFilename); + obj.Copy(dataset.obj); + + obj.ErrorCheck() + obj.ErrorCheckName() + b = true; end % ----------------------------------------------------------------------------------- - function GetDataSet(obj) + function GetDataSet(obj, skipsaveddata) + if ~exist('skipsaveddata','var') + skipsaveddata = ''; + end + if ~strcmp(skipsaveddata, 'skipsaveddata') + if obj.LoadDataSet() + return + end + end + + obj.changed = true; + + mode = ''; + if obj.dirFormats.type == 0 + iFormats = 1:length(obj.dirFormats.choices); + else + iFormats = obj.dirFormats.type; + mode = 'noerrorcheck'; + end if exist(obj.rootdir, 'dir')~=7 error('Invalid subject folder: ''%s''', obj.rootdir); end - - for ii = 1:length(obj.dirFormats.choices) + for ii = iFormats obj.InitLookupTable(); obj.FindDataSet(ii); + if strcmp(mode,'noerrorcheck') + continue + end + if isempty(obj.files) + continue + end + % Remove any files that cannot pass the basic test of loading % its data obj.ErrorCheck(); if ~isempty(obj.files) break end + + % If we're here it means we still have not found a valid dataset in any of the formats looked + % at so far, but obj.dirFormats.type may have been set to non-zero because a invalid dataset + % may have been found by obj.FindDataSet() then flagged as invalid by obj.ErrorCheck(). So make + % sure to reset format type here otherwise wrong format could be associated with saved dataset. + % To reproduce comment out this fix and add ONE error snirf file that does NOT belong in the + % dataset due to being in incompatible format folder structure. BUG FIX: Jay Dubb, May 8, 2024 + obj.dirFormats.type = 0; + end + if strcmp(mode,'noerrorcheck') + return end - obj.ErrorCheckFinal(); end % ----------------------------------------------------------------------------------- - function InitFolderFromats(obj) + function InitFolderFormats(obj) obj.dirFormats.choices = { %%%% 1. Flat #1 { - ['*_run*.', obj.filetype]; + ['*_run*.', obj.filetype]; } %%%% 2. Flat #2 { - ['*.', obj.filetype]; + ['*.', obj.filetype]; } %%%% 3. Subjects { - '*'; - ['*.', obj.filetype]; + '*'; + ['*.', obj.filetype]; } %%%% 4. BIDS #1, sub-