text:
function info = roi_contour_map(dat, varargin)
% Draw a pattern map of one slice (either saggital, axial, or coronal) that
% shows the most voxels, or the slice that you specify (e.g., x = #).
% You can also draw outlines for the significant voxels from a statistical test.
%
% :Usage:
% ::
%
%    info = roi_contour_map(dat, varargin)
%
% ..
%     Author and copyright information:
%
%     Copyright (C) 2014 Choong-Wan (Wani) Woo
%
%     This program is free software: you can redistribute it and/or modify
%     it under the terms of the GNU General Public License as published by
%     the Free Software Foundation, either version 3 of the License, or
%     (at your option) any later version.
%
%     This program is distributed in the hope that it will be useful,
%     but WITHOUT ANY WARRANTY; without even the implied warranty of
%     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
%     GNU General Public License for more details.
%
%     You should have received a copy of the GNU General Public License
%     along with this program.  If not, see <http://www.gnu.org/licenses/>.
% ..
%
% :Inputs:
%
%   **dat:**
%        dat can be fmri_data, statistic_image (to mark significant
%        voxels), and region objects. If your dat is the "region"
%        object, please add 'cluster' as an optional input.
%        You can display two pattern maps for the purpose of
%        comparison by putting additional columns of data in cell array.
%
%        An example of displaying two pattern maps:
%        ::
%
%            dat{1} = region('img1.nii');
%            dat{2} = region('img2.nii'); )
%
% :Optional Inputs:
%
%   **'cluster':**
%        When the data is a region object, you need this option.
%
%   **'sig':**
%        This option outlines significant voxels. To use this
%        option, data in a format of statistic_image with a "sig"
%        field should be given.
%
%   **'colorbar':**
%        display colorbar under the plot.
%
%   **'use_same_range':**
%        When you display two pattern maps, this option uses the
%        same color range for the two maps.
%
%   **'surf':**
%        surface plot rather than voxel-by-voxel mapping.
%
%   **'xyz':**
%        When you want a specific view and slice, you can use this
%        option with 'coord'. (1:x - saggital view, 2:y - coronal
%        view, 3:z - axial view)
%
%   **'coord':**
%        With 'xyz' option, this specifies the slice displayed.
%
%   **'notfill':**
%        Default is to fill in the blank voxels using a black
%        color. With this option, you can color the blank voxels
%        with the white color.
%
%   **'whole':**
%        Default is dividing the data into contiguous regions and
%        show only one region that has the most voxels. This option
%        akes this function not to divide into contiguous regions.
%
%   **'colors' or 'color':**
%        you can specify your own colormap.
%
%   **'contour':**
%        Not fully implemented yet.
%
% :Outputs:
%
%   **info:**
%        information about the display with the following fields.
%
%   info.dat:
%     - [2x30 double] (xyz mesh)
%     - Z: [1x30 double] (z values)
%     - xyz: 3 (1:x-saggital, 2:y-coronal, 3:z-axial)
%     - xyz_coord: -2 (slice coordinate; in this case, z = 3)
%     - region_idx: 1
%
%
% :Examples: you can also see the same example and output in
% http://wagerlab.colorado.edu/wiki/doku.php/help/core/figure_gallery
% ::
%
%    mask{1} = 'dACC_hw_pattern_sl6mm.nii';
%    mask{2} = 'dACC_rf_pattern_sl6mm.nii';
%    for i = 1:2, cl{i} = region(mask{i}); end
%    info = roi_contour_map([cl{1} cl{2}], 'cluster', 'use_same_range', 'colorbar');
%
% ..
%    Copyright (C) 2014  Wani Woo
% ..

docolor = 0; % parsing optional inputs
dosig = 0;
docolorbar = 0;
usesamerange = 0;
usecluster = 0;
dosurf = 0;
docontour = 0; % needs to be implemented.
doxyz = 0;
docoord = 0;
donotfill = 0;
dowhole = 0;
outline_color_pos = 'r';
outline_color_neg = 'b';

for i = 1:length(varargin)
    if ischar(varargin{i})
        switch varargin{i}
            % functional commands
            case {'colors', 'color'}
                docolor = 1; colors = varargin{i+1};
            case {'sig'}
                dosig = 1;
            case {'sigcolor'}
                outline_color_pos = varargin{i+1};
                outline_color_neg = varargin{i+2};
            case {'colorbar'}
                docolorbar = 1;
            case {'use_same_range'}
                usesamerange = 1;
            case {'cluster'}
                usecluster = 1;
            case {'surf'}
                dosurf = 1;
            case {'contour'}
                docontour = 1;
            case {'xyz'}
                doxyz = 1;
                xyz_idx = varargin{i+1};
            case {'coord'}
                docoord = 1;
                coord_idx = varargin{i+1};
            case {'notfill'}
                donotfill = 1;
            case {'whole'}
                dowhole = 1;
        end
    end
end

% data inputs
if ~usecluster
    if ~iscell(dat)
        dat_temp{1} = dat;
        dat = dat_temp;
    end
end

rnum  = numel(dat);

