[Master Index] [Index for Toolbox]

# zscore_gui

## (Toolbox/zscore_gui.m in BrainStorm 2.0 (Alpha))

### Function Synopsis

    varargin = zscore_gui(action,varargin)

### Help Text

    ZSCORE_GUI - : management of zscore visualization
     function varargin = zscore_gui(action,varargin)

### Cross-Reference Information

This function calls

    bst_color_scheme                   C:\BrainStorm_2001\Toolbox\bst_color_scheme.m
    bst_imaging_display                C:\BrainStorm_2001\Toolbox\bst_imaging_display.m
    bst_layout                         C:\BrainStorm_2001\Toolbox\bst_layout.m
    bst_message_window                 C:\BrainStorm_2001\Toolbox\bst_message_window.m
    dataplot_cb                        C:\BrainStorm_2001\Toolbox\dataplot_cb.m
    findclosest                        C:\BrainStorm_2001\Toolbox\findclosest.m
    inorcol                            C:\BrainStorm_2001\Toolbox\inorcol.m
    patch_swell                        C:\BrainStorm_2001\Toolbox\patch_swell.m
    save_fieldnames                    C:\BrainStorm_2001\Toolbox\save_fieldnames.m
    togglebuttoncolor                  C:\BrainStorm_2001\Toolbox\togglebuttoncolor.m
    vertices_connectivity              C:\BrainStorm_2001\Toolbox\vertices_connectivity.m

This function is called by

    dataplot_cb                        C:\BrainStorm_2001\Toolbox\dataplot_cb.m