for jj = 1:rnum
    if ~usecluster
        if dosig
            if any(strcmp(fields(dat{jj}), 'sig'))
                datsig{jj} = dat{jj};
                datsig{jj}.dat = dat{jj}.sig;
                datsig{jj}.dat = dat{jj}.dat.*~dat{jj}.sig*10^-6 + dat{jj}.sig;
                dat{jj}.sig = ones(size(dat{jj}.dat)); datsig{jj}.sig = ones(size(datsig{jj}.dat));
            else
                error('There is no sig data. Please check your data.');
            end
        else
            if any(strcmp(fields(dat{jj}), 'sig'))
                dat{jj}.sig = ones(size(dat{jj}.dat));
            end
        end

        ri = region(dat{jj});
        if dosig, risig = region(datsig{jj}); end
        if dowhole, ri = combine_region(ri); risig = combine_region(risig); end
    else
        ri = dat(jj);
    end
    clear max_vox idx idxidx;

    for i = 1:numel(ri)

        if ~docoord
            mm_c{i} = center_of_mass(ri(i).XYZmm, ones(1,size(ri(i).XYZmm, 2)));
        elseif docoord && ~doxyz
            warning('''coord'' option should be used with ''xyz''. This will use center of mass.');
            mm_c{i} = center_of_mass(ri(i).XYZmm, ones(1,size(ri(i).XYZmm, 2)));
        else
            mm_c{i} = center_of_mass(ri(i).XYZmm, ones(1,size(ri(i).XYZmm, 2)));
            [~, idx_coord_idx] = min(abs(ri(i).XYZmm(xyz_idx,:)-coord_idx));
            new_coord_idx = ri(i).XYZmm(xyz_idx,idx_coord_idx);
            mm_c{i}(xyz_idx) = new_coord_idx;
        end

        if ~doxyz
            [max_vox(i), idx(i)] = max([sum(ri(i).XYZmm(1,:) == mm_c{i}(1)), sum(ri(i).XYZmm(2,:) == mm_c{i}(2)), sum(ri(i).XYZmm(3,:) == mm_c{i}(3))]);
        else
            voxs = [sum(ri(i).XYZmm(1,:) == mm_c{i}(1)), sum(ri(i).XYZmm(2,:) == mm_c{i}(2)), sum(ri(i).XYZmm(3,:) == mm_c{i}(3))];
            idx(i) = xyz_idx;
            max_vox(i) = voxs(xyz_idx);
        end

    end

    [~, idxidx] = max(max_vox);

    xyz = 1:3;
    xyz(idx(idxidx)) = [];

    coord{jj}.dat = ri(idxidx).XYZmm(xyz,ri(idxidx).XYZmm(idx(idxidx),:) == mm_c{idxidx}(idx(idxidx)));
    coord{jj}.Z = ri(idxidx).Z(ri(idxidx).XYZmm(idx(idxidx),:) == mm_c{idxidx}(idx(idxidx)));
    if dosig
        coord{jj}.sig = risig(idxidx).Z(risig(idxidx).XYZmm(idx(idxidx),:) == mm_c{idxidx}(idx(idxidx)));
    end
    coord{jj}.xyz = idx(idxidx);
    coord{jj}.xyz_coord = mm_c{idxidx}(idx(idxidx));
    coord{jj}.region_idx = idxidx;

end

info = coord;

%% basic settingfor drawing

close all;
clear xyz;
figure;
if docolorbar
    set(gcf, 'color', 'w', 'position', [ 6   221   500*rnum   485]); % with colorbar
else
    set(gcf, 'color', 'w', 'position', [ 6   322   500*rnum   384]);
end

col = [0.0461    0.3833    0.5912
    0.2461    0.5833    0.7912
    0.4000    0.7608    0.6471
    0.6706    0.8667    0.6431
    0.9020    0.9608    0.5961
    1.0000    1.0000    0.7490
    0.9961    0.8784    0.5451
    0.9922    0.6824    0.3804
    0.9569    0.4275    0.2627
    0.8353    0.2431    0.3098
    0.6196    0.0039    0.2588];

colormap(col);
if docolor, colormap(colors); end

for jj = 1:rnum

    if dosig
        xyz = [coord{jj}.dat(1,:)', coord{jj}.dat(2,:)', coord{jj}.Z', sign(coord{jj}.Z') .* coord{jj}.sig'];
    else
        xyz = [coord{jj}.dat(1,:)', coord{jj}.dat(2,:)', coord{jj}.Z'];
    end

    [vX, vY] = meshgrid(unique(xyz(:,1)), unique(xyz(:,2)));

    vZ{jj} = NaN(size(vX));

    for i = 1:size(vX,2)
        for j = 1:size(vY,1)
            for k = 1:length(xyz)
                if isequal([vX(1,i) vY(j,1)], xyz(k,1:2))
                    vZ{jj}(j,i) = xyz(k,3);
                    if dosig
                        vsig{jj}(j,i) = xyz(k,4);
                    end
                end
            end
        end
    end

    vZ{jj} = flipud(vZ{jj});
    if dosig
        vsig{jj} = flipud(vsig{jj});
    end
    if docontour
        vZ{jj} = interp2(vZ{jj},5);
    end

    vox_size = vX(1,2)-vX(1,1);
    newX = vX/vox_size-vX(1,1)/vox_size+1;
    newY = vY/vox_size-vY(1,1)/vox_size+1;

    if docontour
        newX = repmat(1:size(vZ{jj},2), size(vZ{jj},1),1);
        newY = repmat((1:size(vZ{jj},1))', 1, size(vZ{jj},2));
    end

    %newX = (vX+2)/3-(vX(1,1)+2)/3+1;
    %newY = (vY+2)/3-(vY(1,1)+2)/3+1;

    xx{jj} = [newX(isnan(vZ{jj}))-.5 newX(isnan(vZ{jj}))-.5 newX(isnan(vZ{jj}))+.5 newX(isnan(vZ{jj}))+.5];
    yy{jj} = [newY(isnan(vZ{jj}))-.5 newY(isnan(vZ{jj}))+.5 newY(isnan(vZ{jj}))+.5 newY(isnan(vZ{jj}))-.5];

    if donotfill
        xx_outline{jj} = [newX(~isnan(vZ{jj}))-.5 newX(~isnan(vZ{jj}))-.5 newX(~isnan(vZ{jj}))+.5 newX(~isnan(vZ{jj}))+.5 newX(~isnan(vZ{jj}))-.5];
        yy_outline{jj} = [newY(~isnan(vZ{jj}))-.5 newY(~isnan(vZ{jj}))+.5 newY(~isnan(vZ{jj}))+.5 newY(~isnan(vZ{jj}))-.5 newY(~isnan(vZ{jj}))-.5];
    end

    if dosig
        x_out_pos{jj} = [newX(vsig{jj}==1)-.5 newX(vsig{jj}==1)-.5 newX(vsig{jj}==1)+.5 newX(vsig{jj}==1)+.5 newX(vsig{jj}==1)-.5];
        y_out_pos{jj} = [newY(vsig{jj}==1)-.5 newY(vsig{jj}==1)+.5 newY(vsig{jj}==1)+.5 newY(vsig{jj}==1)-.5 newY(vsig{jj}==1)-.5];
        x_out_neg{jj} = [newX(vsig{jj}==-1)-.5 newX(vsig{jj}==-1)-.5 newX(vsig{jj}==-1)+.5 newX(vsig{jj}==-1)+.5 newX(vsig{jj}==-1)-.5];
        y_out_neg{jj} = [newY(vsig{jj}==-1)-.5 newY(vsig{jj}==-1)+.5 newY(vsig{jj}==-1)+.5 newY(vsig{jj}==-1)-.5 newY(vsig{jj}==-1)-.5];
    end

    lim_max(jj) = max(max(vZ{jj}));
    lim_min(jj) = min(min(vZ{jj}));
end

if usesamerange
    lim_max = max(lim_max); lim_min = min(lim_min);
end

for jj = 1:rnum

    subplot(1,rnum,jj);

    if dosurf % not work that well yet, but kind of work
        % surf(newX, newY, vZ{jj}, 'EdgeColor','k');
        surf(interp2(vZ{jj}, 5), 'EdgeColor','none');
        set(gca, 'view', [-21.5 70])
        % axis off;
        % camlight left; lighting phong;
        set(gca, 'linewidth', 1, 'xtick', [], 'ytick', [], 'ztick', [], 'fontsize', 15);
        if usesamerange
            set(gca, 'zlim', [lim_min lim_max]);
        end
    else
        if ~usesamerange
            if docontour
                imagesc(vZ{jj}, [lim_min(jj) lim_max(jj)]);
            else
                imagesc(vZ{jj}, [lim_min(jj) lim_max(jj)]);
                if ~donotfill
                    set(gca, 'linewidth', 5, 'xtick', [], 'ytick', []);
                else
                    set(gca, 'linewidth', 0.0001, 'xtick', [], 'ytick', [], 'xcolor', 'w', 'ycolor', 'w', 'zcolor', 'w');
                end
            end
        else
            if docontour
                imagesc(vZ{jj}, [lim_min lim_max]);
                set(gca, 'linewidth', 0.0001, 'xtick', [], 'ytick', [], 'xcolor', 'w', 'ycolor', 'w', 'zcolor', 'w');
            else
                imagesc(vZ{jj}, [lim_min lim_max]);
                if ~donotfill
                    set(gca, 'linewidth', 5, 'xtick', [], 'ytick', []);
                else
                    set(gca, 'linewidth', 0.0001, 'xtick', [], 'ytick', [], 'xcolor', 'w', 'ycolor', 'w', 'zcolor', 'w');
                end
            end
        end

        if ~docontour
            y = size(vZ{jj},1); x = size(vZ{jj},2);
            if ~donotfill
                %for i = 1:x, line([i+.5,i+.5], [0 y+.5], 'color', [.2 .2 .2], 'linewidth', 1.5, 'linestyle', '-'); end
                %for i = 1:y, line([0 x+.5], [i+.5,i+.5], 'color', [.2 .2 .2], 'linewidth', 1.5, 'linestyle', '-'); end

                % black lines
                for i = 1:x, line([i+.5,i+.5], [0 y+.5], 'color', [.2 .2 .2], 'linewidth', 1, 'linestyle', '-'); end
                for i = 1:y, line([0 x+.5], [i+.5,i+.5], 'color', [.2 .2 .2], 'linewidth', 1, 'linestyle', '-'); end

                % white lines
%                 for i = 1:x, line([i+.5,i+.5], [0 y+.5], 'color', 'w', 'linewidth', 1, 'linestyle', '-'); end
%                 for i = 1:y, line([0 x+.5], [i+.5,i+.5], 'color', 'w', 'linewidth', 1, 'linestyle', '-'); end
            else
                for i = 1:x, line([i+.5,i+.5], [0 y+.5], 'color', repmat(.2, 1, 3), 'linewidth', 1.5, 'linestyle', '-'); end
                for i = 1:y, line([0 x+.5], [i+.5,i+.5], 'color', repmat(.2, 1, 3), 'linewidth', 1.5, 'linestyle', '-'); end
            end
        end

        hold on;
        if docontour
            for i = 1:size(xx{jj},1), fill(xx{jj}(i,:), yy{jj}(i,:), 'w'); end % this is wrong.
        else
            if ~donotfill
                % fill black
                for i = 1:size(xx{jj},1), fill(xx{jj}(i,:), yy{jj}(i,:), 'k'); end

                % fill white
%                 for i = 1:size(xx{jj},1), fill(xx{jj}(i,:), yy{jj}(i,:), 'w', 'edgecolor', 'w'); end
            else
                for i = 1:size(xx{jj},1), fill(xx{jj}(i,:), yy{jj}(i,:), 'w', 'edgecolor', 'w'); end
                xy_outline_idx{jj} = outline(xx_outline{jj}, yy_outline{jj});
                draw_outline(xx_outline{jj},yy_outline{jj},xy_outline_idx{jj},'k',5);
            end
        end

        if dosig
            xy_idx_pos{jj} = outline(x_out_pos{jj}, y_out_pos{jj});
            draw_outline(x_out_pos{jj},y_out_pos{jj},xy_idx_pos{jj},outline_color_pos,7)

            xy_idx_neg{jj} = outline(x_out_neg{jj}, y_out_neg{jj});
            draw_outline(x_out_neg{jj},y_out_neg{jj},xy_idx_neg{jj},outline_color_neg,7)
%             for i = 1:size(x_out{jj},1)
%                 for i2 = 1:4
%                     if xy_idx{jj}(i,i2)
%                         line([x_out{jj}(i,i2) x_out{jj}(i,i2+1)], [y_out{jj}(i,i2) y_out{jj}(i,i2+1)], 'color', 'y', 'linewidth', 5);
%                     end
%                 end
%             end
        end
    end
    if docolorbar
        hh = colorbar('southoutside');
        set(hh, 'fontSize', 25, 'lineWidth', 3);
    end

    info{jj}.vZ = vZ{jj};
end

end

% subfunctions
function xy_idx = outline(x, y)

xy_idx = true(size(x,1),4);

for i = 1:size(x,1)
    for ii = 1:4
        for j = i+1:size(x,1)
            for jj = 1:4
                if sum(abs([sort(x(i,(ii:ii+1))) sort(y(i,(ii:ii+1)))] - [sort(x(j,(jj:jj+1))), sort(y(j,(jj:jj+1)))]) < .0001) == 4
                    xy_idx(i,ii) = false; xy_idx(j,jj) = false;
                end
            end
        end
    end
end

end

function draw_outline(x,y,idx,color, linewidth)
for i = 1:size(x,1)
    for j = 1:4
        if idx(i,j)
            line([x(i,j) x(i,j+1)], [y(i,j) y(i,j+1)], 'color', color, 'linewidth', linewidth);
        end
    end
end

end

function r = combine_region(r)
rr = r(1);

if numel(r) > 1
    for i = 2:numel(r)
        rr.XYZ = [rr.XYZ r(i).XYZ];
        rr.XYZmm = [rr.XYZmm r(i).XYZmm];
        rr.val = [rr.val; r(i).val];
        rr.Z = [rr.Z r(i).Z];
    end
end

r = rr;

end

author: canlab

repo: CanlabCore

save_path: github-repos/MATLAB/canlab-CanlabCore


text:
function [P1,P2]=mg_ellblock(nelx,nely,x,y)
%mg_ellblock   prolongation for part of L-shaped and step domains
%   [P1,P2] = mg_ellblock(nelx,nely,x,y);
%   input
%          nelx    number of elements in x-direction
%          nely    number of elements in y-direction
%          x       x coordinate vector for coarse grid
%          y       y coordinate vector for coarse grid
%   output
%          P1      component of prolongation operator for upper left of
%                  step [-1,h]x[0,1]
%          P2      component of prolngation operator for border between
%                  upper left of step and right part of step
%
%   IFISS function: AR; 19 November, 2001.
% Copyright (c) 2005 D.J. Silvester, H.C. Elman, A. Ramage
ecx=nelx/2; ecy=nely/2;
tx = x(2:nelx+1) - x(1:nelx);
bx = zeros(nelx,1);
bx(1:2:nelx-1) = x(3:2:nelx+1) - x(1:2:nelx-1);
bx(2:2:nelx)   = bx(1:2:nelx-1);
dx = tx./bx;
%
ty = y(2:nely+1) - y(1:nely);
by = zeros(nely,1);
by(1:2:nely-1) = y(3:2:nely+1) - y(1:2:nely-1);
by(2:2:nely)   = by(1:2:nely-1);
dy = ty./by;
%
p = spalloc(2*ecx+1,ecx+1,3*ecx+1);
for j=2:ecx
   p(2*j-2,j) = dx(2*j-3);
   p(2*j-1,j) = 1.0;
   p(2*j  ,j) = dx(2*j);
end
j=1;
   p(2*j-1,1) = 1.0;
   p(2*j  ,1) = dx(2*j);
j=ecx+1;
   p(2*j-2,j) = dx(2*j-3);
   p(2*j-1,j) = 1.0;
%
% main block (last row/column omitted)
P1=spalloc( (2*ecy+1)*(2*ecx), (ecy+1)*(ecx), (3*ecy+1)*(3*ecx+1) );
%
sp=p(1:end-1,1:end-1);
for j=2:ecy,
    cols = (j-1)*(ecx)+[1:ecx];
    P1((2*j-3)*(2*ecx)+[1:2*ecx],cols) = dy(2*j-3)*sp;
    P1((2*j-2)*(2*ecx)+[1:2*ecx],cols) = sp;
    P1((2*j-1)*(2*ecx)+[1:2*ecx],cols) = dy(2*j)*sp;
end
j=1;
    cols = (j-1)*(ecx)+[1:ecx];
    P1((2*j-2)*(2*ecx)+[1:2*ecx],cols) = sp;
    P1((2*j-1)*(2*ecx)+[1:2*ecx],cols) = dy(2*j)*sp;
j=ecy+1;
    cols = (j-1)*(ecx)+[1:ecx];
    P1((2*j-3)*(2*ecx)+[1:2*ecx],cols) = dy(2*j-3)*sp;
    P1((2*j-2)*(2*ecx)+[1:2*ecx],cols) = sp;
%
% coupling with large block
P2=spalloc((2*ecy+1)*(2*ecx),(ecy+1)*(ecx+1),(3*ecy+1)*(3*ecx+1) );
entry=p(2*ecx,ecx+1);
for j=2:ecy,
    cols = (j-1)*(ecx+1)+1;
    P2((2*j-3)*(2*ecx)+[2*ecx],cols) = dy(2*j-3)*entry;
    P2((2*j-2)*(2*ecx)+[2*ecx],cols) = entry;
    P2((2*j-1)*(2*ecx)+[2*ecx],cols) = dy(2*j)*entry;
end
j=1;
    cols = (j-1)*(ecx+1)+1;
    P2((2*j-2)*(2*ecx)+[2*ecx],cols) = entry;
    P2((2*j-1)*(2*ecx)+[2*ecx],cols) = dy(2*j)*entry;
j=ecy+1;
    cols = (j-1)*(ecx+1)+1;
    P2((2*j-3)*(2*ecx)+[2*ecx],cols) = dy(2*j-3)*entry;
    P2((2*j-2)*(2*ecx)+[2*ecx],cols) = entry;
return

author: johannesgerer

repo: jburkardt-m


text:
function varargout = gui_producevisualization(varargin)
% GUI_PRODUCEVISUALIZATION MATLAB code for gui_producevisualization.fig
%      GUI_PRODUCEVISUALIZATION, by itself, creates a new GUI_PRODUCEVISUALIZATION or raises the existing
%      singleton*.
%
%      H = GUI_PRODUCEVISUALIZATION returns the handle to a new GUI_PRODUCEVISUALIZATION or the handle to
%      the existing singleton*.
%
%      GUI_PRODUCEVISUALIZATION('CALLBACK',hObject,eventData,handles,...) calls the local
%      function named CALLBACK in GUI_PRODUCEVISUALIZATION.M with the given input arguments.
%
%      GUI_PRODUCEVISUALIZATION('Property','Value',...) creates a new GUI_PRODUCEVISUALIZATION or raises the
%      existing singleton*.  Starting from the left, property value pairs are
%      applied to the GUI before gui_producevisualization_OpeningFcn gets called.  An
%      unrecognized property name or invalid value makes property application
%      stop.  All inputs are passed to gui_producevisualization_OpeningFcn via varargin.
%
%      *See GUI Options on GUIDE's Tools menu.  Choose "GUI allows only one
%      instance to run (singleton)".
%
% See also: GUIDE, GUIDATA, GUIHANDLES

% Edit the above text to modify the response to help gui_producevisualization

% Last Modified by GUIDE v2.5 17-Nov-2010 20:40:05

% Begin initialization code - DO NOT EDIT
gui_Singleton = 1;
gui_State = struct('gui_Name',       mfilename, ...
                   'gui_Singleton',  gui_Singleton, ...
                   'gui_OpeningFcn', @gui_producevisualization_OpeningFcn, ...
                   'gui_OutputFcn',  @gui_producevisualization_OutputFcn, ...
                   'gui_LayoutFcn',  [] , ...
                   'gui_Callback',   []);
if nargin && ischar(varargin{1})
    gui_State.gui_Callback = str2func(varargin{1});
end

if nargout
    [varargout{1:nargout}] = gui_mainfcn(gui_State, varargin{:});
else
    gui_mainfcn(gui_State, varargin{:});
end
% End initialization code - DO NOT EDIT

% --- Executes just before gui_producevisualization is made visible.
function gui_producevisualization_OpeningFcn(hObject, eventdata, handles, varargin)
% This function has no output args, see OutputFcn.
% hObject    handle to figure
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)
% varargin   command line arguments to gui_producevisualization (see VARARGIN)

% populate model selector listbox with admissible models
vars = evalin('base','whos');
vars = {vars(strcmp({vars.class},'struct')).name};
model_vars = {};
model_times = [];
% for all struct vars...
for v=vars
    % select those that are results
    var = evalin('base',v{1});
    if isfield(var,'tracking') && ~isfield(var,'pipelines') && isfield(var.tracking,'prediction_function')
        model_vars{end+1} = v{1};
        if isfield(var,'timestamp')
            model_times(end+1) = var.timestamp;
        else
            model_times(end+1) = 0;
        end
    end
end

if isempty(model_vars)
    % disable the OK button
    set(handles.popupmenu1,'String','(no model available)');
    set(handles.pushbutton1,'Enable','off');
else
    % add to pulldown menu...
    set(handles.popupmenu1,'String',model_vars);
end
guidata(hObject, handles);

% populate stream selector listbox with admissible streams
vars = evalin('base','whos');
vars = {vars(strcmp({vars.class},'struct')).name};
stream_vars = {};
% for all struct vars...
for v=vars
    % select those that are results
    var = evalin('base',v{1});
    if isfield(var,{'buffer','smax'})
        stream_vars{end+1} = v{1}; end
end

if isempty(stream_vars)
    % disable the OK button
    set(handles.popupmenu2,'String','(no stream available)');
    set(handles.pushbutton1,'Enable','off');
else
    % add to pulldown menu...
    set(handles.popupmenu2,'String',stream_vars);
end
guidata(hObject, handles);

% Choose default command line output for gui_producevisualization
handles.output = hObject;

% Update handles structure
guidata(hObject, handles);

% UIWAIT makes gui_producevisualization wait for user response (see UIRESUME)
uiwait(handles.figure1);

% --- Outputs from this function are returned to the command line.
function varargout = gui_producevisualization_OutputFcn(hObject, eventdata, handles)
% varargout  cell array for returning output args (see VARARGOUT);
% hObject    handle to figure
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

% Get default command line output from handles structure
try
    varargout{1} = handles.output;
    close(handles.figure1);
catch
    varargout = {};
end

function pushbutton1_Callback(hObject, eventdata, handles)
% ok pressed
% get model
choices = get(handles.popupmenu1,'String');
selection = get(handles.popupmenu1,'Value');
model = evalin('base',strtok(choices{selection}));
% get stream
choices = get(handles.popupmenu2,'String');
selection = get(handles.popupmenu2,'Value');
stream = strtok(choices{selection});
% get format
choices = get(handles.popupmenu3,'String');
selection = get(handles.popupmenu3,'Value');
fmt = choices{selection};
% start visualizer
run_producevisualization('pred_model',model,'pred_name',get(handles.edit2,'String'),'in_stream',stream, ...
    'vis_func',get(handles.edit3,'String'), 'freq',str2num(get(handles.edit1,'String')), 'out_form',fmt, ...
    'create_fig',true,'start_delay',1);
% resume gui
uiresume(handles.figure1);

function pushbutton5_Callback(hObject, eventdata, handles)
% help pressed
disp('coming soon...');

function pushbutton4_Callback(hObject, eventdata, handles)
% cancel pressed
uiresume(handles.figure1);

function edit1_Callback(hObject, eventdata, handles)
gui_check_scalar(hObject,'10');

function edit2_Callback(hObject, eventdata, handles)
gui_check_varname(hObject,'lastpredictor');

% --- auto-generated junk ---

function popupmenu1_Callback(hObject, eventdata, handles)
function popupmenu1_CreateFcn(hObject, eventdata, handles)
function popupmenu2_Callback(hObject, eventdata, handles)
function popupmenu2_CreateFcn(hObject, eventdata, handles)
function edit1_CreateFcn(hObject, eventdata, handles)
function popupmenu3_Callback(hObject, eventdata, handles)
function popupmenu3_CreateFcn(hObject, eventdata, handles)
function edit2_CreateFcn(hObject, eventdata, handles)
function edit3_Callback(hObject, eventdata, handles)
function edit3_CreateFcn(hObject, eventdata, handles)

repo: REST


text:
function [O,lo] = update_order_greedy(A)
    % greedy method for updating order
    K = size(A,2);
    if ~(size(A,1) == K)
        A = A'*A;
    end

    A = A>0;
    A(1:K+1:K^2) = 0;

    O{1} = 1;
    lo(1) = 1;
    for i = 2:K
       cont = true;
       for j = 1:length(O)
           if ~any(A(i,O{j}))
               O{j} = [O{j},i];
               lo(j) = lo(j) + 1;
               cont = false;
               break;
           end
       end
       if cont
           O{length(O)+1} = j;
           lo = [lo,1];
       end
    end
end

author: flatironinstitute

repo: CaImAn-MATLAB


text:
function value = sphere_shell_volume_nd ( n, r1, r2 )

%*****************************************************************************80
%
%% SPHERE_SHELL_VOLUME_ND computes the volume of a spherical shell in ND.
%
%  Integration region:
%
%    Points X(1:N) such that:
%
%      R1**2 <= Sum ( X(1:N) - XC(1:N) )**2 <= R2**2.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    20 May 2004
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, integer N, the dimension of the space.
%
%    Input, real R1, R2, the radiuses of the inner and
%    outer spheres.
%
%    Output, real SPHERE_SHELL_VOLUME_ND, the volume of the
%    spherical shell.
%
  value = ball_volume_nd ( n, r2 ) - ball_volume_nd ( n, r1 );

  return
end


text:
function pass = test_revolution(~)

%% Check results for a cone
f = chebfun(@(x) x, [0,1]);
result = cheb.revolution(f);

% Known exact results
exact.surfaceArea = pi*sqrt(2);
exact.volume = pi/3;
exact.centroidZ = 3/4;
exact.momentOfInertia = pi/10;

% Happy?
pass(1) = abs(result.surfaceArea - exact.surfaceArea) < 1e-15;
pass(2) = abs(result.volume - exact.volume) < 1e-15;
pass(3) = abs(result.centroidZ - exact.centroidZ) < 1e-15;
pass(4) = abs(result.momentOfInertia - exact.momentOfInertia) < 1e-15;

%% Check results for a surface of revolution on an unbounded domain
f = chebfun(@(x) exp(-x), [0,Inf]);
warnstate = warning;
warning('off', 'CHEBFUN:UNBNDFUN:sum:slowDecay')
result = cheb.revolution(f);
warning(warnstate);

% Known exact results
exact.surfaceArea = pi*(sqrt(2) + asinh(1));
exact.volume = pi/2;
exact.centroidZ = 1/2;
exact.momentOfInertia = pi/8;

% Happy?
pass(5) = abs(result.surfaceArea - exact.surfaceArea) < 1e-12;
pass(6) = abs(result.volume - exact.volume) < 1e-12;
pass(7) = abs(result.centroidZ - exact.centroidZ) < 1e-6;
pass(8) = abs(result.momentOfInertia - exact.momentOfInertia) < 1e-12;

end

author: chebfun

repo: chebfun

save_path: github-repos/MATLAB/chebfun-chebfun


text:
function rules1
global x y z alpha

xn=str2num(get(x,'String'));
yn=str2num(get(y,'String'));
zn=str2num(get(z,'String'));

alphan=str2num(get(alpha,'String'));

rules(xn,yn,zn,alphan);

author: Sable

repo: mcbench-benchmarks


text:
%DEMOSAICING  Demosaicing algorithm
%
%     dst = cv.demosaicing(src, code)
%     dst = cv.demosaicing(..., 'OptionName', optionValue, ...)
%
% ## Input
% * __src__ source image, 1-channel, `uint8` or `uint16` depth.
% * __code__ color conversion code. The following strings are supported:
%   * For Bayer->Gray demosaicing:
%     * __BayerBG2GRAY__
%     * __BayerGB2GRAY__
%     * __BayerRG2GRAY__
%     * __BayerGR2GRAY__
%   * For Bayer->RGB demosaicing:
%     * __BayerBG2BGR__, __BayerBG2RGB__
%     * __BayerGB2BGR__, __BayerGB2RGB__
%     * __BayerRG2BGR__, __BayerRG2RGB__
%     * __BayerGR2BGR__, __BayerGR2RGB__
%   * For Bayer->RGB demosaicing using Variable Number of Gradients:
%     * **BayerBG2BGR_VNG**, **BayerBG2RGB_VNG**
%     * **BayerGB2BGR_VNG**, **BayerGB2RGB_VNG**
%     * **BayerRG2BGR_VNG**, **BayerRG2RGB_VNG**
%     * **BayerGR2BGR_VNG**, **BayerGR2RGB_VNG**
%   * For Bayer->RGB Edge-Aware demosaicing:
%     * **BayerBG2BGR_EA**, **BayerBG2RGB_EA**
%     * **BayerGB2BGR_EA**, **BayerGB2RGB_EA**
%     * **BayerRG2BGR_EA**, **BayerRG2RGB_EA**
%     * **BayerGR2BGR_EA**, **BayerGR2RGB_EA**
%
% ## Output
% * __dst__ output image of same row-/col- size and depth as `src`, and of
%   specified number of channels (see `Channels` option).
%
% ## Options
% * __Channels__ Number of channels. If <= 0, automatically determined based
%   on `code`. default 0
%
% See also: cv.cvtColor, demosaic, vision.DemosaicInterpolator
%

author: kyamagu

repo: mexopencv

save_path: github-repos/MATLAB/kyamagu-mexopencv


text:
function [overlaps, score, class_ids, avg_cube, avg_cube_per_class] = detector_wiggle_from_iou_score(filename)
  ld = load(filename, 'ious', 'scores', 'classes');

  avg_cube = [];
  avg_cube_per_class = [];
  score = ld.scores;
  overlaps = ld.ious;
  class_ids = ld.classes;
end

author: hosang

repo: detection-proposals


text: puntos=draw('UMH');

repo: ARTE


text:
function count = compute_counts(data, sz)
% COMPUTE_COUNTS Count the number of times each combination of discrete assignments occurs
% count = compute_counts(data, sz)
%
% data(i,t) is the value of variable i in case t
% sz(i) : values for variable i are assumed to be in [1:sz(i)]
%
% Example: to compute a transition matrix for an HMM from a sequence of labeled states:
% transmat = mk_stochastic(compute_counts([seq(1:end-1); seq(2:end)], [nstates nstates]));

assert(length(sz) == size(data, 1));
P = prod(sz);
indices = subv2ind(sz, data'); % each row of data' is a case
%count = histc(indices, 1:P);
count = hist(indices, 1:P);
count = myreshape(count, sz);

author: fieldtrip

repo: fieldtrip

save_path: github-repos/MATLAB/fieldtrip-fieldtrip


text:
%  create auxiliary parameters to be used in the perturbation approach by
%  Foerster, Rubio-Ramirez, Waggoner and Zha(2016)
%

author: jmaih


text:
function [Y, optinf] = bpdn(D, S, lambda, opt)

% bpdn -- Basis Pursuit DeNoising
%
%         argmin_x (1/2)||D*x - s||_2^2 + lambda*||x||_1
%
%         The solution of the BPDN problem (see chen-1998-atomic) is
%         computed using the ADMM approach (see boyd-2010-distributed).
%
% Usage:
%       [Y, optinf] = bpdn(D, S, lambda, opt)
%
% Input:
%       D           Dictionary matrix
%       S           Signal vector (or matrix)
%       lambda      Regularization parameter
%       opt         Options/algorithm parameters structure (see below)
%
% Output:
%       Y           Dictionary coefficient vector (or matrix)
%       optinf      Details of optimisation
%
%
% Options structure fields:
%   Verbose           Flag determining whether iteration status is displayed.
%                     Fields are iteration number, functional value,
%                     data fidelity term, l1 regularisation term, and
%                     primal and dual residuals (see Sec. 3.3 of
%                     boyd-2010-distributed). The value of rho is also
%                     displayed if options request that it is automatically
%                     adjusted.
%   MaxMainIter       Maximum main iterations
%   AbsStopTol        Absolute convergence tolerance (see Sec. 3.3.1 of
%                     boyd-2010-distributed)
%   RelStopTol        Relative convergence tolerance (see Sec. 3.3.1 of
%                     boyd-2010-distributed)
%   L1Weight          Weighting array for coefficients in l1 norm of X
%   Y0                Initial value for Y
%   U0                Initial value for U
%   rho               ADMM penalty parameter
%   AutoRho           Flag determining whether rho is automatically updated
%                     (see Sec. 3.4.1 of boyd-2010-distributed)
%   AutoRhoPeriod     Iteration period on which rho is updated
%   RhoRsdlRatio      Primal/dual residual ratio in rho update test
%   RhoScaling        Multiplier applied to rho when updated
%   AutoRhoScaling    Flag determining whether RhoScaling value is
%                     adaptively determined (see wohlberg-2015-adaptive). If
%                     enabled, RhoScaling specifies a maximum allowed
%                     multiplier instead of a fixed multiplier.
%   RhoRsdlTarget     Residual ratio targeted by auto rho update policy.
%   StdResiduals      Flag determining whether standard residual definitions
%                     (see Sec 3.3 of boyd-2010-distributed) are used instead
%                     of normalised residuals (see wohlberg-2015-adaptive)
%   RelaxParam        Relaxation parameter (see Sec. 3.4.3 of
%                     boyd-2010-distributed)
%   NonNegCoef        Flag indicating whether solution should be forced to
%                     be non-negative
%   AuxVarObj         Flag determining whether objective function is computed
%                     using the auxiliary (split) variable
%
%
% Author: Brendt Wohlberg <brendt@lanl.gov>  Modified: 2015-07-23
%
% This file is part of the SPORCO library. Details of the copyright
% and user license can be found in the 'License' file distributed with
% the library.

if nargin < 4,
  opt = [];
end
checkopt(opt, defaultopts([]));
opt = defaultopts(opt);

% Default lambda is 1/10 times the lambda value beyond which the
% solution is a zero vector
if nargin < 3 | isempty(lambda),
  lambda = 0.1*max(vec(abs(D'*S)));
end

% Set up status display for verbose operation
hstr = 'Itn   Fnc       DFid      l1        r         s      ';
sfms = '%4d %9.2e %9.2e %9.2e %9.2e %9.2e';
nsep = 54;
if opt.AutoRho,
  hstr = [hstr '   rho   '];
  sfms = [sfms ' %9.2e'];
  nsep = nsep + 10;
end
if opt.Verbose && opt.MaxMainIter > 0,
  disp(hstr);
  disp(char('-' * ones(1,nsep)));
end

% Start timer
tstart = tic;

% Set up algorithm parameters and initialise variables
rho = opt.rho;
if isempty(rho), rho = 50*lambda+1; end;
[Nr, Nc] = size(D);
Nm = size(S,2);
Nx = Nc*Nm;
DTS = D'*S;
[luL, luU] = factorise(D, rho);
optinf = struct('itstat', [], 'opt', opt);
r = Inf;
s = Inf;
epri = 0;
edua = 0;

% Initialise main working variables
X = [];
if isempty(opt.Y0),
  Y = zeros(Nc,Nm);
else
  Y = opt.Y0;
end
Yprv = Y;
if isempty(opt.U0),
  if isempty(opt.Y0),
    U = zeros(Nc,Nm);
  else
    U = (lambda/rho)*sign(Y);
  end
else
  U = opt.U0;
end

% Main loop
k = 1;
while k <= opt.MaxMainIter && (r > epri | s > edua),

  % Solve X subproblem
  X = linsolve(D, rho, luL, luU, DTS + rho*(Y - U));

  % See pg. 21 of boyd-2010-distributed
  if opt.RelaxParam == 1,
    Xr = X;
  else
    Xr = opt.RelaxParam*X + (1-opt.RelaxParam)*Y;
  end

  % Solve Y subproblem
  Y = shrink(Xr + U, (lambda/rho)*opt.L1Weight);
  if opt.NonNegCoef,
    Y(Y < 0) = 0;
  end

  % Update dual variable
  U = U + Xr - Y;

  % Objective function and convergence measures
  if opt.AuxVarObj,
    Jdf = sum(vec(abs(D*Y - S).^2))/2;
    Jl1 = sum(abs(vec(bsxfun(@times, opt.L1Weight, Y))));
  else
    Jdf = sum(vec(abs(D*X - S).^2))/2;
    Jl1 = sum(abs(vec(bsxfun(@times, opt.L1Weight, X))));
  end
  Jfn = Jdf + lambda*Jl1;

  nX = norm(X(:)); nY = norm(Y(:)); nU = norm(U(:));
  if opt.StdResiduals,
    % See pp. 19-20 of boyd-2010-distributed
    r = norm(vec(X - Y));
    s = norm(vec(rho*(Yprv - Y)));
    epri = sqrt(Nx)*opt.AbsStopTol+max(nX,nY)*opt.RelStopTol;
    edua = sqrt(Nx)*opt.AbsStopTol+rho*nU*opt.RelStopTol;
  else
    % See wohlberg-2015-adaptive
    r = norm(vec(X - Y))/max(nX,nY);
    s = norm(vec(Yprv - Y))/nU;
    epri = sqrt(Nx)*opt.AbsStopTol/max(nX,nY)+opt.RelStopTol;
    edua = sqrt(Nx)*opt.AbsStopTol/(rho*nU)+opt.RelStopTol;
  end

  % Record and display iteration details
  tk = toc(tstart);
  optinf.itstat = [optinf.itstat; [k Jfn Jdf Jl1 r s epri edua rho tk]];
  if opt.Verbose,
    if opt.AutoRho,
      disp(sprintf(sfms, k, Jfn, Jdf, Jl1, r, s, rho));
    else
      disp(sprintf(sfms, k, Jfn, Jdf, Jl1, r, s));
    end
  end

  % See wohlberg-2015-adaptive and pp. 20-21 of boyd-2010-distributed
  if opt.AutoRho,
    if k ~= 1 && mod(k, opt.AutoRhoPeriod) == 0,
      if opt.AutoRhoScaling,
        rhomlt = sqrt(r/(s*opt.RhoRsdlTarget));
        if rhomlt < 1, rhomlt = 1/rhomlt; end
        if rhomlt > opt.RhoScaling, rhomlt = opt.RhoScaling; end
      else
        rhomlt = opt.RhoScaling;
      end
      rsf = 1;
      if r > opt.RhoRsdlTarget*opt.RhoRsdlRatio*s, rsf = rhomlt; end
      if s > (opt.RhoRsdlRatio/opt.RhoRsdlTarget)*r, rsf = 1/rhomlt; end
      rho = rsf*rho;
      U = U/rsf;
      if rsf ~= 1,
        [luL, luU] = factorise(D, rho);
      end
    end
  end

  Yprv = Y;
  k = k + 1;

end

% Record run time and working variables
optinf.runtime = toc(tstart);
optinf.X = X;
optinf.Y = Y;
optinf.U = U;
optinf.lambda = lambda;
optinf.rho = rho;

% End status display for verbose operation
if opt.Verbose && opt.MaxMainIter > 0,
  disp(char('-' * ones(1,nsep)));
end

return

function u = vec(v)

  u = v(:);

return

function u = shrink(v, lambda)

  if isscalar(lambda),
    u = sign(v).*max(0, abs(v) - lambda);
  else
    u = sign(v).*max(0, bsxfun(@minus, abs(v), lambda));
  end

return

function [L,U] = factorise(A, c)

  [N,M] = size(A);
  % If N < M it is cheaper to factorise A*A' + cI and then use the
  % matrix inversion lemma to compute the inverse of A'*A + cI
  if N >= M,
    [L,U] = lu(A'*A + c*eye(M,M));
  else
    [L,U] = lu(A*A' + c*eye(N,N));
  end

return

function x = linsolve(A, c, L, U, b)

  [N,M] = size(A);
  if N >= M,
    x = U \ (L \ b);
  else
    x = (b - A'*(U \ (L \ (A*b))))/c;
  end

return

function opt = defaultopts(opt)

  if ~isfield(opt,'Verbose'),
    opt.Verbose = 0;
  end
  if ~isfield(opt,'MaxMainIter'),
    opt.MaxMainIter = 1000;
  end
  if ~isfield(opt,'AbsStopTol'),
    opt.AbsStopTol = 0;
  end
  if ~isfield(opt,'RelStopTol'),
    opt.RelStopTol = 1e-4;
  end
  if ~isfield(opt,'L1Weight'),
    opt.L1Weight = 1;
  end
  if ~isfield(opt,'Y0'),
    opt.Y0 = [];
  end
  if ~isfield(opt,'U0'),
    opt.U0 = [];
  end
  if ~isfield(opt,'rho'),
    opt.rho = [];
  end
  if ~isfield(opt,'AutoRho'),
    opt.AutoRho = 1;
  end
  if ~isfield(opt,'AutoRhoPeriod'),
    opt.AutoRhoPeriod = 10;
  end
  if ~isfield(opt,'RhoRsdlRatio'),
    opt.RhoRsdlRatio = 1.2;
  end
  if ~isfield(opt,'RhoScaling'),
    opt.RhoScaling = 100;
  end
  if ~isfield(opt,'AutoRhoScaling'),
    opt.AutoRhoScaling = 1;
  end
  if ~isfield(opt,'RhoRsdlTarget'),
    opt.RhoRsdlTarget = 1;
  end
  if ~isfield(opt,'StdResiduals'),
    opt.StdResiduals = 0;
  end
  if ~isfield(opt,'RelaxParam'),
    opt.RelaxParam = 1.8;
  end
  if ~isfield(opt,'NonNegCoef'),
    opt.NonNegCoef = 0;
  end
  if ~isfield(opt,'AuxVarObj'),
    opt.AuxVarObj = 1;
  end

return

author: thfylsty


text:
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%   PARAMETERS Returns a data structure containing the parameters of the
%   ABB IRB1600iD.
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Authors:Daniel Vivancos Unica
%        Jose David Martinez Exposito
%        Maria Jose Martinez Liza
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% Copyright (C) 2012, by Arturo Gil Aparicio
%
% This file is part of ARTE (A Robotics Toolbox for Education).
%
% ARTE is free software: you can redistribute it and/or modify
% it under the terms of the GNU Lesser General Public License as published by
% the Free Software Foundation, either version 3 of the License, or
% (at your option) any later version.
%
% ARTE is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
% GNU Lesser General Public License for more details.
%
% You should have received a copy of the GNU Leser General Public License
% along with ARTE.  If not, see <http://www.gnu.org/licenses/>.

function robot = parameters()

robot.name= 'ABB_IRB1600ID';

%Path where everything is stored for this robot
robot.path = 'robots/abb/IRB1600ID';

robot.DH.theta= '[q(1) q(2)-pi/2 q(3) q(4) q(5) q(6)+pi]';
robot.DH.d='[0.4865 0 0 0.640 0 0.2]';
robot.DH.a='[0.15 0.7 0.11 0 0 0]';
robot.DH.alpha= '[-pi/2 0 -pi/2 pi/2 -pi/2 0]';
robot.J=[];

robot.inversekinematic_fn = 'inversekinematic_irb1600id(robot, T)';

%number of degrees of freedom
robot.DOF = 6;

%rotational: 0, translational: 1
robot.kind=['R' 'R' 'R' 'R' 'R' 'R'];

%minimum and maximum rotation angle in rad
robot.maxangle =[-pi pi; %Axis 1, minimum, maximum
                deg2rad(-90) deg2rad(150); %Axis 2, minimum, maximum
                deg2rad(-238) deg2rad(79); %Axis 3
                deg2rad(-155) deg2rad(155); %Axis 4
                deg2rad(-90) deg2rad(135); %Axis 5
                deg2rad(-200) deg2rad(200)]; %Axis 6

%maximum absolute speed of each joint rad/s or m/s
robot.velmax = [deg2rad(180); %Axis 1, rad/s
                deg2rad(180); %Axis 2, rad/s
                deg2rad(180); %Axis 3, rad/s
                deg2rad(320); %Axis 4, rad/s
                deg2rad(380); %Axis 5, rad/s
                deg2rad(460)];%Axis 6, rad/s
robot.accelmax=robot.velmax/0.1; % 0.1 is here an acceleration time

            % end effectors maximum velocity
robot.linear_velmax = 2.5; %m/s

%base reference system
robot.T0 = eye(4);

%INITIALIZATION OF VARIABLES REQUIRED FOR THE SIMULATION
%position, velocity and acceleration
robot=init_sim_variables(robot);

% GRAPHICS
robot.graphical.has_graphics=1;
robot.graphical.color = [255 102 51]./255;
%for transparency
robot.graphical.draw_transparent=0;
%draw DH systems
robot.graphical.draw_axes=1;
%DH system length and Font size, standard is 1/10. Select 2/20, 3/30 for
%bigger robots
robot.graphical.axes_scale=1;
%adjust for a default view of the robot
robot.axis=[-2 2 -2 2 0 2.2];
%read graphics files
robot = read_graphics(robot);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%DYNAMICS

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
robot.has_dynamics=1;

%consider friction in the computations
robot.dynamics.friction=0;

%link masses (kg)
%Tenemos que repartir 250kg (peso total robot IRB1600ID entre los 6 eslabones)
%Utilizando como material Aluminio, los datos obtenidos del programa Inventor seran:

robot.dynamics.masses=[0 109.440 28.242 31.079 6.936 1.430 0.249];
%COM of each link with respect to own reference system
robot.dynamics.r_com=[0.052       -0.012    0.343; %(rx, ry, rz) link 1
         0.150       0.791     -0.183;%(rx, ry, rz) link 2
         1.172	     0.140	   0.017; %(rx, ry, rz) link 3
         1.296       0.063     0.561; %(rx, ry, rz) link 4
         1.296       0.870     0.027; %(rx, ry, rz) link 5
         0           1.296     0.966];%(rx, ry, rz) link 6

%Inertia matrices of each link with respect to its D-H reference system.
% Ixx	Iyy	Izz	Ixy	Iyz	Ixz, for each row
robot.dynamics.Inertia=[0      0	    0   	0	0	0;
         2.305  2.696	3.802	0	0	0;
         .0700	1.765	1.784	0	0	0;
         0.267  0.324	0.332	0	0	0;
         0.297	0.290	.023	0	0	0;
         .008	.002	.008	0   0	0;
         0      0       0       0   0   0];
%Los momentos de Inercia y los centros de gravedad quedan corregidos de los
%obtenidos del programa Inventor para ser ajustados al sistema de
%coordenadas designado en nuestro robot.

robot.motors=load_motors([5 5 5 4 4 4]);
%Speed reductor at each joint
robot.motors.G=[300 300 300 300 300 300];


text:
function varargout = importerpGUI(varargin)

% Begin initialization code - DO NOT EDIT
gui_Singleton = 1;
gui_State = struct('gui_Name',       mfilename, ...
      'gui_Singleton',  gui_Singleton, ...
      'gui_OpeningFcn', @importerpGUI_OpeningFcn, ...
      'gui_OutputFcn',  @importerpGUI_OutputFcn, ...
      'gui_LayoutFcn',  [] , ...
      'gui_Callback',   []);
if nargin && ischar(varargin{1})
      gui_State.gui_Callback = str2func(varargin{1});
end

if nargout
      [varargout{1:nargout}] = gui_mainfcn(gui_State, varargin{:});
else
      gui_mainfcn(gui_State, varargin{:});
end
% End initialization code - DO NOT EDIT

% -------------------------------------------------------------------------
function importerpGUI_OpeningFcn(hObject, eventdata, handles, varargin)

handles.output   = [];
handles.fulltext = [];
handles.iniftype = [];

%{filename, pathname, ftype,includetime,timeunit,elabel,transpose,fs,xlim});
try
      def = varargin{1};
catch
      def = {'','','',0,1,0,0,1000,[-200 800]};
end

filename     = def{1};
pathname     = def{2};
ftype        = def{3};
includetime  = def{4};
timeunit     = def{5};
elabel       = def{6};
transpose    = def{7};
fs           = def{8};
xlim         = def{9};

%
% Set objects
%
nfile = length(filename);
fulltext = cell(1);

if nfile>0
      for i=1:nfile
            fulltext{i} = fullfile(pathname{i},filename{i});
      end
      fulltext  = cat(1, fulltext', {'new erpset'});
      set(handles.listbox_erpnames, 'String', fulltext');
else
      set(handles.listbox_erpnames, 'String', {'new erpset'});
end

set(handles.checkbox_time, 'Value', includetime) % time exists by default
set(handles.edit_srate,'String', num2str(fs))
set(handles.edit_xlim,'String', num2str(xlim))
if includetime && timeunit==1
      set(handles.radiobutton_time_sec,'Value',1)
      set(handles.radiobutton_time_msec,'Value',0)
elseif includetime && timeunit==1E-3
      set(handles.radiobutton_time_sec,'Value',0)
      set(handles.radiobutton_time_msec,'Value',1)
elseif ~includetime
      set(handles.radiobutton_time_sec,'Value',0)
      set(handles.radiobutton_time_msec,'Value',0)
      set(handles.radiobutton_time_sec,'Enable', 'off')
      set(handles.radiobutton_time_msec,'Enable', 'off')
else
      msgboxText =  ['No valid time unit.\n'...
            'Only 1 (sec) and 1E-3 (msec) are allowed.'];
      title = 'ERPLAB: importerpGUI() -> wrong input';
      errorfound(sprintf(msgboxText), title);
      return
end

set(handles.checkbox_elabels, 'Value', elabel) % time exists by default

if transpose==1
      set(handles.radiobutton_er_pc,'Value',0)
      set(handles.radiobutton_pr_ec,'Value',1)
else
      set(handles.radiobutton_er_pc,'Value',1)
      set(handles.radiobutton_pr_ec,'Value',0)
end

%
% Name & version
%
version = geterplabversion;
set(handles.gui_chassis,'Name', ['ERPLAB ' version '   -   IMPORT AVERAGED ERP DATA GUI'])

%
% Color GUI
%
handles = painterplab(handles);

%
% Set font size
%
handles = setfonterplab(handles);

% Update handles structure
guidata(hObject, handles);

% help
helpbutton

% UIWAIT makes grandaveragerGUI wait for user response (see UIRESUME)
uiwait(handles.gui_chassis);

% -------------------------------------------------------------------------
function varargout = importerpGUI_OutputFcn(hObject, eventdata, handles)

% Get default command line output from handles structure
varargout{1} = handles.output;

% The figure can be deleted now
delete(handles.gui_chassis);
pause(0.1)

% -------------------------------------------------------------------------
function listbox_erpnames_Callback(hObject, eventdata, handles)

% -------------------------------------------------------------------------
function listbox_erpnames_CreateFcn(hObject, eventdata, handles)

if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
      set(hObject,'BackgroundColor','white');
end

% -------------------------------------------------------------------------
function pushbutton_adderpfile_Callback(hObject, eventdata, handles)

% [erpfname, erppathname, findex] = uigetfile({'*.txt','Universal text format(*.txt)'; ...
%       '*.avg','Neuroscan AVG format (*.avg)'; ...
%       '*.txt','ERPSS text format (*.txt)'; ...
%       '*.*',  'All Files (*.*)'}, ...
%       'Select an edited file', ...
%       'MultiSelect', 'on');

[erpfname, erppathname, findex] = uigetfile({'*.txt','Universal text format(*.txt)'; ...
      '*.*',  'All Files (*.*)'}, ...
      'Select an edited file', ...
      'MultiSelect', 'on');

if isequal(erpfname,0)
      disp('User selected Cancel')
      return
else
      iniftype = handles.iniftype;
      if ~isempty(iniftype)
            if findex~= iniftype
                  fprintf('\nWARNING: You are putting together different formats of averaged files.\n\n')
            end
      else
            handles.iniftype = findex;
      end

      if ~iscell(erpfname)
            erpfname = {erpfname};
      end

      nerpn = length(erpfname);

      for i=1:nerpn

            newline  = fullfile(erppathname, erpfname{i});      % new member
            currline = get(handles.listbox_erpnames, 'Value');  % current position at the list
            fulltext = get(handles.listbox_erpnames, 'String'); % current list
            indxline = length(fulltext); % last line (virtual)

            if currline==indxline % is the current line the last line?
                  % extra line forward
                  fulltext  = cat(1, fulltext, {'new erpset'});
                  set(handles.listbox_erpnames, 'Value', currline+1)
            else
                  set(handles.listbox_erpnames, 'Value', currline)
                  resto = fulltext(currline:indxline);
                  fulltext  = cat(1, fulltext, {'new erpset'});
                  set(handles.listbox_erpnames, 'Value', currline+1)
                  [fulltext{currline+1:indxline+1}] = resto{:};
            end

            fulltext{currline} = sprintf('Bin %g\t@\t%s',currline, newline);
            set(handles.listbox_erpnames, 'String', fulltext)
      end

      indxline = length(fulltext);
      handles.indxline = indxline;
      handles.fulltext = fulltext;
      handles.listname = [];
      %set(handles.edit_filelist,'String','');
      % Update handles structure
      guidata(hObject, handles);
end

% -------------------------------------------------------------------------
function pushbutton_delerpset_Callback(hObject, eventdata, handles)

fulltext = get(handles.listbox_erpnames, 'String'); % current list
indxline = length(fulltext); % last line
fulltext = char(fulltext); % string matrix
currline = get(handles.listbox_erpnames, 'Value'); % current pos

if currline>=1 && currline<indxline

      fulltext(currline,:) = [];
      fulltext = cellstr(fulltext); % cell string
      nfile = length(fulltext);

      if nfile>1
            fulltext = regexprep(fulltext,'Bin.*?@','', 'ignorecase');
            fulltext = strtrim(fulltext);

            for i=1:nfile-1;
                  fulltext{i} = sprintf('Bin %g\t@\t%s',i, fulltext{i});
            end
      else
            handles.iniftype = [];
      end

      set(handles.listbox_erpnames, 'String', fulltext);
      listbox_erpnames_Callback(hObject, eventdata, handles)

      handles.fulltext = fulltext;
      %indxline = length(fulltext);
      %handles.listname = [];
      %set(handles.edit_filelist,'String','');

      % Update handles structure
      %guidata(hObject, handles);
end

% -------------------------------------------------------------------------
function pushbutton_up_Callback(hObject, eventdata, handles)

fulltext = get(handles.listbox_erpnames, 'String'); % current list
indxline = length(fulltext); % last line
currline = get(handles.listbox_erpnames, 'Value'); % current pos

if currline>1 && currline<indxline

      auxft = fulltext{currline -1};
      fulltext{currline-1} = fulltext{currline};
      fulltext{currline}   = auxft;
      fulltext = regexprep(fulltext,'Bin.*?@','', 'ignorecase');% remove bin prefix
      fulltext = strtrim(fulltext);

      for i=1:length(fulltext)-1;
            fulltext{i} = sprintf('Bin %g\t@\t%s',i, fulltext{i}); % add bin prefix
      end

      set(handles.listbox_erpnames, 'String', fulltext);
      listbox_erpnames_Callback(hObject, eventdata, handles)
      set(handles.listbox_erpnames, 'Value', currline-1);
      handles.fulltext = fulltext;
      %indxline = length(fulltext);
      %handles.listname = [];
      %set(handles.edit_filelist,'String','');

      % Update handles structure
      %guidata(hObject, handles);
end

% -------------------------------------------------------------------------
function pushbutton_down_Callback(hObject, eventdata, handles)

fulltext = get(handles.listbox_erpnames, 'String'); % current list
indxline = length(fulltext); % last line
currline = get(handles.listbox_erpnames, 'Value'); % current pos

if currline>=1 && currline<indxline-1

      auxft = fulltext{currline + 1};
      fulltext{currline+1} = fulltext{currline};
      fulltext{currline}   = auxft;
      fulltext = regexprep(fulltext,'Bin.*?@','', 'ignorecase');% remove bin prefix
      fulltext = strtrim(fulltext);

      for i=1:length(fulltext)-1;
            fulltext{i} = sprintf('Bin %g\t@\t%s',i, fulltext{i});% add bin prefix
      end

      set(handles.listbox_erpnames, 'String', fulltext);
      listbox_erpnames_Callback(hObject, eventdata, handles)
      set(handles.listbox_erpnames, 'Value', currline+1);
      handles.fulltext = fulltext;
      %indxline = length(fulltext);
      %handles.listname = [];
      %set(handles.edit_filelist,'String','');

      % Update handles structure
      %guidata(hObject, handles);
end

% -------------------------------------------------------------------------
function ok_Callback(hObject, eventdata, handles)

fulltext = get(handles.listbox_erpnames, 'String'); % current list
fulltext = fulltext(1:end-1);
fulltext = regexprep(fulltext,'Bin.*?@','', 'ignorecase'); % remove bin prefix
fulltext = strtrim(fulltext);
filename = cell(1);
pathname = cell(1);
filetype = zeros(1, length(fulltext));

for i=1:length(fulltext)
      [pathstr, name, ext] = fileparts(fulltext{i});
      pathname{1,i} = pathstr;
      filename{1,i} = [name ext];

      if strcmpi(ext,'.txt')
            filetype(1,i) = 1; % text
      elseif strcmpi(ext,'.avg')
            filetype(1,i) = 2; % neuroscan
      else
            filetype(1,i) = 1; % text by default
      end
end

istime = get(handles.checkbox_time,'Value');

if get(handles.radiobutton_time_sec,'Value') && ~get(handles.radiobutton_time_msec,'Value')
      tunit = 1;
elseif ~get(handles.radiobutton_time_sec,'Value') && get(handles.radiobutton_time_msec,'Value')
      tunit = 1E-3;
elseif ~get(handles.radiobutton_time_sec,'Value') && ~get(handles.radiobutton_time_msec,'Value')
      tunit = [];
else
      msgboxText =  'Please, choose one option.';
      title = 'ERPLAB: importerpGUI() -> missing input';
      errorfound(sprintf(msgboxText), title);
      return
end

iselabel = get(handles.checkbox_elabels,'Value');

if get(handles.radiobutton_er_pc,'Value');
      transpose = 0;
end
if get(handles.radiobutton_pr_ec,'Value');
      transpose = 1;
end
srate = str2num(get(handles.edit_srate, 'String'));
xlim  = str2num(get(handles.edit_xlim, 'String'));

if (isempty(srate) || isempty(xlim)) && istime==0
      msgboxText =  ['Since there is not time info, you must specify\n'...
            'sampling rate in Hz, and time range (xmin xmax) in msec.'];
      title = 'ERPLAB: importerpGUI() -> missing input';
      errorfound(sprintf(msgboxText), title);
      return
end
if ~isempty(srate)
      if numel(srate)>1
            msgboxText =  'You must specify a unique value for sampling rate in Hz';
            title = 'ERPLAB: importerpGUI() -> wrong input';
            errorfound(sprintf(msgboxText), title);
            return
      end
      if srate<=0
            msgboxText =  'You must specify a positive integer (no zero) value for sampling rate in Hz';
            title = 'ERPLAB: importerpGUI() -> wrong input';
            errorfound(sprintf(msgboxText), title);
            return
      end
end
if ~isempty(xlim)
      if numel(xlim)~=2
            msgboxText =  'You must specify two values for time range [xmin xmax] in msec';
            title = 'ERPLAB: importerpGUI() -> wrong input';
            errorfound(sprintf(msgboxText), title);
            return
      end
      if xlim(1)>=xlim(2)
            msgboxText =  ['Time range is incorrect.\n'...
                  'xmin must be lesser than xmax'];
            title = 'ERPLAB: importerpGUI() -> wrong input';
            errorfound(sprintf(msgboxText), title);
            return
      end
      if xlim(1)>0 || xlim(2)<0
            msgboxText =  'time zero must be included within the time range.';
            title = 'ERPLAB: importerpGUI() -> wrong input';
            errorfound(sprintf(msgboxText), title);
            return
      end
end

handles.output = {filename pathname filetype istime tunit iselabel transpose srate xlim};
% Update handles structure
guidata(hObject, handles);
uiresume(handles.gui_chassis);

% -------------------------------------------------------------------------
function cancel_Callback(hObject, eventdata, handles)
handles.output = [];
% Update handles structure
guidata(hObject, handles);
uiresume(handles.gui_chassis);

% -------------------------------------------------------------------------
function edit1_Callback(hObject, eventdata, handles)

% -------------------------------------------------------------------------
function edit1_CreateFcn(hObject, eventdata, handles)

if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
      set(hObject,'BackgroundColor','white');
end

% -------------------------------------------------------------------------
function edit2_Callback(hObject, eventdata, handles)

% -------------------------------------------------------------------------
function edit2_CreateFcn(hObject, eventdata, handles)

if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
      set(hObject,'BackgroundColor','white');
end

% -------------------------------------------------------------------------
function pushbutton_CL_Callback(hObject, eventdata, handles)

fulltext = get(handles.listbox_erpnames, 'String'); % current list
indxline = length(fulltext); % last line

if indxline>1
      BackERPLABcolor = [1 0.9 0.3];    % yellow
      question = ['Current list will be deleted.\n\n'...
                  'Are you sure?'];
      title = 'List of averaged files for importing';
      oldcolor = get(0,'DefaultUicontrolBackgroundColor');
      set(0,'DefaultUicontrolBackgroundColor',BackERPLABcolor)
      button = questdlg(sprintf(question), title,'Cancel','Yes', 'Cancel');
      set(0,'DefaultUicontrolBackgroundColor',oldcolor)

      if ~strcmpi(button,'Yes')
            return
      end

      set(handles.listbox_erpnames, 'Value', 1);
      set(handles.listbox_erpnames, 'String', {'new erpset'});
      handles.fulltext = [];
      handles.iniftype = [];

      % Update handles structure
      guidata(hObject, handles);
end

% -------------------------------------------------------------------------
function radiobutton_er_pc_Callback(hObject, eventdata, handles)
if get(hObject,'Value')
      set(handles.radiobutton_pr_ec,'value', 0)
else
      set(handles.radiobutton_er_pc,'value', 1)
end

% -------------------------------------------------------------------------
function radiobutton_pr_ec_Callback(hObject, eventdata, handles)
if get(hObject,'Value')
      set(handles.radiobutton_er_pc,'value', 0)
else
      set(handles.radiobutton_pr_ec,'value', 1)
end

% -------------------------------------------------------------------------
function checkbox_elabels_Callback(hObject, eventdata, handles)

% -------------------------------------------------------------------------
function checkbox_time_Callback(hObject, eventdata, handles)

if get(hObject,'Value')
      %set(handles.edit_srate,'Enable','off')
      %set(handles.edit_xlim,'Enable','off')
      set(handles.radiobutton_time_sec,'Enable','on')
      set(handles.radiobutton_time_msec,'Enable','on')
      set(handles.edit_srate, 'Enable', 'off')
      set(handles.edit_xlim, 'Enable', 'off')
else
      %set(handles.edit_srate,'Enable','on')
      %set(handles.edit_xlim,'Enable','on')
      set(handles.radiobutton_time_sec,'Enable','off')
      set(handles.radiobutton_time_msec,'Enable','off')
      set(handles.edit_srate, 'Enable', 'on')
      set(handles.edit_xlim, 'Enable', 'on')
end

% -------------------------------------------------------------------------
function edit_srate_Callback(hObject, eventdata, handles)

% -------------------------------------------------------------------------
function edit_srate_CreateFcn(hObject, eventdata, handles)

if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
    set(hObject,'BackgroundColor','white');
end

% -------------------------------------------------------------------------
function edit_xlim_Callback(hObject, eventdata, handles)

% -------------------------------------------------------------------------
function edit_xlim_CreateFcn(hObject, eventdata, handles)

if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
    set(hObject,'BackgroundColor','white');
end

% -------------------------------------------------------------------------
function radiobutton_time_sec_Callback(hObject, eventdata, handles)

if get(hObject,'Value')
      set(handles.radiobutton_time_msec,'value',0)
else
      set(handles.radiobutton_time_sec,'value',1)
end

% -------------------------------------------------------------------------
function radiobutton_time_msec_Callback(hObject, eventdata, handles)
if get(hObject,'Value')
      set(handles.radiobutton_time_sec,'value',0)
else
      set(handles.radiobutton_time_msec,'value',1)
end

% -------------------------------------------------------------------------
function pushbutton_take_a_look_Callback(hObject, eventdata, handles)

fulltext = get(handles.listbox_erpnames, 'String'); % current list
indxline = length(fulltext); % last line
currline = get(handles.listbox_erpnames, 'Value'); % current pos

if currline>=1 && currline<indxline
      filename = fulltext{currline};
      filename = regexprep(filename,'Bin.*?@','', 'ignorecase');% remove bin prefix
      filename = strtrim(filename);

      [pathstr, name, ext] = fileparts(filename) ;

      if ismember_bc2({ext},{'.txt' '.asc' '.ascii' '.csv' '.csa' '.dat' ''})
            uiopen(filename,1)
      else
            msgboxText =  'Current line does not refer a text file.';
            title = 'ERPLAB: importerpGUI() -> wrong input';
            errorfound(sprintf(msgboxText), title);
            return
      end
else
      msgboxText =  'Current line does not contain a valid filename.';
      title = 'ERPLAB: importerpGUI() -> wrong input';
      errorfound(sprintf(msgboxText), title);
      return
end

% -------------------------------------------------------------------------
function pushbutton_help_Callback(hObject, eventdata, handles)
web https://github.com/lucklab/erplab/wiki/Manual -browser % pending

% -------------------------------------------------------------------------
function button_clearfile_Callback(hObject, eventdata, handles)
set(handles.edit5,'String','');
set(handles.button_savelist, 'Enable', 'off')
handles.listname = [];
% Update handles structure
guidata(hObject, handles);

function button_savelistas_Callback(hObject, eventdata, handles)
fulltext = char(get(handles.listbox_erpnames,'String'));

if length(fulltext)>1

        fullname = savelist(hObject, eventdata, handles);

        if isempty(fullname)
                return
        end

        set(handles.edit5, 'String', fullname )
        set(handles.button_savelist, 'Enable', 'on')
        handles.listname = fullname;
        % Update handles structure
        guidata(hObject, handles);
else
        set(handles.button_savelistas,'Enable','off')
        msgboxText =  'You have not specified any erpset!';
        title = 'ERPLAB: averager GUI few inputs';
        errorfound(msgboxText, title);
        set(handles.button_savelistas,'Enable','on')
        return
end

%--------------------------------------------------------------------------
function fullname = savelist(hObject, eventdata, handles)

fulltext  = char(strtrim(get(handles.listbox_erpnames,'String')));
pre_fname = char(strtrim(get(handles.edit5,'String')));
if isempty(pre_fname)
    pre_fname = '';
end

%
% Save OUTPUT file
%
[filename, filepath, filterindex] = uiputfile({'*.txt';'*.dat';'*.*'},'Save erpset list as', pre_fname);

if isequal(filename,0)
        disp('User selected Cancel')
        fullname =[];
        return
else

        [px, fname, ext] = fileparts(filename);

        if strcmp(ext,'')

                if filterindex==1 || filterindex==3
                        ext   = '.txt';
                else
                        ext   = '.dat';
                end
        end

        fname = [ fname ext];
        fullname = fullfile(filepath, fname);
        disp(['For saving erpset list, user selected <a href="matlab: open(''' fullname ''')">' fullname '</a>'])
        fid_list   = fopen( fullname , 'w');

        for i=1:size(fulltext,1)-1
                fprintf(fid_list,'%s\n', fulltext(i,:));
        end

        fclose(fid_list);
end

% -------------------------------------------------------------------------
function button_savelist_Callback(hObject, eventdata, handles)
fulltext = char(strtrim(get(handles.listbox_erpnames,'String')));

if length(fulltext)>1

        fullname = get(handles.edit5, 'String');

        if ~strcmp(fullname,'')

                fid_list   = fopen( fullname , 'w');

                for i=1:size(fulltext,1)
                        fprintf(fid_list,'%s\n', fulltext(i,:));
                end

                fclose(fid_list);
                handles.listname = fullname;

                % Update handles structure
                guidata(hObject, handles);
                disp(['Saving equation list at <a href="matlab: open(''' fullname ''')">' fullname '</a>'])
        else
                button_savelistas_Callback(hObject, eventdata, handles)
                return
        end
else
        set(handles.button_savelistas,'Enable','off')
        msgboxText =  'You have not written any formula yet!';
        title = 'ERPLAB: chanoperGUI few inputs';
        errorfound(msgboxText, title);
        set(handles.button_savelistas,'Enable','on')
        return
end

% -------------------------------------------------------------------------
function button_loadlist_Callback(hObject, eventdata, handles, optionx)
if nargin<4
        optionx=0;
end
if optionx==0
        [listname, lispath] = uigetfile({  '*.txt','Text File (*.txt)'; ...
                '*.*',  'All Files (*.*)'}, ...
                'Select an edited list', ...
                'MultiSelect', 'off');

        if isequal(listname,0)
                disp('User selected Cancel')
                return
        else
                fullname = fullfile(lispath, listname);
                disp(['For erpset list user selected  <a href="matlab: open(''' fullname ''')">' fullname '</a>'])
        end
else
        fullname = handles.erpset;
        if isnumeric(fullname)
                fullname = '';
        end
end

fid_list = fopen( fullname );
formcell = textscan(fid_list, '%[^\n]','CommentStyle','#', 'whitespace', '');
lista = formcell{:};
% extra line forward
lista   = cat(1, lista, {'new erpset'});
lentext = length(lista);
fclose(fid_list);

if lentext>1
        set(handles.listbox_erpnames,'String',lista);
        set(handles.edit5,'String',fullname);
        handles.listname = fullname;
        set(handles.button_savelist, 'Enable','on')

        % Update handles structure
        guidata(hObject, handles);
else
        msgboxText =  'This list is empty!';
        title = 'ERPLAB: geterpvaluesGUI inputs';
        errorfound(msgboxText, title);
        handles.listname = [];
        set(handles.button_savelist, 'Enable','off')

        % Update handles structure
        guidata(hObject, handles);
end
set(handles.listbox_erpnames,'String',lista);

% -------------------------------------------------------------------------
function edit5_Callback(hObject, eventdata, handles)

% -------------------------------------------------------------------------
function edit5_CreateFcn(hObject, eventdata, handles)

% Hint: edit controls usually have a white background on Windows.
%       See ISPC and COMPUTER.
if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
    set(hObject,'BackgroundColor','white');
end

%--------------------------------------------------------------------------
function gui_chassis_CloseRequestFcn(hObject, eventdata, handles)

if isequal(get(handles.gui_chassis, 'waitstatus'), 'waiting')
      %The GUI is still in UIWAIT, us UIRESUME
      handles.output = [];
      %Update handles structure
      guidata(hObject, handles);
      uiresume(handles.gui_chassis);
else
      % The GUI is no longer waiting, just close it
      delete(handles.gui_chassis);
end

author: ucdavis

repo: erplab

save_path: github-repos/MATLAB/ucdavis-erplab


text:
function DVSet = CorrelationAnalysis(Problem,Population,DV,nCor)
% Detect the group of each distance variable

%------------------------------- Copyright --------------------------------
% Copyright (c) 2023 BIMK Group. You are free to use the PlatEMO for
% research purposes. All publications which use this platform or any code
% in the platform should acknowledge the use of "PlatEMO" and reference "Ye
% Tian, Ran Cheng, Xingyi Zhang, and Yaochu Jin, PlatEMO: A MATLAB platform
% for evolutionary multi-objective optimization [educational forum], IEEE
% Computational Intelligence Magazine, 2017, 12(4): 73-87".
%--------------------------------------------------------------------------

    DVSet = {};
    for v = DV
        RelatedSet = [];
        for d = 1 : length(DVSet)
            for u = DVSet{d}
                drawnow('limitrate');
                sign = false;
                for i = 1 : nCor
                    p    = Population(randi(length(Population)));
                    a2   = unifrnd(Problem.lower(v),Problem.upper(v));
                    b2   = unifrnd(Problem.lower(u),Problem.upper(u));
                    decs = repmat(p.dec,3,1);
                    decs(1,v)     = a2;
                    decs(2,u)     = b2;
                    decs(3,[v,u]) = [a2,b2];
                    F = Problem.Evaluation(decs);
                    delta1 = F(1).obj - p.obj;
                    delta2 = F(3).obj - F(2).obj;
                    if any(delta1.*delta2<0)
                        sign = true;
                        RelatedSet = [RelatedSet,d];
                        break;
                    end
                end
                if sign
                    break;
                end
            end
        end
        if isempty(RelatedSet)
            DVSet = [DVSet,v];
        else
            DVSet = [DVSet,[cell2mat(DVSet(RelatedSet)),v]];
            DVSet(RelatedSet) = [];
        end
    end
end

author: BIMK

repo: PlatEMO

save_path: github-repos/MATLAB/BIMK-PlatEMO


text:
function navego_print_version()

versionstr = 'NaveGo, release v1.3';

fprintf('\n%s.\n', versionstr)

end

author: rodralez

repo: NaveGo

save_path: github-repos/MATLAB/rodralez-NaveGo


text:
function mesh = prepare_mesh_tetrahedral(cfg, mri)

% PREPARE_MESH_TETRAHEDRAL
%
% See also PREPARE_MESH_MANUAL, PREPARE_MESH_HEADSHAPE,
% PREPARE_MESH_HEXAHEDRAL, PREPARE_MESH_SEGMENTATION

% Copyrights (C) 2016, Robert Oostenveld
%
% This file is part of FieldTrip, see http://www.fieldtriptoolbox.org
% for the documentation and details.
%
%    FieldTrip is free software: you can redistribute it and/or modify
%    it under the terms of the GNU General Public License as published by
%    the Free Software Foundation, either version 3 of the License, or
%    (at your option) any later version.
%
%    FieldTrip is distributed in the hope that it will be useful,
%    but WITHOUT ANY WARRANTY; without even the implied warranty of
%    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
%    GNU General Public License for more details.
%
%    You should have received a copy of the GNU General Public License
%    along with FieldTrip. If not, see <http://www.gnu.org/licenses/>.
%
% $Id$

% ensure that the input is consistent with what this function expects
mri = ft_checkdata(mri, 'datatype', {'volume', 'segmentation'}, 'hasunit', 'yes');

% get the default options
cfg.tissue = ft_getopt(cfg, 'tissue');

if isempty(cfg.tissue)
  mri = ft_datatype_segmentation(mri, 'segmentationstyle', 'indexed');
  fn = fieldnames(mri);
  for i = 1:numel(fn)
    if (numel(mri.(fn{i})) == prod(mri.dim)) && (~strcmp(fn{i}, 'inside'))
      segfield = fn{i};
    end
  end
  cfg.tissue = setdiff(unique(mri.(segfield)(:)), 0);
end

if ischar(cfg.tissue)
  % it should either be something like {'brain', 'skull', 'scalp'}, or something like [1 2 3]
  cfg.tissue = {cfg.tissue};
end

if iscell(cfg.tissue)
  % the code below assumes that it is a probabilistic representation
  if any(strcmp(cfg.tissue, 'brain'))
    mri = ft_datatype_segmentation(mri, 'segmentationstyle', 'probabilistic', 'hasbrain', 'yes');
  else
    mri = ft_datatype_segmentation(mri, 'segmentationstyle', 'probabilistic');
  end
  % combine all tissue types
  seg = false(mri.dim);
  for i=1:numel(cfg.tissue)
    seg = seg | mri.(cfg.tissue{i});
  end
else
  % the code below assumes that it is an indexed representation
  mri = ft_datatype_segmentation(mri, 'segmentationstyle', 'indexed');
  % combine all tissue types
  seg = (mri.seg>0);
end

% this requires the external iso2mesh toolbox
ft_hastoolbox('iso2mesh', 1);

[node, elem, face] = vol2mesh(seg, 1:mri.dim(1), 1:mri.dim(2), 1:mri.dim(3), 2, 2, 0,'cgalsurf');

mesh = keepfields(mri, {'coordsys', 'unit'});
mesh.pos = ft_warp_apply(mri.transform, node);
mesh.tet = elem(:,1:4);

author: spm


text:
function db_update(LatestDbVersion, sProtocol, saveMetadata)
% DB_UPDATE: Updates any existing Brainstorm database to the current database structure.
%
% USAGE:  db_update(LatestDbVersion);              % Updates all the protocols of the database loaded in GlobalData
%         db_update(LatestDbVersion, iProtocols);  % Updates specific protocols
%         db_update(LatestDbVersion, sProtocol);   % Updates a new protocol
%
% INPUT:
%     - LatestDbVersion : The latest version number of the database structure
%     - sProtocol : Protocol structure, with the following fields (some fields are ignored, depending on the action)
%          |- Comment  : Name of the protocol
%          |- SUBJECTS : Directory that contains the anatomies of the subjects (MRI + surfaces)
%          |- STUDIES  : Directory that contains the functional data (recordings, sensors, sources...)
%          |- iStudy   : Ignored
%          |- UseDefaultAnat    : Ignored
%          |- UseDefaultChannel : Ignored
%     - saveMetadata : Whether to save the modifications of the protocol metadata (default: yes)
%
% @=============================================================================
% This function is part of the Brainstorm software:
% https://neuroimage.usc.edu/brainstorm
%
% Copyright (c) University of Southern California & McGill University
% This software is distributed under the terms of the GNU General Public License
% as published by the Free Software Foundation. Further details on the GPLv3
% 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.
%
% For more information type "brainstorm license" at command prompt.
% =============================================================================@
%
% Authors: Francois Tadel, 2009-2013; Martin Cousineau, 2019

global GlobalData;

%% Parse inputs
if nargin < 2 || isempty(sProtocol)
    iProtocols = 1:length(GlobalData.DataBase.ProtocolInfo);
    sProtocol = [];
elseif isnumeric(sProtocol)
    iProtocols = sProtocol;
    sProtocol = [];
end
if nargin < 3 || isempty(saveMetadata)
    saveMetadata = 1;
end
if isempty(sProtocol)
    % We are updating the GlobalData database
    isGlobal = 1;
    ProtocolsListInfos    = GlobalData.DataBase.ProtocolInfo(iProtocols);
    ProtocolsListStudies  = GlobalData.DataBase.ProtocolStudies(iProtocols);
    ProtocolsListSubjects = GlobalData.DataBase.ProtocolSubjects(iProtocols);
    CurrentDbVersion      = GlobalData.DataBase.DbVersion;
else
    % We are updating a specific protocol
    isGlobal = 0;
    [ProtocolMat, ProtocolFile] = GetProtocolMat(sProtocol);
    if isempty(ProtocolMat)
        % If not enough information to update protocol, return
        return;
    end
    ProtocolsListInfos    = sProtocol;
    ProtocolsListStudies  = ProtocolMat.ProtocolStudies;
    ProtocolsListSubjects = ProtocolMat.ProtocolSubjects;
    CurrentDbVersion      = ProtocolMat.DbVersion;
end

%% Check whether the update is necessary
if abs(CurrentDbVersion - LatestDbVersion) < 1e-8
    % Already up to date, return
    return;
elseif CurrentDbVersion > LatestDbVersion
    % The database is more up to date than the software, prompt for BST update
    if isGlobal
        bst_splash('hide');
    end
    errMsg = 'You are trying to load a protocol that was created with a newer version of Brainstorm.';
    [res, isCancel] = java_dialog('question', [errMsg 10 'Would you like to update Brainstorm now to continue?'], 'Update database');
    if ~isCancel && strcmpi(res, 'yes')
        bst_update();
        return;
    else
        error([errMsg 10 'Please update Brainstorm to continue.']);
    end
end

isFirstWarning = 1;

%% ===== UPDATE 31-Jul-2009 =====
% Modification: Add 'Image' field to the studies
function sStudy = AddImage(sStudy)
    sStudy.Image = repmat(db_template('Image'), 0);
end
% Perform update
if ~isempty(ProtocolsListStudies) && ~isempty(ProtocolsListStudies(1).Study) && ~isfield(ProtocolsListStudies(1).Study, 'Image')
    disp('BST> Database structure: Adding image support...');
    ApplyToDatabase(@AddImage);
    disp('BST> Database structure: Done.');
end

%% ===== UPDATE 16-Oct-2009 =====
% Modification: Add 'NoiseCov' field to the studies
function sStudy = AddNoiseCov(sStudy)
    sStudy.NoiseCov = repmat(db_template('NoiseCov'), 0);
end
% Perform update
if ~isempty(ProtocolsListStudies) && ~isempty(ProtocolsListStudies(1).Study) && ~isfield(ProtocolsListStudies(1).Study, 'NoiseCov')
    disp('BST> Database structure: Adding noise covariance matrix...');
    ApplyToDatabase(@AddNoiseCov);
    disp('BST> Database structure: Done.');
end

%% ===== UPDATE 20-Mar-2010 =====
% Modification: Add 'Dipoles' field to the studies
function sStudy = AddDipoles(sStudy)
    sStudy.Dipoles = repmat(db_template('Dipoles'), 0);
end
% Perform update
if ~isempty(ProtocolsListStudies) && ~isempty(ProtocolsListStudies(1).Study) && ~isfield(ProtocolsListStudies(1).Study, 'Dipoles')
    disp('BST> Database structure: Adding dipoles...');
    ApplyToDatabase(@AddDipoles);
    disp('BST> Database structure: Done.');
end

%% ===== UPDATE 15-Apr-2010 =====
% Modification: Add 'Timefreq' field to the studies
function sStudy = AddTimefreq(sStudy)
    sStudy.Timefreq = repmat(db_template('Timefreq'), 0);
end
% Perform update
if ~isempty(ProtocolsListStudies) && ~isempty(ProtocolsListStudies(1).Study) && ~isfield(ProtocolsListStudies(1).Study, 'Timefreq')
    disp('BST> Database structure: Adding time-frequency...');
    ApplyToDatabase(@AddTimefreq);
    disp('BST> Database structure: Done.');
end

%% ===== UPDATE 25-Jun-2010 =====
% Modifications: - Rename 'HeadModelName' field in 'Comment'
%                - Remove GridName field
function sStudy = UpdateHeadmodel(sStudy)
    % If no head models
    if isempty(sStudy.HeadModel)
        sStudy.HeadModel = repmat(db_template('headmodel'), 0);
    else
        if isfield(sStudy.HeadModel, 'HeadModelName')
            % For each headmodel in this study
            for iHeadModel = 1:length(sStudy.HeadModel)
                % Copy HeadModelName to Comment
                sStudy.HeadModel(iHeadModel).Comment = sStudy.HeadModel(iHeadModel).HeadModelName;
            end
            % Remove HeadModelName field
            sStudy.HeadModel = rmfield(sStudy.HeadModel, 'HeadModelName');
        end
        % Remove GridName field
        if isfield(sStudy.HeadModel, 'GridName')
            sStudy.HeadModel = rmfield(sStudy.HeadModel, 'GridName');
        end
    end
end
% If update was not performed yet
if ~isempty(ProtocolsListStudies) && ~isempty(ProtocolsListStudies(1).DefaultStudy) && isfield(ProtocolsListStudies(1).DefaultStudy.HeadModel, 'GridName')
    disp('BST> Database structure: Updating headmodel definition...');
    ApplyToDatabase(@UpdateHeadmodel);
    disp('BST> Database structure: Done.');
end

%% ===== UPDATE 25-Jun-2010 =====
% Modification: Add 'HeadModelType' to the Results structure
function sStudy = AddHeadModelType(sStudy)
    if ~isempty(sStudy.Result)
        [sStudy.Result.HeadModelType] = deal('surface');
    else
        sStudy.Result = repmat(db_template('results'), 0);
    end
end
% If update was not performed yet
if ~isempty(ProtocolsListStudies) && ~isempty(ProtocolsListStudies(1).AnalysisStudy) && ~isfield(ProtocolsListStudies(1).AnalysisStudy.Result, 'HeadModelType')
    disp('BST> Database structure: Adding HeadModelType to the Results structures...');
    ApplyToDatabase(@AddHeadModelType);
    disp('BST> Database structure: Done.');
end

%% ===== UPDATE 12-Jul-2010 =====
% Modification: Add 'BadTrial' and 'DataType' to the Data structure
% Function to update structure
function sStudy = AddBadTrial(sStudy)
    if ~isempty(sStudy.Data)
        [sStudy.Data.DataType] = deal('recordings');
        [sStudy.Data.BadTrial] = deal(0);
    else
        sStudy.Data = repmat(db_template('data'), 0);
    end
end
% If update was not performed yet
if ~isempty(ProtocolsListStudies) && ~isempty(ProtocolsListStudies(1).AnalysisStudy) && ~isfield(ProtocolsListStudies(1).AnalysisStudy.Data, 'DataType')
    disp('BST> Database structure: Adding BadTrial definition...');
    ApplyToDatabase(@AddBadTrial);
    disp('BST> Database structure: Done.');
end

%% ===== UPDATE 14-Jul-2010 =====
% Modification: Remove the "linkresults" files, now links are represented by filenames of the type "link|resultsfile.mat|datafile.mat"
if (CurrentDbVersion < 3)
    % Hide splash screen
    bst_splash('hide');
    % A full reload of the database is necessary to update all the links
    java_dialog('msgbox', ['Brainstorm was updated: the database structure changed.' 10 10 ...
                           'Please click on "Ok" and wait while all the database is being reloaded.' 10 ...
                           'Be patient, it might take several minutes...' 10 10], ...
                           'Database update');
    disp('BST> Reloading database... This might take several minutes...');
    % Reload database
    db_reload_database();
    % Get again the list of protocols for successive updates
    SaveProtocolStudies();
end

%% ===== UPDATE 07-October-2010 =====
% Modification: Add 'Matrix' field to the studies
function sStudy = AddMatrix(sStudy)
    sStudy.Matrix = repmat(db_template('Matrix'), 0);
end
% Perform update
if ~isempty(ProtocolsListStudies) && ~isempty(ProtocolsListStudies(1).Study) && ~isfield(ProtocolsListStudies(1).Study, 'Matrix')
    disp('BST> Database structure: Adding custom "matrix" file type...');
    ApplyToDatabase(@AddMatrix);
    disp('BST> Database structure: Done.');
end

%% ===== UPDATE 03-March-2011 =====
% Modification: Setting the subject file of "inter-subject" data to the default anatomy
if (CurrentDbVersion < 3.2)
    disp('BST> Database update: Setting the inter-subject nodes anatomy to the default anatomy.');
    % Set it protocol by protocol (if not set)
    for i = 1:length(ProtocolsListStudies)
        defSubjFile = bst_fullfile(bst_get('DirDefaultSubject'), 'brainstormsubject.mat');
        defSubjFileFull = bst_fullfile(ProtocolsListInfos(i).SUBJECTS, defSubjFile);
        if ~isempty(ProtocolsListStudies(i).AnalysisStudy) && isempty(ProtocolsListStudies(i).AnalysisStudy.BrainStormSubject) && file_exist(defSubjFileFull)
            ProtocolsListStudies(i).AnalysisStudy.BrainStormSubject = file_win2unix(defSubjFile);
        end
    end
    % Save database
    SaveProtocolStudies();
end

%% ===== UPDATE 15-March-2011 =====
% Database update: Changing the time reference in RAW FIF files from relative to absolute (the first sample is no longer t=0).
function sStudy = UpdateFifTimeRef(sStudy)
    for iData = 1:length(sStudy.Data)
        if strcmpi(sStudy.Data(iData).DataType, 'raw')
            try
                DataFile = file_fullpath(sStudy.Data(iData).FileName);
                disp(['UPDATE> Checking file: ' DataFile]);
                % Load file header
                DataMat = in_bst_data(DataFile);
                % If it is a RAW FIF file: Update time definition
                if strcmpi(DataMat.F.format, 'FIF') && isfield(DataMat.F.header,'raw') && (DataMat.F.header.raw.first_samp ~= 0) && (DataMat.Time(1) == 0)
                    disp('UPDATE> Updating file...');
                    addSamples = double(DataMat.F.header.raw.first_samp);
                    addTime    = addSamples / DataMat.F.prop.sfreq;
                    DataMat.Time         = DataMat.Time + addTime;
                    DataMat.F.prop.times = DataMat.F.prop.times + addTime;
                    for iEvt = 1:length(DataMat.F.events)
                        DataMat.F.events(iEvt).times = DataMat.F.events(iEvt).times + addTime;
                    end
                    % Save file back
                    bst_save(DataFile, DataMat, 'v6');
                end
            catch
                disp('UPDATE> Error updating file...');
            end
        end
    end
end
if (CurrentDbVersion < 3.3)
    disp('BST> Database update: Changing the time reference in RAW FIF files from relative to absolute (the first sample is no longer t=0).');
    ApplyToDatabase(@UpdateFifTimeRef);
    disp('BST> Database update: Done.');
end

%% ===== UPDATE 07-Dec-2011 =====
% Modification: Subject name is now defined based on the folder name (previously: "Name" folder)
if (CurrentDbVersion < 3.4)
    disp('BST> Database update: Forcing the subjects names to match the folder name.');
    % Update protocol by protocol
    for i = 1:length(ProtocolsListSubjects)
        % Process each subject
        for iSubj = 1:length(ProtocolsListSubjects(i).Subject)
            % Update name based on folder name
            ProtocolsListSubjects(i).Subject(iSubj).Name = bst_fileparts(ProtocolsListSubjects(i).Subject(iSubj).FileName);
        end
    end
    % Save database
    SaveProtocolSubjects();
end

%% ===== UPDATE 31-Jan-2013 =====
% Modified lots of things in the pipeline editor: resetting the options
if (CurrentDbVersion < 3.51)
    disp('BST> Software update: Resetting the process options...');
    % Reset all the saved options
    bst_set('ProcessOptions', []);
end

%% ===== UPDATE 21-Fev-2013 =====
% Function to update structure
function sStudy = AddIntraElectrodes(sStudy)
    if ~isempty(sStudy.HeadModel)
        [sStudy.HeadModel.ECOGMethod] = deal('');
        [sStudy.HeadModel.SEEGMethod] = deal('');
    else
        sStudy.Data = repmat(db_template('headmodel'), 0);
    end
end
% Modification: Add 'SEEGMethod' and 'ECOGMethod' to the Headmodel structure
if (CurrentDbVersion < 3.6) && ~isempty(ProtocolsListStudies) && ~isfield(ProtocolsListStudies(1).AnalysisStudy.HeadModel, 'ECOGMethod')
    disp('BST> Database structure: Adding intra-cranial electrodes definition...');
    ApplyToDatabase(@AddIntraElectrodes);
    disp('BST> Database structure: Done.');
    % Reset all the saved options
    bst_set('ProcessOptions', []);
end

%% ===== UPDATES: 7-May-2019, 22-May-2019, 31-Aug-2019 =====
% Modification: add fibers and FEM objects
if (CurrentDbVersion < 5.02)
    disp('BST> Database update: Adding support for fibers and FEM...');
    for iProt = 1:length(ProtocolsListSubjects)
        subjFields = fieldnames(ProtocolsListSubjects(iProt));
        for iField = 1:length(subjFields)
            subjField = subjFields{iField};
            nSubjects = length(ProtocolsListSubjects(iProt).(subjField));
            sSubjects = repmat(db_template('subject'), 1, nSubjects);
            for iSubj = 1:nSubjects
                % Add iFibers to subject
                sSubjects(iSubj) = struct_copy_fields(db_template('subject'), ProtocolsListSubjects(iProt).(subjField)(iSubj));
                % Add Fibers to subjectmat
                SubjectFile = bst_fullfile(ProtocolsListInfos(iProt).SUBJECTS, sSubjects(iSubj).FileName);
                if file_exist(SubjectFile)
                    SubjectMat = load(SubjectFile);
                    SubjectMat = struct_copy_fields(db_template('subjectmat'), SubjectMat);
                    bst_save(SubjectFile, SubjectMat, 'v7');
                end
            end
            ProtocolsListSubjects(iProt).(subjField) = sSubjects;
        end
    end
    disp('BST> Database update: Done.');
end

%% ===== JUST BEFORE RETURNING TO STARTUP FUNCTION =====
% Save the new database version
if saveMetadata
    if isGlobal
        GlobalData.DataBase.ProtocolInfo(iProtocols)     = ProtocolsListInfos;
        GlobalData.DataBase.ProtocolSubjects(iProtocols) = ProtocolsListSubjects;
        GlobalData.DataBase.ProtocolStudies(iProtocols)  = ProtocolsListStudies;
        GlobalData.DataBase.DbVersion        = LatestDbVersion;
        db_save();
    else
        ProtocolMat.ProtocolInfo     = ProtocolsListInfos;
        ProtocolMat.ProtocolSubjects = ProtocolsListSubjects;
        ProtocolMat.ProtocolStudies  = ProtocolsListStudies;
        ProtocolMat.DbVersion        = LatestDbVersion;
        bst_save(ProtocolFile, ProtocolMat, 'v7');
    end
end

%% ===================================================================================================
%  ===== HELPER FUNCTIONS ============================================================================
%  ===================================================================================================
%% ===== APPLY FUNCTION TO DATABASE =====
    function ApplyToDatabase(UpdateFcn)
        % ===== WARNING/BACKUP =====
        % Display warning before updating database
        if isFirstWarning
            % Hide splash screen
            bst_splash('hide');
            % Display warning message
            msgWarning = ['If after the update you cannot start Brainstorm anymore: ' 10 ...
                          '  - Do not try to run a previous version of the software, it would not help.' 10 ...
                          '  - Type "brainstorm reset", start Brainstorm, then re-import your protocols.' 10 ...
                          '  - Report the bug on the user forum'];
            disp(['*********************************************************************' 10 ...
                  msgWarning 10 ...
                  '*********************************************************************']);
            isOk = java_dialog('confirm', ...
                ['Warning: You are running a newer version of Brainstorm.' 10 ...
                 'Some modifications have to be performed on your database.' 10 10 ...
                 msgWarning  10 10 ...
                 'Update database now?' 10 10], 'Database update');
            if ~isOk
                error(['User cancelled database update.' 10 'This version of Brainstorm cannot be started, please try your previous version.']);
            end
            % Make a backup copy of the current database
            c = clock;
            DbFile = bst_get('BrainstormDbFile');
            BakDbFile = file_unique(strrep(DbFile, '.mat', sprintf('_%02.0f%02.0f%02.0f_%02.0f%02.0f.bak', c(1)-2000, c(2:5))));
            try
                file_copy(DbFile, BakDbFile);
            catch
                disp(['UPDATE> Cannot write backup file "' BakDbFile '".']);
            end
            % Do not display the warning again
            isFirstWarning = 0;
        end

        % ===== PROCESSING ======
        % Progress bar
        bst_progress('start', 'Update database', 'Updating database structure...');
        % For each protocol
        for ip = 1:length(ProtocolsListStudies)
            % If no study in the protocol
            if isempty(ProtocolsListStudies(ip).Study)
                %ProtocolsListStudies(ip).Study = repmat(db_template('Study'), 0);
                ProtocolsListStudies(ip).Study(1).Comment = 'update';
                ProtocolsListStudies(ip).Study = repmat(UpdateFcn(ProtocolsListStudies(ip).Study), 0);
            else
                % Update each study in this protocol
                ProtocolStudies = repmat(ProtocolsListStudies(ip).Study, 0);
                for is = 1:length(ProtocolsListStudies(ip).Study)
                    ProtocolStudies(is) = UpdateFcn(ProtocolsListStudies(ip).Study(is));
                end
                ProtocolsListStudies(ip).Study = ProtocolStudies;
            end
            if ~isempty(ProtocolsListStudies(ip).DefaultStudy)
                ProtocolsListStudies(ip).DefaultStudy = UpdateFcn(ProtocolsListStudies(ip).DefaultStudy);
            end
            if ~isempty(ProtocolsListStudies(ip).AnalysisStudy)
                ProtocolsListStudies(ip).AnalysisStudy = UpdateFcn(ProtocolsListStudies(ip).AnalysisStudy);
            end
        end
        % Save database
        SaveProtocolStudies();
        bst_progress('stop');
    end

    %% ===== SAVE CHANGES TO DATABASE =====
    function SaveProtocolStudies()
        if saveMetadata
            if isGlobal
                GlobalData.DataBase.ProtocolStudies(iProtocols) = ProtocolsListStudies;
            else
                ProtocolMat.ProtocolStudies = ProtocolsListStudies;
                bst_save(ProtocolFile, ProtocolMat, 'v7');
            end
        end
    end
    function SaveProtocolSubjects()
        if saveMetadata
            if isGlobal
                GlobalData.DataBase.ProtocolSubjects(iProtocols) = ProtocolsListSubjects;
            else
                ProtocolMat.ProtocolSubjects = ProtocolsListSubjects;
                bst_save(ProtocolFile, ProtocolMat, 'v7');
            end
        end
    end

end

function [ProtocolMat, ProtocolFile] = GetProtocolMat(sProtocol)
    ProtocolMat  = [];
    ProtocolFile = [];

    if isstruct(sProtocol) && isfield(sProtocol, 'STUDIES') && ~isempty(sProtocol.STUDIES)
        ProtocolFile = bst_fullfile(sProtocol.STUDIES, 'protocol.mat');
        if file_exist(ProtocolFile)
            ProtMat = load(ProtocolFile);
            if isfield(ProtMat, 'DbVersion')
                ProtocolMat = ProtMat;
            end
        end
    end
end

author: brainstorm-tools


text:
function [suspicious_index lof] = LOF(A, k)
%
% Local Outlier Factor
% Authors: Markus M. Breunig, Hans-Peter Kriegel,
%          Raymond T. Ng, J?rg Sander
% Original paper :
% LOF: Identifying Density-Based Local Outliers
% e-mail : { breunig | kriegel | sander }
%          @dbs.informatik.uni-muenchen.de
%          rng@cs.ubc.ca
% Programmer: Yi-Ren Yeh(yirenyeh@gmail.com)
% modified by: Zi-Wen Gui(evan176@hotmail.com)
%
%
% Inputs
%   A: the data matrix, each row represents an instance
%   k: the number of nearest neighbors, specified as an integer or
%      as a fraction of the total number of data points
%
% Outputs
%   lof: the local outlier factor for each instance
%   suspicious_index: the ranking of instances according to their
%                     suspicious score
%                     For example, suspicious_index(i)=j means the
%                     ith instance is in jth position in the ranking
%

if k < 1
    [numrows ~] = size(A);
    k = round(k*numrows);
end

try
    %Find the nearest neighbors by "KDTree" for each elements
    [k_index, k_dist] = knnsearch(A,A,'k',k+1,'nsmethod','kdtree','IncludeTies',true);
    %Ignore first element(itself) at nearest neighbors
    k_index = cellfun(@(x) x(2:end),k_index,'UniformOutput',false);
    numneigh = cellfun('length',k_index);
    %Get k-distance
    k_dist1 = cell2mat(cellfun(@(x) x(end),k_dist,'UniformOutput',false));
    %Get row length of matrix A
    n = length(A(:,1));
    %Initialize lrd_value vector
    lrd_value = zeros(n,1);
    %Calculate lrd for each elements
    for i = 1:n
        lrd_value(i) = lrd(A, i, k_dist1, k_index, numneigh(i));
    end
    %Initialize lof vector
    lof = zeros(n,1);
    %Calculate LOF
    for i = 1:n
        lof(i) = sum(lrd_value(k_index{i})/lrd_value(i))/numneigh(i);
    end
    %Indices from sorting lof are the suspicious score rankings
    [~,suspicious_index]=sort(lof,'descend');

catch err
    if (strcmp(err.message, 'Invalid parameter name: IncludeTies.'))
        warning('MATLAB:LOF', 'Matlab not newest version? Falling back to old version.')
        [suspicious_index lof] = LOF_old(A, k);
    else
        rethrow(err)
    end

end

%=========================================================================
function lrd_value = lrd(A, index_p, k_dist,k_index, numneighbors)
%Calculate the reachability distance for nearest neighbors
Temp = repmat(A(index_p,:), numneighbors, 1) - A(k_index{index_p}, :);
Temp = sqrt(sum(Temp.^2,2));
reach_dist = max([Temp k_dist(k_index{index_p})],[],2);
%Calculate the local reachability density for each elements
lrd_value = numneighbors/sum(reach_dist);

author: dsmi-lab-ntust

repo: AnomalyDetectionToolbox


text:
function [ h, hmin, hmax ] = p16_stepsize ( option )

%*****************************************************************************80
%
%% P16_STEPSIZE returns step sizes for problem 16.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    24 September 2008
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, integer OPTION, the option index.
%
%    Output, real H, HMIN, HMAX, suggested values for the
%    initial step, the minimum step, and the maximum step.
%
  h =    0.200;
  hmin = 0.001;
  hmax = 2.000;

  return
end