### Listing of function C:\BrainStorm_2001\Toolbox\zscore_gui.m

    function varargin = zscore_gui(action,varargin)
    %ZSCORE_GUI - : management of zscore visualization
    % function varargin = zscore_gui(action,varargin)

    %<autobegin> ---------------------- 26-May-2004 11:34:45 -----------------------
    % --------- Automatically Generated Comments Block Using AUTO_COMMENTS ---------
    %
    % CATEGORY: GUI and Related
    %
    % Alphabetical list of external functions (non-Matlab):
    %   toolbox\bst_color_scheme.m
    %   toolbox\bst_imaging_display.m
    %   toolbox\bst_layout.m
    %   toolbox\bst_message_window.m
    %   toolbox\dataplot_cb.m
    %   toolbox\findclosest.m
    %   toolbox\inorcol.m
    %   toolbox\patch_swell.m
    %   toolbox\save_fieldnames.m
    %   toolbox\togglebuttoncolor.m
    %   toolbox\vertices_connectivity.m
    %   toolbox\zscore_gui.m  NOTE: Routine calls itself explicitly
    %
    % Application data and their calls in this file:
    %   'TileType'
    %
    %   setappdata(zscore_fig,'TileType','T')
    %
    %
    % Figure Files opened by this function:
    %   'zscore_gui.fig'
    %
    %   Format of strings below: Type:Style:Tag, "String", CallBack Type and Call
    %   <automatic> callback is <Tag>_Callback by Matlab default
    %
    % Callbacks by figure zscore_gui.fig
    %   uicontrol:edit:Baseline "" uses Callback for zscore_gui DefineBaseline
    %   uicontrol:edit:CorrectionFactor "" uses Callback for zscore_gui ZScoreThreshold
    %   uicontrol:edit:ZScoreThreshold "5" uses Callback for zscore_gui ZScoreThreshold
    %   uicontrol:pushbutton:AnalyzeCorticalMap "Spatiotemp. Segment" uses Callback for
    %     zscore_gui AnalyzeCorticalMap
    %   uicontrol:pushbutton:AverageMaps "Avrg" uses Callback for zscore_gui AverageMaps
    %   uicontrol:pushbutton:pushbutton55 "Load" uses Callback for zscore_gui LoadBaseline
    %   uicontrol:pushbutton:pushbutton56 "Fetch files" uses Callback for zscore_gui ZScoreBatch
    %   uicontrol:pushbutton:Qpushbutton7 "Close" uses Callback for close(gcbf)
    %   uicontrol:pushbutton:STHistogram "Histo" uses Callback for zscore_gui STHistogram
    %   uicontrol:togglebutton:IntSTH "Integ." uses Callback for zscore_gui IntegSTHistogram
    %   uicontrol:togglebutton:ZScore "Z-score" uses Callback for togglebuttoncolor ; zscore_gui ZScore
    %   uicontrol:togglebutton:ZScoreThresholdApply "Apply" uses Callback for
    %     togglebuttoncolor; zscore_gui ZScoreThreshold
    %
    % At Check-in: $Author: Mosher $  $Revision: 12 $  $Date: 5/26/04 10:02a $
    %
    % This software is part of BrainStorm Toolbox Version 2.0 (Alpha) 24-May-2004
    %
    % Principal Investigators and Developers:
    % ** Richard M. Leahy, PhD, Signal & Image Processing Institute,
    %    University of Southern California, Los Angeles, CA
    % ** John C. Mosher, PhD, Biophysics Group,
    %    Los Alamos National Laboratory, Los Alamos, NM
    % ** Sylvain Baillet, PhD, Cognitive Neuroscience & Brain Imaging Laboratory,
    %    CNRS, Hopital de la Salpetriere, Paris, France
    %
    % See BrainStorm website at http://neuroimage.usc.edu for further information.
    %
    % Copyright (c) 2004 BrainStorm by the University of Southern California
    % This software distributed  under the terms of the GNU General Public License
    % as published by the Free Software Foundation. Further details on the GPL
    % license can be found at http://www.gnu.org/copyleft/gpl.html .
    %
    % FOR RESEARCH PURPOSES ONLY. THE SOFTWARE IS PROVIDED "AS IS," AND THE
    % UNIVERSITY OF SOUTHERN CALIFORNIA AND ITS COLLABORATORS DO NOT MAKE ANY
    % WARRANTY, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO WARRANTIES OF
    % MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE, NOR DO THEY ASSUME ANY
    % LIABILITY OR RESPONSIBILITY FOR THE USE OF THIS SOFTWARE.
    %<autoend> ------------------------ 26-May-2004 11:34:45 -----------------------

    global data

    ImageGridTime = data.Display.Time(1):data.Display.SamplingRate:data.Display.Time(end);

    switch(action)

    case 'create'
        zscore_fig = openfig('zscore_gui.fig','reuse');
        setappdata(zscore_fig,'TileType','T')
        bst_color_scheme(zscore_fig)
        bst_layout('align',zscore_fig,2,2,4)
        set(zscore_fig,'visible','on')

        %----------------------------------------------------------------------------

    case 'DefineBaseline'
        hZScore_fig = guihandles(gcbf);
        set(hZScore_fig.ZScore,'Value',1)
        dataplot_cb('ToggleButtonColor',hZScore_fig.ZScore)

        if ~strcmp(get(hZScore_fig.Baseline,'String'),'File')
            set(hZScore_fig.Baseline,'Userdata',0) % FLAG: Baseline defined on current file
        end

        togglebuttoncolor
        zscore_gui('ZScore',1) % Force ZScore computation

        %----------------------------------------------------------------------------

    case 'LoadBaseline' % Define Baseline from a different result file
        cd(Users.STUDIES)
        [resfile, respath] = uigetfile('*result*.mat','Please Select a Result File Containing the Activation Baseline');
        if resfile == 0, return, end
        Zscore_fig = findobj(0,'Tag','tesselation_select');
        hZScore_fig = guihandles(Zscore_fig);
        set(hZScore_fig.LoadBaseline,'userdata',fullfile(respath, resfile));
        set(hZScore_fig.Baseline,'String','File');
        set(hZScore_fig.Baseline,'Userdata',1) % FLAG: Baseline defined on a pre-loaded file

        dataplot_cb DefineBaseline

        %----------------------------------------------------------------------------

    case 'ZScore' % Switches from Z-Score to Absolute current density mapping

        ZScore_fig = gcbf;
        hZScore_fig = guihandles(ZScore_fig);

        switch get(hZScore_fig.ZScore,'Value')
        case 1 % Switch to ZScore map

            % Are we switching from ZScore to Absolute mapping or vice-versa ?
            if get(hZScore_fig.Baseline,'Userdata') == 1  % Baseline is defined from an existing result file, load it
                data.Results = load(get(hZScore_fig.LoadBaseline,'Userdata'));
                FlagBaselineFile = 1;
            else % else use current loaded results
                FlagBaselineFile = 0;
            end

            if isempty(data.Results), return, end % No results were loaded

            test = 0; % By default, consider not computing a new ZScore - check if a previously computation is available -
            % if test is set to 1 anytime below, this will conduct to a new ZScore computation either because the baseline has changed
            % or the computation is conducted on a new set of data.
            if nargin == 1
                varargin{1} = 0;
            end

            if (~isfield(data.Results,'ZScore') | varargin{1} == 1) % No ZScore was defined beforehand or computation is forced
                data.Results.ZScore.Baseline = [-Inf Inf]; % Dummy Score
                test = 1; % Force ZScore computation
            end

            %Get baseline time limits
            if FlagBaselineFile == 1
                set(hZScore_fig.Baseline,'String',...
                    [num2str(ImageGridTime(1)*1000,'%3.1f'),' ',num2str(ImageGridTime(end)*1000,'%3.1f')]);

                baseline = [ImageGridTime(1) ImageGridTime(end)];
            else
                baseline = str2num(get(hZScore_fig.Baseline,'String'))/1000;
            end

            if length(baseline) < 2
                errordlg('Please enter 2 time values to define the baseline interval')
                return
            end

            if isfield(data.Results.ZScore,'FlagBaselineFile') % Was this ZScore computed from a stored baseline ?
                if data.Results.ZScore.FlagBaselineFile == 0 % No - check new baseline limits
                    % Get exact time samples closest to these values
                    if round(1000*baseline(1)) < round(1000*ImageGridTime(1)) | round(1000*baseline(2)) > round(1000*ImageGridTime(end)) % Accuracy set to the ms
                        errordlg('Please set baseline time extremas within the available time window for this result file')
                        return
                    end

                    [tmp, tmin] = min(abs(ImageGridTime-baseline(1)));
                    [tmp, tmax] = min(abs(ImageGridTime-baseline(2)));

                    baseline_new = [tmin tmax];

                    test =  (sum(baseline_new == ...
                        data.Results.ZScore.Baseline) ~= 2); % New baseline requested - recompute Z-Score

                elseif test ~= 1 % Yes - Baseline was stored in a file - do not recompute the ZScore

                    [tmp, tmin] = min(abs(ImageGridTime-baseline(1)));
                    [tmp, tmax] = min(abs(ImageGridTime-baseline(2)));

                    test = 0;

                elseif varargin{1} == 1 % Force computation
                    [tmp, tmin] = min(abs(ImageGridTime-baseline(1)));
                    [tmp, tmax] = min(abs(ImageGridTime-baseline(2)));

                    test = 1;

                end

            else
                [tmp, tmin] = min(abs(ImageGridTime-baseline(1)));
                [tmp, tmax] = min(abs(ImageGridTime-baseline(2)));

                test = 1;
            end

            if 1%test

                set(ZScore_fig,'Pointer','watch')

                bst_message_window('wrap',...
                    {'Computing a new Z-Score map from this baseline. . .',...
                        sprintf('Baseline is %d time samples large (out of %d samples available)',...
                        tmax-tmin+1, length(ImageGridTime)),...
                    })

                % Store original results patter in memory
                if 1%(~isfield(data.Results,'ImageGridAmpOrig') & FlagBaselineFile == 0) | varargin{1} == 1
                    data.Results.ImageGridAmpOrig = data.Results.ImageGridAmp;
                end

                if FlagBaselineFile == 0 % remove the zero here !!!
                    data.Results.ImageGridAmp = data.Results.ImageGridAmpOrig;
                end

                %data.Results.ImageGridAmp = abs(data.Results.ImageGridAmp);

                % Compute ZScores on that baseline
                data.Results.ZScore.Baseline = [tmin tmax];
                data.Results.ZScore.BaselineTime = [baseline(1) baseline(2)];

                data.Results.ImageGridAmp = double(data.Results.ImageGridAmp);

                % Average activity and its standard deviation over the baseline, for each source
                data.Results.ZScore.ImageGridZ.mean = ...
                    mean((data.Results.ImageGridAmp(:,tmin:tmax)),2);
                data.Results.ZScore.ImageGridZ.std = ...
                    std((data.Results.ImageGridAmp(:,tmin:tmax)),0,2);
                data.Results.ZScore.ImageGridZ.std = ...
                    data.Results.ZScore.ImageGridZ.std + max(data.Results.ZScore.ImageGridZ.std)*eps; % Avoid devide-by-zero issues

                if ~isempty(get(hZScore_fig.CorrectionFactor,'String'))
                    data.Results.ZScore.ImageGridZ.std = str2num(get(hZScore_fig.CorrectionFactor,'String'))*data.Results.ZScore.ImageGridZ.std;
                end

                % Compute Z-Scores
                if FlagBaselineFile == 1 % Now load current results
                    ZScore = data.Results.ZScore;
                    data.Results = get(hZScore_fig.ResultFiles,'Userdata');
                    data.Results.ZScore = ZScore; clear ZScore
                    if ~isfield(data.Results,'ImageGridAmpOrig')
                        data.Results.ImageGridAmpOrig = [];
                    end
                    if isempty(data.Results.ImageGridAmpOrig)
                        data.Results.ImageGridAmpOrig = data.Results.ImageGridAmp; % Store original current density map
                    end
                    data.Results.ImageGridAmp = data.Results.ImageGridAmpOrig;
                    data.Results.ZScore.BaselineFile = get(hZScore_fig.LoadBaseline,'Userdata');
                    data.Results.ZScore.FlagBaselineFile = 1;
                else
                    data.Results.ZScore.FlagBaselineFile = 0;
                    data.Results.ZScore.BaselineFile = '';
                end

                %             bst_message_window('Updating Result File...'), drawnow
                %             %             data.Results.ImageGridAmp = single(data.Results.ImageGridAmpOrig);
                %             %             tmp = data.Results.ImageGridAmpOrig;
                %             %             data.Results.ImageGridAmpOrig = [];
                %             %             save_fieldnames(data.Results,get(hZScore_fig.Refresh,'Userdata'));
                %             ZScore = data.Results.ZScore;
                %             Users = get_user_directory;
                %             save(fullfile(Users.STUDIES,data.Results.FileName{1}),'ZScore','-append');
                %             clear ZScore
                %             %            data.Results.ImageGridAmpOrig = tmp; clear tmp
                %            data.Results.ImageGridAmp = double(data.Results.ImageGridAmp);

                data.Results.ImageGridAmp = (data.Results.ImageGridAmp)-repmat(data.Results.ZScore.ImageGridZ.mean,1,size(data.Results.ImageGridAmp,2));
                iStd = spdiags(1./data.Results.ZScore.ImageGridZ.std, 0, length(data.Results.ZScore.ImageGridZ.std), length(data.Results.ZScore.ImageGridZ.std)) ;

                %data.Results.ZScore.ImageGridZ.Amp = data.Results.ImageGridAmp;

                bst_message_window('Computing a new Z-Score map from this baseline...-> DONE'), drawnow

            else % Compute ZScore from saved mean and std maps

                data.Results.ImageGridAmpOrig = double(data.Results.ImageGridAmp);
                %data.Results.ImageGridAmp= data.Results.ZScore.ImageGridZ.Amp ;

                data.Results.ImageGridAmp = (data.Results.ImageGridAmp)-repmat(data.Results.ZScore.ImageGridZ.mean,1,size(data.Results.ImageGridAmp,2));
                iStd = spdiags(1./data.Results.ZScore.ImageGridZ.std, 0, length(data.Results.ZScore.ImageGridZ.std), length(data.Results.ZScore.ImageGridZ.std)) ;

            end

            %data.Results.ImageGridAmp = rand(size(iStd*data.Results.ImageGridAmp));
            data.Results.ImageGridAmp = single(abs(iStd*data.Results.ImageGridAmp));

            set(ZScore_fig,'Pointer','arrow')

        case 0 % Switch to absolute current density mapping

            bst_message_window('Previous Z-Score map on same baseline is used'), drawnow
            data.Results.ImageGridAmp = data.Results.ImageGridAmpOrig;
            set(hZScore_fig.ZScoreThresholdApply,'Value',0)
            dataplot_cb('ToggleButtonColor',hZScore_fig.ZScoreThresholdApply)
        end

        % Update display
        global hDataplot
        ctime = get(hDataplot.TimeSlider,'Value');
        tindx = findclosest(ctime,[data.Display.Time(1):data.Display.SamplingRate:data.Display.Time(end)]');

        if ishandle(data.Results.Display.Handles.hPatch) % surface handle for result display is valid
            data.Results.Display.OPTIONS.Time = tindx;
            [tmp,OPTIONS] = bst_imaging_display(data.Results.ImageGridAmp,data.Results.Display.OPTIONS); clear tmp
        end

        %----------------------------------------------------------------------------
    case 'ZScoreThreshold'
        %Thresholded ZScore map in terms of multiples of sigmas in the baseline

        Zscore_fig = findobj(0,'Tag','tesselation_select');
        hZScore_fig = guihandles(Zscore_fig);
        TessWin  = findobj(0,'tag','tessellation_window');
        %if isempty(TessWin), return, end
        %hTessWin = guihandles(TessWin);

        %TessWinUserData = get(TessWin,'Userdata');

        % Are we switching from ZScore to Absolute mapping or vice-versa ?
        data.Results = get(hZScore_fig.ResultFiles,'Userdata');

        switch get(hZScore_fig.ZScoreThresholdApply,'Value')
        case 1 % Switch to ZScore map and threshold it
            if isempty(data.Results), return, end % No results were loaded

            if ~isfield(data.Results,'ZScore') % No ZScore was defined beforehand
                dataplot_cb ZScore % Compute the ZScore map
                data.Results = get(hZScore_fig.ResultFiles,'Userdata');
            end

            set(hZScore_fig.ZScore,'Value',1);
            dataplot_cb('ToggleButtonColor',hZScore_fig.ZScore)

            ZThres = str2num(get(hZScore_fig.ZScoreThreshold,'String'));
            if isempty(ZThres)
                set(hZScore_fig.ZScoreThreshold,'String',2); % Apply Default
                ZThres = str2num(get(hZScore_fig.ZScoreThreshold,'String'));
            end

            data.Results.ZScore.ImageGridZ.ZThres = ZThres;

            set(hZScore_fig.ZScoreThreshold,'Userdata',ZThres / max(data.Results.ImageGridAmp(:)))

        case 0
            % Don't apply thresholded map
            set(hZScore_fig.ZScoreThreshold,'Userdata',0)

        end

        %data.Results.ImageGridAmp(iZeroed) = eps ; %Thresholded Map

        % Update display
        set(hZScore_fig.ResultFiles,'Userdata',data.Results);

        dataplot_cb ScaleColormap
        %dataplot_cb tesselation_select_done
        %----------------------------------------------------------------------------
    case 'IntegSTHistogram' % Integrate Saptio-temporal histogram overtime
        Zscore_fig = findobj(0,'Tag','tesselation_select');
        hZScore_fig = guihandles(Zscore_fig);

        avResults = get(hZScore_fig.ResultFiles,'Userdata');
        if ~isfield(avResults,'HistoSrc'), return, end % This is not an histogram result set

        % Integrate over time
        figure, imagesc(avResults.ImageGridTime, 1:size(avResults.ImageGridAmp,1),avResults.ImageGridAmp)
        colorbar
        avResults.ImageGridAmp = max(avResults.ImageGridAmp,[],2) *ones(1,size(avResults.ImageGridAmp,2));

        set(hZScore_fig.ResultFiles,'Userdata',avResults);

        %----------------------------------------------------------------------------

    case 'STHistogram' % Spatio-temporal histogram of activations (in terms of ZScore) on a series of result files

        Zscore_fig = findobj(0,'Tag','tesselation_select');
        hZScore_fig = guihandles(Zscore_fig);

        % Look for data.Results file in the current study folder
        % Get the path
        [path,file,ext] = fileparts(Users.CurrentData.StudyFile);
        % Dir list and look for result files
        cd(fullfile(Users.STUDIES,path))
        ResFiles = dir('*result*.mat');
        for k = 1:length(ResFiles)
            ResFiles(k).name = strrep(ResFiles(k).name,'.mat','');
        end
        [ResSelect,ok] = listdlg('Liststring',{ResFiles.name},'Selectionmode','multiple','Name','Please select one or sevral result files',...
            'listsize',[400 400]);

        if ok == 0
            return
        end

        hw = waitbar(0,['Computing the spatio-temporal histogram of activation maps through a set of ',int2str(length(ResSelect)),' result files...']);
        nDataFiles = 0;

        % Get Value for Z threshold
        ZThres = str2num(get(hZScore_fig.ZScoreThreshold,'String'));
        if isempty(ZThres)
            set(hZScore_fig.ZScoreThreshold,'String',2); % Apply Default
            ZThres = str2num(get(hZScore_fig.ZScoreThreshold,'String'));
        end

        for k = 1:length(ResSelect)
            dataplot_cb('LoadResultFile',ResFiles(ResSelect(k)).name)
            ResFiles(ResSelect(k)).name
            data.Results = get(hZScore_fig.ResultFiles,'Userdata');

            if k == 1
                avResults = data.Results;
                avResults.ImageGridAmp = 0* avResults.ImageGridAmp; % ImageGridAmp becomes a hit-count table -
                % cell (i,j) will indicate how many times source i had a ZScore > ZThres at time j

                avResults.HistoSrc = cell(length(ResSelect));
            end

            if isfield(data.Results,'ZScore')

                %             % Adapt length to the length of shortest result file
                %             cropMin = 0;
                %             cropMax = 0;
                %             if data.Results.Time(1) > avResults.Time(1)
                %                 cropMin = data.Results.Time(1)-avResults.Time(1);
                %                 avResults.ImageGridAmp = avResults.ImageGridAmp(:,cropMin+1:end);
                %                 avResults.ImageGridTime = avResults.ImageGridTime(:,cropMin+1:end);
                %                 avResults.Fsynth = avResults.Fsynth(:,cropMin+1:end);
                %                 avResults.ZScore.ImageGridZ.Amp = avResults.ZScore.ImageGridZ.Amp(:,cropMin+1:end);
                %
                %             end
                %             if data.Results.Time(1) < avResults.Time(1)
                %                 cropMin = -data.Results.Time(1)+avResults.Time(1);
                %                 %                 data.Results.ImageGridAmp = data.Results.ImageGridAmp(:,cropMin+1:end);
                %                 data.Results.ZScore.ImageGridZ.Amp = data.Results.ZScore.ImageGridZ.Amp(:,cropMin+1:end);
                %                 data.Results.Time = data.Results.Time(cropMin+1:end);
                %             end
                %             if data.Results.Time(end) > avResults.Time(end)
                %                 cropMax = data.Results.Time(end)-avResults.Time(end);
                %                 %                 data.Results.ImageGridAmp = data.Results.ImageGridAmp(:,1:end-cropMax);
                %                 data.Results.ZScore.ImageGridZ.Amp = data.Results.ZScore.ImageGridZ.Amp(:,1:end-cropMax);
                %                 data.Results.Time = data.Results.Time(1:end-cropMax);
                %             end
                %             if data.Results.Time(end) < avResults.Time(end)
                %                 cropMax = -data.Results.Time(end)+avResults.Time(end);
                %                 avResults.ImageGridAmp = avResults.ImageGridAmp(:,1:end-cropMax);
                %                 avResults.ImageGridTime = avResults.ImageGridTime(:,1:end-cropMax);
                %                 avResults.Fsynth = avResults.Fsynth(:,1:end-cropMax);
                %                 avResults.ZScore.ImageGridZ.Amp = avResults.ZScore.ImageGridZ.Amp(:,1:end-cropMax);
                %             end
                %
                avResults.Time = data.Results.Time;

                data.Results.ZScore.ImageGridZ.Amp = (data.Results.ImageGridAmp)-repmat(data.Results.ZScore.ImageGridZ.mean,1,size(data.Results.ImageGridAmp,2));
                iStd = spdiags(1./data.Results.ZScore.ImageGridZ.std, 0, length(data.Results.ZScore.ImageGridZ.std), length(data.Results.ZScore.ImageGridZ.std)) ;
                data.Results.ZScore.ImageGridZ.Amp = abs(iStd*data.Results.ZScore.ImageGridZ.Amp);

                avResults.HistoSrc{k} = sparse(data.Results.ZScore.ImageGridZ.Amp > ZThres);
                % avResults.HistoSrc{k} stores the source indices and time where sources had greater ZScore that ZThres
                avResults.ImageGridAmp = avResults.ImageGridAmp + avResults.HistoSrc{k} ;

                if k == 1
                    prev = figure;
                end

                nDataFiles = nDataFiles + 1;
                k
                figure(prev)
                plot(avResults.ImageGridTime, sum(avResults.HistoSrc{k},1))

            else % No ZScore available - skip file

                bst_message_window([ResFiles(ResSelect(k)).name,' does not contain any ZScore map - skip file'])

            end

            waitbar(k/length(ResSelect))

        end

        avResults.ImageGridAmp = 100*avResults.ImageGridAmp/nDataFiles; % Scale to 100% hits

        cd(Users.STUDIES);
        [path,file,ext] = fileparts(Users.CurrentData.StudyFile);
        cd(path)
        c = clock;
        ResFiletmp = strrep(ResFiles(ResSelect(end)).name,Users.STUDIES,'');
        ResFiletmp = strrep(ResFiletmp,'.mat','');
        Iunder = findstr(ResFiletmp,'_');
        DataFile= [ResFiletmp(1:Iunder(end)-1)];

        newname = ([DataFile sprintf('_HISTO_%02.0f%02.0f',c(4:5)) ext]);
        i = 0;
        while(exist(newname,'file')),
            i = i+1; % subtract another minute
            c(5) = mod(c(5) - 1,60);
            newname = ([DataFile sprintf('_HISTO_%02.0f%02.0f',MethodeCode,c(4),c(5)) ext]);
        end

        save_fieldnames(avResults,newname)

        delete(hw)

        msgbox(['Average activation maps saved in ', newname])
        set(hZScore_fig.ResultFiles,'Userdata',avResults)

        % Refresh results list
        dataplot_cb mesh_rendering
        ResFiles = get(hZScore_fig.ResultFiles,'String');
        iFile = find(strcmp(ResFiles,strrep(newname,'.mat','')));
        set(hZScore_fig.ResultFiles,'Value',iFile)

        %----------------------------------------------------------------------------

    case 'AverageMaps' % Average activations across multiple result files

        Zscore_fig = findobj(0,'Tag','tesselation_select');
        hZScore_fig = guihandles(Zscore_fig);

        % Look for data.Results file in the current study folder
        % Get the path
        [path,file,ext] = fileparts(Users.CurrentData.StudyFile);
        % Dir list and look for result files
        cd(fullfile(Users.STUDIES,path))
        ResFiles = dir('*result*.mat');
        for k = 1:length(ResFiles)
            ResFiles(k).name = strrep(ResFiles(k).name,'.mat','');
        end
        [ResSelect,ok] = listdlg('Liststring',{ResFiles.name},'Selectionmode','multiple','Name','Please select one or sevral result files',...
            'listsize',[400 400]);

        if ok == 0
            return
        end

        % Now compute ZScore for each of the selected files
        hw = waitbar(0,['Averaging activation maps through a set of ',int2str(length(ResSelect)),' result files...']);
        for k = 1:length(ResSelect)
            dataplot_cb('LoadResultFile',ResFiles(ResSelect(k)).name)
            ResFiles(ResSelect(k)).name
            data.Results = get(hZScore_fig.ResultFiles,'Userdata');

            if k == 1
                avResults = data.Results;
                avResults.ImageGridAmp = avResults.ImageGridAmp/length(ResSelect);
            end

            % !!!!!!Suppose all result files were computed on the same time window - discard what's commented below

            %         % Adapt length to the length of shortes result file
            %         cropMin = 0;
            %         cropMax = 0;
            %         if data.Results.Time(1) > avResults.Time(1)
            %             cropMin = data.Results.Time(1)-avResults.Time(1);
            %             avResults.ImageGridAmp = avResults.ImageGridAmp(:,cropMin+1:end);
            %             avResults.ImageGridTime = avResults.ImageGridTime(:,cropMin+1:end);
            %             avResults.Time = avResults.Time(:,cropMin+1:end);
            %             avResults.Fsynth = avResults.Fsynth(:,cropMin+1:end);
            %
            %             if isfield(data.Results,'ZScore')
            %                 avResults.ZScore.ImageGridZ.Amp = avResults.ZScore.ImageGridZ.Amp(:,cropMin+1:end);
            %             end
            %         end
            %         if data.Results.Time(1) < avResults.Time(1)
            %             cropMin = -data.Results.Time(1)+avResults.Time(1);
            %             data.Results.ImageGridAmp = data.Results.ImageGridAmp(:,cropMin+1:end);
            %             ImageGridTime = ImageGridTime(:,cropMin+1:end);
            %             data.Results.Time = data.Results.Time(:,cropMin+1:end);
            %             if isfield(data.Results,'ZScore')
            %                 data.Results.ZScore.ImageGridZ.Amp = data.Results.ZScore.ImageGridZ.Amp(:,cropMin+1:end);
            %             end
            %         end
            %         if data.Results.Time(end) > avResults.Time(end)
            %             cropMax = data.Results.Time(end)-avResults.Time(end);
            %             data.Results.ImageGridAmp = data.Results.ImageGridAmp(:,1:end-cropMax);
            %             ImageGridTime = ImageGridTime(:,1:end-cropMax);
            %             data.Results.Time = data.Results.Time(:,1:end-cropMax);
            %             if isfield(data.Results,'ZScore')
            %                 data.Results.ZScore.ImageGridZ.Amp = data.Results.ZScore.ImageGridZ.Amp(:,1:end-cropMax);
            %             end
            %         end
            %         if data.Results.Time(end) < avResults.Time(end)
            %             cropMax = -data.Results.Time(end)+avResults.Time(end);
            %             avResults.ImageGridAmp = avResults.ImageGridAmp(:,1:end-cropMax);
            %             avResults.ImageGridTime = avResults.ImageGridTime(:,1:end-cropMax);
            %             avResults.Fsynth = avResults.Fsynth(:,1:end-cropMax);
            %             avResults.Time = avResults.Time(:,1:end-cropMax);
            %             if isfield(data.Results,'ZScore')
            %                 avResults.ZScore.ImageGridZ.Amp = avResults.ZScore.ImageGridZ.Amp(:,1:end-cropMax);
            %             end
            %         end

            if k == 1
                prev = figure;
                if isfield(data.Results,'ZScore')
                    data.Results.ZScore.ImageGridZ.Amp = (data.Results.ImageGridAmp)-repmat(data.Results.ZScore.ImageGridZ.mean,1,size(data.Results.ImageGridAmp,2));
                    iStd = spdiags(1./data.Results.ZScore.ImageGridZ.std, 0, length(data.Results.ZScore.ImageGridZ.std), length(data.Results.ZScore.ImageGridZ.std)) ;
                    data.Results.ZScore.ImageGridZ.Amp = abs(iStd*data.Results.ZScore.ImageGridZ.Amp);
                end
            end

            if k >1
                avResults.ImageGridAmp = avResults.ImageGridAmp + (data.Results.ImageGridAmp)/length(ResSelect);
                if isfield(data.Results,'ZScore')
                    avResults.ZScore.ImageGridZ.Amp = avResults.ZScore.ImageGridZ.Amp + (data.Results.ZScore.ImageGridZ.Amp)/length(ResSelect);
                end
            end
            figure(prev)
            plot(ImageGridTime,mean(abs(data.Results.ImageGridAmp),1),'b'), hold on
            plot(avResults.ImageGridTime,mean(abs(avResults.ImageGridAmp),1),'r'), hold off

            waitbar(k/length(ResSelect))

        end

        %avResults.ImageGridAmp = avResults.ImageGridAmp/length(ResSelect);
        %     if isfield(avResults,'ZScore')
        %         avResults.ZScore.ImageGridZ.Amp = avResults.ZScore.ImageGridZ.Amp/length(ResSelect);
        %     end

        cd(Users.STUDIES);
        [path,file,ext] = fileparts(Users.CurrentData.StudyFile);
        cd(path)
        c = clock;
        ResFiletmp = strrep(ResFiles(ResSelect(end)).name,Users.STUDIES,'');
        ResFiletmp = strrep(ResFiletmp,'.mat','');
        Iunder = findstr(ResFiletmp,'_');
        DataFile= [ResFiletmp(1:Iunder(end-1)-1)];

        newname = ([DataFile sprintf('_AVRresults_%02.0f%02.0f',c(4:5)) ext]);
        i = 0;
        while(exist(newname,'file')),
            i = i+1; % subtract another minute
            c(5) = mod(c(5) - 1,60);
            newname = ([DataFile sprintf('_AVRresults_%02.0f%02.0f',MethodeCode,c(4),c(5)) ext]);
        end

        save_fieldnames(avResults,newname)

        delete(hw)

        msgbox(['Average activation maps saved in ', newname])
        set(hZScore_fig.ResultFiles,'Userdata',avResults)

        % Refresh results list
        dataplot_cb mesh_rendering
        ResFiles = get(hZScore_fig.ResultFiles,'String');
        iFile = find(strcmp(ResFiles,strrep(newname,'.mat','')));
        set(hZScore_fig.ResultFiles,'Value',iFile)

        %----------------------------------------------------------------------------
    case 'ZScoreBatch' % Compute ZScore on a set of Result files

        % Look for data.Results file in the current study folder
        % Get the path
        [path,file,ext] = fileparts(Users.CurrentData.StudyFile);
        % Dir list and look for result files
        cd(fullfile(Users.STUDIES,path))
        ResFiles = dir('*result*.mat');
        for k = 1:length(ResFiles)
            ResFiles(k).name = strrep(ResFiles(k).name,'.mat','');
        end
        [ResSelect,ok] = listdlg('Liststring',{ResFiles.name},'Selectionmode','multiple','Name','Please select one or sevral result files',...
            'listsize',[400 400]);

        if ok == 0
            return
        end

        % Now compute ZScore for each of the selected files
        hw = waitbar(0,['Computing ZScore for the selected set of ',int2str(length(ResSelect)),' result files...']);
        for k = 1:length(ResSelect)
            dataplot_cb('LoadResultFile',ResFiles(ResSelect(k)).name)
            ResFiles(ResSelect(k)).name
            Zscore_fig = findobj(0,'Tag','tesselation_select');
            hZScore_fig = guihandles(Zscore_fig);
            set(hZScore_fig.ZScore,'Value',1)
            dataplot_cb('ToggleButtonColor',hZScore_fig.ZScore)
            dataplot_cb('ZScore',1) % Force the computation of the ZScore for each of the file (ie force test to 1 in dataplot_cb ZScore)
            waitbar(k/length(ResSelect))
        end
        set(hZScore_fig.Baseline,'Userdata',0)

        delete(hw)

        %----------------------------------------------------------------------------

    case 'AnalyzeCorticalMap' % Spatio-temporal segmentation of the cortical current density maps

        Zscore_fig = findobj(0,'Tag','tesselation_select');
        hTessSelect = guihandles(Zscore_fig);
        available = findobj(Zscore_fig,'Tag','available');
        removed =  findobj(Zscore_fig,'Tag','removed');
        availableID = get(removed,'Value');
        IDs = get(removed,'String');
        if isempty(IDs), return, end

        % Identification of the active mesh surfaces
        Visu = get(DATAPLOT,'Userdata');

        try
            load(Visu.Tesselation,'Comment')
        catch
            cd(Users.SUBJECTS)
            load(Visu.Tesselation,'Comment')
        end

        if ~iscell(IDs)
            imesh = find(strcmp(IDs,Comment));
        end

        nsurf = imesh(availableID); % Selected surface

        TessWin  = findobj(0,'tag','tessellation_window');
        hTessWin = guihandles(TessWin);
        TessWinUserData = get(TessWin,'Userdata');

        if ~isfield(TessWinUserData,'VertConn') % Vertex connectivity not available here - load or compute
            nsurfaces = length(Comment);

            % ----------------- Load current surface parameters -
            load(Visu.Tesselation,'Faces','Vertices')
            FV.faces = Faces{nsurf}; clear Faces
            FV.vertices = Vertices{nsurf}'; clear Vertices
            % ----------------- Load current surface parameters - DONE

            load(fullfile(Users.SUBJECTS,Users.CurrentData.SubjectFile),'VertConn'); % Is the vertex connectivity available ?

            if isempty(VertConn) % Vertex connectivity was not defined before
                VertConn = cell(nsurfaces,1);
            else
                if exist(VertConn,'file')
                    load(VertConn)
                else
                    VertConn = cell(nsurfaces,1);
                end
            end

            if isempty(VertConn{nsurf}) % Compute the vertex connectivity
                [pathname, filename] = fileparts(Visu.Tesselation);
                VertConn{nsurf} = vertices_connectivity(FV);
                clear FV
                VertConnFile = [filename,'_vertconn.mat'];

                % Save this file in the current subject folder
                [path,name,ext] = fileparts(fullfile(Users.SUBJECTS,Users.CurrentData.SubjectFile));
                if ~isempty(path)
                    cd(path)
                else
                    cd(Users.SUBJECTS)
                end

                save(VertConnFile,'VertConn');

                VertConn = VertConnFile;
                save(Users.CurrentData.SubjectFile,'VertConn','-append')
            else % Load VC from the file
                load(Users.CurrentData.SubjectFile,'VertConn');
            end
            load(VertConn);
            VertConn = VertConn{nsurf};
            TessWinUserData.VertConn = VertConn; clear VertConn;
        end

        % Get Z threshold
        ZThres = str2num(get(hTessSelect.ZScoreThreshold,'String'))  ;

        % Vertices which maximum Z score is above ZThres
        data.Results = get(hTessSelect.ResultFiles,'Userdata');
        ImageGridAmp =  data.Results.ImageGridAmp; clear data.Results

        % Get baseline info
        data.Results = get(hTessSelect.ResultFiles,'userdata');
        AnalyzeWin = [30 90]/1000; % in sec
        [tmp, tmin] = min(abs(ImageGridTime-AnalyzeWin(1)));
        [tmp, tmax] = min(abs(ImageGridTime-AnalyzeWin(2)));
        clear data.Results

        %Baseline = data.Results.ZScore.Baseline; clear data.Results
        %Time = Baseline(end)+1:size(ImageGridAmp,2); %setdiff([1:size(ImageGridAmp,2)],[Baseline(1):Baseline(end)]);

        Time = [tmin:tmax];

        MaxVerts = max(ImageGridAmp(:,Time)');
        iThres = find(MaxVerts > ZThres);
        % Threshold in terms of correlation
        CorrThres = .9;

        Scout{1} = [];
        iScout = 1;
        imax = 1;

        % Sort MaxVerts
        [s, iMaxVertsSort] = sort(MaxVerts); clear s

        CONNEX = 1; % if 1; enforce connexity of each clusters

        while ~isempty(imax)
            %%%% Segmentation process
            if iScout == 1 % We are starting from scratch
                % Initialization
                % Start region growing from the vertex with maximum ZScore
                imax = iMaxVertsSort(end);
            else
                tmp = setdiff(iMaxVertsSort,[Scout{:}]);
                [mm, imax] = max(MaxVerts(tmp)); clear mm
                imax = tmp(imax);
            end

            % Does imax checks iThres
            if isempty(find(iThres == imax))
                break
            end

            if CONNEX == 1
                % Look for his neighbors
                voiz = patch_swell(imax,TessWinUserData.VertConn);
                % Select those that check the Zthres condition
                voiz = intersect(iThres,voiz);

            else % Just look for all the points above ZThres
                voiz = setdiff(iThres,imax);
            end

            % Remove those previously selected
            voiz = setdiff(voiz,[Scout{:}]);

            % Compute correlation coefficients of their time series
            %CorrCoeff = abs(corrcoef(ImageGridAmp([imax,voiz],Time)'));
            if ~isempty(voiz)
                CorrCoeff = abs(ImageGridAmp(imax,Time)/norm(ImageGridAmp(imax,Time)))...
                    *abs(ImageGridAmp(voiz,Time)'*inorcol(ImageGridAmp(voiz,Time)'));
            else
                CorrCoeff = CorrThres/10;
            end

            % Apply Correlation Threshold
            %icorr = find(CorrCoeff(1,2:end)>CorrThres);
            icorr = find(CorrCoeff>CorrThres);

            Scout{iScout} = [imax,voiz(icorr)];

            while ~isempty(icorr)

                if CONNEX == 1
                    voiz = patch_swell(voiz(icorr),TessWinUserData.VertConn);
                    % Select those that check the Zthres condition
                    voiz = intersect(iThres,voiz);
                else
                    voiz = setdiff(iThres,imax);
                end

                % Remove those previously selected
                voiz = setdiff(voiz,[Scout{:}]);

                % Compute correlation coefficients of their time series
                %CorrCoeff = abs(corrcoef(ImageGridAmp([imax,voiz],Time)'));
                if ~isempty(voiz)
                    CorrCoeff = abs(ImageGridAmp(imax,Time)/norm(ImageGridAmp(imax,Time)))...
                        *abs(ImageGridAmp(voiz,Time)'*inorcol(ImageGridAmp(voiz,Time)'));

                    CorrCoeffBaseLine = abs(ImageGridAmp(imax,125:500)/norm(ImageGridAmp(imax,125:500)))...
                        *abs(ImageGridAmp(voiz,125:500)'*inorcol(ImageGridAmp(voiz,125:500)'));
                else

                    CorrCoeff = CorrThres/10;
                end

                % Apply Correlation Threshold
                %icorr = find(CorrCoeff(1,2:end)>CorrThres);
                icorr = find(CorrCoeff>CorrThres);

                Scout{iScout} = [Scout{iScout},voiz(icorr)];

            end

            if length([Scout{:}]) ~= length(unique([Scout{:}]))
                keyboard
            end

            length(iThres) - length([Scout{:}])

            CorticalScouts.CorticalProbePatches(iScout) = Scout(iScout);
            CorticalScouts.CorticalMarkersLabels{iScout} = int2str(iScout);
            CorticalScouts.CorticalSpots(iScout) = imax;
            if iScout == 1
                load(Visu.Tesselation,'Vertices')
            end
            CorticalScouts.CorticalMarkers(iScout,:) = Vertices{nsurf}(:,Scout{iScout}(1));

            iScout = iScout+1;
            bst_message_window('')

        end

        newScout = 1;

        for iScout = 1:length(Scout)
            if length(Scout{iScout}) == 1 % Single point scout: do we keep it ?
                % Check if its neighbors with activation < ZThres have time series that correlate with his with corrcoeff > CorrThres
                voiz = patch_swell(Scout{iScout},TessWinUserData.VertConn);

                if ~isempty(voiz)
                    CorrCoeff = abs(ImageGridAmp(imax,Time)/norm(ImageGridAmp(imax,Time)))...
                        *abs(ImageGridAmp(voiz,Time)'*inorcol(ImageGridAmp(voiz,Time)'));
                else
                    CorrCoeff = CorrThres/10;
                end

                icorr = find(CorrCoeff>CorrThres);

                if isempty(icorr)
                    bst_message_window('Remove spurious Scout')
                    %               Scout{iScout} = [];
                    %               iScout = iScout - 1;
                else
                    newCorticalScouts.CorticalProbePatches(newScout) = Scout(iScout);
                    newCorticalScouts.CorticalMarkersLabels{newScout} = int2str(newScout);
                    newCorticalScouts.CorticalSpots(newScout) = Scout{iScout}(1);
                    if newScout == 1
                        load(Visu.Tesselation,'Vertices')
                    end
                    newCorticalScouts.CorticalMarkers(newScout,:) = Vertices{nsurf}(:,Scout{iScout}(1));
                    newScout = newScout + 1;
                end

            else
                newCorticalScouts.CorticalProbePatches(newScout) = Scout(iScout);
                newCorticalScouts.CorticalMarkersLabels{newScout} = int2str(newScout);
                newCorticalScouts.CorticalSpots(newScout) = Scout{iScout}(1);
                if newScout == 1
                    load(Visu.Tesselation,'Vertices')
                end
                newCorticalScouts.CorticalMarkers(newScout,:) = Vertices{nsurf}(:,Scout{iScout}(1));
                newScout = newScout + 1;
            end

            clear CorticalScouts
            CorticalScouts = newCorticalScouts;

        end

        CorticalScouts.VertConn = TessWinUserData.VertConn;
        CorticalScouts.CorticalProbeDepth = ones(1,length(Scout));

        cd(Users.STUDIES)

        length(CorticalScouts.CorticalProbePatches)

        save test CorticalScouts
        %----------------------------------------------------------------------------

    end

Produced by color_mat2html, a customized BrainStorm 2.0 (Alpha) version of mat2html on Tue Nov 2 16:17:59 2004
Cross-Directory links are: ON