text:
function [index,distance]=near(x,x0,n);
% NEAR  finds the indices of x that are closest to the point x0.
%function [index,distance]=near(x,x0,[n]);
%     x is an array, x0 is a point, n is the number of closest points to get
%     (in order of increasing distance).  Distance is the abs(x-x0)
% rsignell@usgs.gov

if nargin==2
     n=1;
end
[distance,index]=sort(abs(x-x0));
distance=distance(1:n);
index=index(1:n);

end

author: nctoolbox

repo: nctoolbox

save_path: github-repos/MATLAB/nctoolbox-nctoolbox


text:
function data=locsmooth(data,Fs,Tw,Ts)
%  Running line fit (using local linear regression) - 1d only, continuous
%  processes
%  Usage: data=locsmooth(data,Fs,Tw,Ts)
%  Inputs:
% Note that units of Fs, movinwin have to be consistent.
%  data  (single vector)
%  Fs    (sampling frequency) - optional. Default 1
%  Tw    (length of moving window) - optional.  Default. full length of data (global detrend)
%  Ts    (step size) - optional. Default Tw/2.
%
% Output:
% data   (locally smoothed data).
data=change_row_to_column(data);
N=size(data,1);
if nargin < 2; Fs=1; end;
if nargin < 3; Tw=N/Fs; end;
if nargin < 4; Ts=Tw/2; end;

n=round(Fs*Tw);
dn=round(Fs*Ts);
if ~isreal(data)
  yr=real(data);
  yi=imag(data);
  tmp=runline(yr,n,dn);
  yr=tmp';
  tmp=runline(yi,n,dn);
  yi=tmp';
  data=yr+i*yi;
else
  tmp=runline(data,n,dn);
  data=tmp';
end

author: buzsakilab

repo: buzcode

save_path: github-repos/MATLAB/buzsakilab-buzcode


text:
function [y,fs]=readflac(filename,mode)
%READWAV  Read a .FLAC format sound file [Y,FS]=(FILENAME,MODE)
%
% Input Parameters:
%
%	FILENAME gives the name of the file (with optional .WAV extension) or alternatively
%                 can be the FIDX output from a previous call to READWAV
%	MODE		specifies the following (*=default):
%
%    Scaling: 's'    Auto scale to make data peak = +-1
%             'r'    Raw unscaled data (integer values)
%             'q'    Scaled to make 0dBm0 be unity mean square
%             'p' *	 Scaled to make +-1 equal full scale
%             'o'    Scale to bin centre rather than bin edge (e.g. 127 rather than 127.5 for 8 bit values)
%                     (can be combined with n+p,r,s modes)
%             'n'    Scale to negative peak rather than positive peak (e.g. 128.5 rather than 127.5 for 8 bit values)
%                     (can be combined with o+p,r,s modes)
%
% FLAC (Free lossless audio codec) is a compressed audio file format described here:
% http://flac.sourceforge.net/

%      Copyright (C) Mike Brookes 2008
%      Version: $Id: readflac.m 713 2011-10-16 14:45:43Z dmb $
%
%   VOICEBOX is a MATLAB toolbox for speech processing.
%   Home page: http://www.ee.ic.ac.uk/hp/staff/dmb/voicebox/voicebox.html
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%   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 2 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 can obtain a copy of the GNU General Public License from
%   http://www.gnu.org/copyleft/gpl.html or by writing to
%   Free Software Foundation, Inc.,675 Mass Ave, Cambridge, MA 02139, USA.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if nargin<2
    mode='p';
else
    mode = [mode(:).' 'p'];
end
dirt=voicebox('dir_temp');
[fnp,fnn,fne,fnv]=fileparts(filename);
filetemp=fullfile(dirt,[fnn '.wav' fnv]);
doscom=['"' voicebox('flac') '"' ' -d -f -o "' filetemp '" "' filename '"'];
%                     fprintf(1,'Executing: %s\n',doscom);
[doss,dosr]=dos(doscom); % run the program
if doss % test for errors
    error(sprintf('Error running DOS command: %s',doscom));
end
if exist(filetemp)~=2
    error(sprintf('No output file from: %s',doscom));
end
[y,fs]=readwav(filetemp,mode);
doscom=['del /f "' filetemp '"'];
if dos(doscom) % run the program
    error(sprintf('Error running DOS command: %s',doscom));
end

author: covarep

repo: covarep

save_path: github-repos/MATLAB/covarep-covarep


text:
% pop_prepPipeline() - runs the early stage pipeline to reference and to
% detect bad channels
%
% Usage:
%   >>   [OUTEEG, com] = pop_prepPipeline(INEEG, params);
%
% Inputs:
%   INEEG   - input EEG dataset
%   params  - (optional) structure with parameters to override defaults
%
% Outputs:
%   OUTEEG  - output dataset
%
% See also:
%   prepPipeline, prepPipelineReport, EEGLAB

% Copyright (C) 2015  Kay Robbins with contributions from Nima
% Bigdely-Shamlo, Christian Kothe, Tim Mullen, and Cassidy Matousek
%
% 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 2 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, write to the Free Software
% Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA

function [EEG, com] = pop_prepPipeline(EEG, params)
com = ''; % Return something if user presses the cancel button
okay = true;
if nargin < 1  %% display help if not enough arguments
    help pop_prepPipeline;
    return;
elseif nargin < 2
    params = struct();
end

%% Add path to prepPipeline subdirectories if not in the list
tmp = which('getPrepDefaults');
if isempty(tmp)
    myPath = fileparts(which('prepPipeline'));
    addpath(genpath(myPath));
end

%% Pop up window
if nargin < 2
    userData = getUserData();
    [params, okay] = MasterGUI([],[],userData, EEG);
end

%% Return if user pressed cancel
if (~okay)
	return;
end

%% Being the pipeline execution
userData = getUserData();
com = sprintf('%s = pop_prepPipeline(%s, %s);', inputname(1), ...
    struct2str(params));

reportMode = userData.report.reportMode.value;
consoleFID = userData.report.consoleFID.value;
publishOn = userData.report.publishOn.value;
summaryFilePath = userData.report.summaryFilePath.value;
sessionFilePath = userData.report.sessionFilePath.value;

if okay
    if strcmpi(reportMode, 'normal') || strcmpi(reportMode, 'skipReport')
        EEG = prepPipeline(EEG, params);
    end
    if (strcmpi(reportMode, 'normal') || ...
            strcmpi(reportMode, 'reportOnly')) %&& publishOn
        publishPrepReport(EEG, summaryFilePath, sessionFilePath, ...
            consoleFID, publishOn);
    end
    if strcmpi(reportMode, 'normal') || strcmpi(reportMode, 'skipReport')
        EEG = prepPostProcess(EEG, params);
    end
end

    function userData = getUserData()
        %% Gets the userData defaults and merges it with the parameters
        userData = struct('boundary', [], 'detrend', [], ...
            'lineNoise', [], 'reference', [], ...
            'report', [],  'postProcess', []);
        stepNames = fieldnames(userData);
        for k = 1:length(stepNames)
            defaults = getPrepDefaults(EEG, stepNames{k});
            [theseValues, errors] = checkStructureDefaults(params, ...
                defaults);
            if ~isempty(errors)
                error('pop_prepPipeline:BadParameters', ['|' ...
                    sprintf('%s|', errors{:})]);
            end
            userData.(stepNames{k}) = theseValues;
        end
    end  % getUserData

end % pop_prepPipeline

author: VisLab

repo: EEG-Clean-Tools


text:
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% setting GPU & Enviroment etc.
addpath('caffe/matlab');
caffe.reset_all();

caffe.set_mode_gpu();
gpu_id = 1;
caffe.set_device(gpu_id);
pp=0;

[foldername, videoname, filename, labelid] = textread('lst/UCF-101_testSplit1_25frames.lst','%s %s %d %d');

net_model = ['prototxt/spatial_cls.prototxt'];
net_weights = ['models/spatial_UCF101_list01.caffemodel'];
phase = 'test';
if ~exist(net_weights, 'file')
  error('Please download CaffeNet from Model Zoo before you run this demo');
end

net = caffe.Net(net_model, net_weights, phase);
im_data_ori = zeros(340,256, 3,25,'single');
mean_data = caffe.io.read_mean('mean/VGG_mean.binaryproto');

input_data = zeros(224,224,3,25);
scores_total=zeros(101,size(foldername,1)/25);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Processing videos
for i=1:25:size(filename)
    ii=1;
    t_io = 0;
    t1 = clock;
    for j=0:24
	for k=filename(i+j):filename(i+j)
	  i_str=sprintf('%04d',k);
      t2 = clock;
      im = caffe.io.load_image(['ucf101_image/',char(foldername{i+j}),'/',char(videoname{i+j}),'/flow_i_',i_str,'.jpg']);
      t3 = clock;
      t_io = etime(t3,t2) + t_io;

	  im_data_ori(:,:,:,j+1) = im(:,:,:);
%      im_data_mirror(:,:,:,j+1) = im(end:-1:1,:,:);
	end
    end

	input_data = im_data_ori(60:60+223,16:16+223,:,:);
	input_data = bsxfun(@minus,input_data,mean_data);

    net.blobs('data').set_data(input_data);
    net.forward_prefilled();
    scores = net.blobs('prob').get_data();
    scores_total(:,(i-1)/25+1) = mean(scores, 2);
    t4 = clock;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    total_time = etime(t4,t1);
    io_time = t_io;
    computational_time = etime(t4,t1)-t_io;
    fprintf(['Time (w/o I/O) of processing 25 rgb frames: ',num2str(computational_time),'\n'])

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%show results:
    [~, maxlabel] = max(mean(scores,2));
    if (maxlabel==labelid(i)+1)
        pp=pp+1;
    end
    fprintf('%d/%d %d %d %f\n',pp,(i-1)/25+1,maxlabel, labelid(i)+1, pp/((i-1)/25+1));
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Cleaning and etc.
save('UCF101_spatial_list01.mat','scores_total');
caffe.reset_all();
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

author: zbwglory

repo: MV-release

save_path: github-repos/MATLAB/zbwglory-MV-release


text:
classdef ROIPooling < dagnn.Layer
  % DAGNN.ROIPOOLING  Region of interest pooling layer

  % Copyright (C) 2016 Hakan Bilen.
  % All rights reserved.
  %
  % This file is part of the VLFeat library and is made available under
  % the terms of the BSD license (see the COPYING file).

  properties
    method = 'max'
    subdivisions = [6 6]
    transform = 1
    flatten = false
  end

  methods
    function outputs = forward(obj, inputs, params)
      numROIs = numel(inputs{2}) / 5 ;
      outputs{1} = vl_nnroipool(...
        inputs{1}, inputs{2}, ...
        'subdivisions', obj.subdivisions, ...
        'transform', obj.transform, ...
        'method', obj.method) ;
      if obj.flatten
        outputs{1} = reshape(outputs{1},1,1,[],numROIs) ;
      end
    end

    function [derInputs, derParams] = backward(obj, inputs, params, derOutputs)
      numROIs = numel(inputs{2}) / 5 ;
      if obj.flatten
        % unflatten
        derOutputs{1} = reshape(...
          derOutputs{1},obj.subdivisions(1),obj.subdivisions(2),[],numROIs) ;
      end
      derInputs{1} = vl_nnroipool(...
        inputs{1}, inputs{2}, derOutputs{1}, ...
        'subdivisions', obj.subdivisions, ...
        'transform', obj.transform, ...
        'method', obj.method) ;
      derInputs{2} = [];
      derParams = {} ;
    end

    function outputSizes = getOutputSizes(obj, inputSizes)
      if isempty(inputSizes{1})
        n = 0 ;
      else
        n = prod(inputSizes{2})/5 ;
      end
      outputSizes{1} = [obj.subdivisions, inputSizes{1}(3), n] ;
    end

    function obj = ROIPooling(varargin)
      obj.load(varargin) ;
    end
  end
end

author: guosheng

repo: refinenet

save_path: github-repos/MATLAB/guosheng-refinenet


text:
function handles = togglePlotSilhouette(hObject, handles, state)

set(handles.btnPlotSilhouette, 'Enable', state);

if isequal(handles.foundExportFig, 1)
    set(handles.btnPlotSilhouetteSave, 'Enable', state);
end

guidata(hObject, handles);

author: beckel

repo: nilm-eval

save_path: github-repos/MATLAB/beckel-nilm-eval


text:
function xlsPasteTo(filename,sheetname,width, height,varargin)

%Paste current figure to selected Excel sheet and cell
%
%
% xlsPasteTo(filename,sheetname,width, height,range)
%Example:
%xlsPasteTo('File.xls','Sheet1',200, 200,'A1')
% this will paset into A1 at Sheet1 at File.xls the current figure with
% width and height of 200
%
% tal.shir@hotmail.com

options = varargin;

    range = varargin{1};

[fpath,file,ext] = fileparts(char(filename));
if isempty(fpath)
    fpath = pwd;
end
Excel = actxserver('Excel.Application');
set(Excel,'Visible',0);
Workbook = invoke(Excel.Workbooks, 'open', [fpath filesep file ext]);

sheet = get(Excel.Worksheets, 'Item',sheetname);
invoke(sheet,'Activate');

    ExAct = Excel.Activesheet;

  ExActRange = get(ExAct,'Range',range);
    ExActRange.Select;
    pos=get(gcf,'Position');
    set(gcf,'Position',[ pos(1:2) width height])
    print -dmeta

invoke(Excel.Selection,'PasteSpecial');

invoke(Workbook, 'Save');
invoke(Excel, 'Quit');
delete(Excel);

author: Sable

repo: mcbench-benchmarks


text:
%% Copyright (C) 2016-2017, 2019 Colin B. Macdonald
%%
%% This file is part of OctSymPy.
%%
%% OctSymPy 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 software 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 software; see the file COPYING.
%% If not, see <http://www.gnu.org/licenses/>.

%% -*- texinfo -*-
%% @documentencoding UTF-8
%% @defmethod  @@sym ezsurf (@var{z})
%% @defmethodx @@sym ezsurf (@var{f1}, @var{f2}, @var{f3})
%% @defmethodx @@sym ezsurf (@dots{}, @var{dom})
%% @defmethodx @@sym ezsurf (@dots{}, @var{N})
%% Simple 3D surface plots of symbolic expressions.
%%
%% Example 3D surface plot:
%% @example
%% @group
%% syms x y
%% z = sin(2*x)*sin(y)
%%   @result{} z = (sym) sin(2⋅x)⋅sin(y)
%% ezsurf(z)                                    % doctest: +SKIP
%% @end group
%% @end example
%%
%% Example parametric surface plot of a Möbius strip:
%% @example
%% @group
%% syms u v
%% x = (1+v*cos(u/2))*cos(u)
%%   @result{} x = (sym)
%%       ⎛     ⎛u⎞    ⎞
%%       ⎜v⋅cos⎜─⎟ + 1⎟⋅cos(u)
%%       ⎝     ⎝2⎠    ⎠
%% y = (1+v*cos(u/2))*sin(u);
%% z = v*sin(u/2);
%%
%% ezsurf(x, y, z, [0 2*pi -0.5 0.5], 32)       % doctest: +SKIP
%% axis equal
%% @end group
%% @end example
%%
%% See help for the (non-symbolic) @code{ezsurf}, which this
%% routine calls after trying to convert sym inputs to
%% anonymous functions.
%%
%% @seealso{ezsurf, @@sym/ezmesh, @@sym/ezplot, @@sym/function_handle}
%% @end defmethod

function varargout = ezsurf(varargin)

  % first input is handle, shift
  if (ishandle(varargin{1}))
    firstpotsym = 2;
  else
    firstpotsym = 1;
  end

  maxnumsym = 3;
  firstsym = [];

  for i = firstpotsym:nargin
    if (isa(varargin{i}, 'sym'))
      if (i < firstpotsym + maxnumsym)
        % one of the fcns to plot, covert to handle fcn

        % Each is function of one var, and its the same var for all
        thissym = symvar(varargin{i});
        assert(length(thissym) <= 2, ...
          'ezsurf: parameterized: functions should have at most two inputs');
        if (isempty(thissym))
          % a number, create a constant function in a dummy variable
          % (0*t works around some Octave oddity on 3.8 and hg Dec 2014)
          thisf = inline(sprintf('%g + 0*t', double(varargin{i})), 't');
          %thisf = @(t) 0*t + double(varargin{i});  % no
        else
          % check variables match (sanity check)
          if (isempty(firstsym))
            firstsym = thissym;
          else
            assert(all(logical(thissym == firstsym)), ...
              'ezsurf: all functions must be in terms of the same variables');
          end
          thisf = function_handle(varargin{i});
        end

        varargin{i} = thisf;

      else
        % plot ranges, etc, convert syms to doubles
        varargin{i} = double(varargin{i});
      end
    end
  end

  h = ezsurf(varargin{:});

  if (nargout)
    varargout{1} = h;
  end

end

%!error <all functions must be in terms of the same variables>
%! syms u v t
%! ezsurf(u*v, 2*u*v, 3*v*t)

%!error <functions should have at most two inputs>
%! syms u v t
%! ezsurf(u*v, 2*u*v, u*v*t)

repo: octsympy


text:
function [ xc, rank, seed ] = comp_random_grlex ( kc, rank1, rank2, seed )

%*****************************************************************************80
%
%% COMP_RANDOM_GRLEX randomly selects a composition in a given rank range.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    06 September 2014
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, integer KC, the number of parts of the composition.
%    1 <= KC.
%
%    Input, integer RANK1, RANK2, the lowest and highest ranks to consider.
%    1 <= RANK1 <= RANK2.
%
%    Input/output, integer SEED, a seed for the random number generator.
%
%    Output, integer XC(KC), a randomly chosen composition with rank
%    between RANK1 and RANK2.
%
%    Output, integer RANK, the rank of the composition.
%

%
%  Ensure that 1 <= KC.
%
  if ( kc < 1 )
    fprintf ( 1, '\n' );
    fprintf ( 1, 'COMP_RANDOM_GRLEX - Fatal error!\n' );
    fprintf ( 1, '  KC < 1\n' );
    error ( 'COMP_RANDOM_GRLEX - Fatal error!' );
  end
%
%  Ensure that 1 <= RANK1.
%
  if ( rank1 < 1 )
    fprintf ( 1, '\n' );
    fprintf ( 1, 'COMP_RANDOM_GRLEX - Fatal error!\n' );
    fprintf ( 1, '  RANK1 < 1\n' );
    error ( 'COMP_RANDOM_GRLEX - Fatal error!' );
  end
%
%  Ensure that RANK1 <= RANK2.
%
  if ( rank2 < rank1 )
    fprintf ( 1, '\n' );
    fprintf ( 1, 'COMP_RANDOM_GRLEX - Fatal error!\n' );
    fprintf ( 1, '  RANK2 < RANK1\n' );
    error ( 'COMP_RANDOM_GRLEX - Fatal error!' );
  end
%
%  Select the rank.
%
  [ rank, seed ] = i4_uniform_ab ( rank1, rank2, seed );
%
%  Recover the composition.
%
  xc = comp_unrank_grlex ( kc, rank );

  return
end

author: johannesgerer

repo: jburkardt-m


text:
% SHOW - Displays an image with the right size and colors and with a title.
%
% Usage:
%         h = show(im)
%         h = show(im, figNo)
%         h = show(im, title)
%         h = show(im, figNo, title)
%
% Arguments:  im    - Either a 2 or 3D array of pixel values or the name
%                     of an image file;
%             figNo - Optional figure number to display image in. If
%                     figNo is 0 the current figure or subplot is
%                     assumed.
%             title - Optional string specifying figure title
%
% Returns:    h     - Handle to the figure.  This allows you to set
%                     additional figure attributes if desired.
%
% The function displays the image, automatically setting the colour map to
% grey if it is a 2D image, or leaving it as colour otherwise, and setting
% the axes to be 'equal'.  The image is also displayed as 'TrueSize', that
% is, pixels on the screen match pixels in the image (if it is possible
% to fit it on the screen, otherwise MATLAB rescales it to fit).
%
% Unless you are doing a subplot (figNo==0) the window is sized to match
% the image, leaving no border, and hence saving desktop real estate.
%
% If figNo is omitted a new figure window is created for the image.  If
% figNo is supplied, and the figure exists, the existing window is reused to
% display the image, otherwise a new window is created. If figNo is 0 the
% current figure or subplot is assumed.
%
% See also: SHOWSURF

% Copyright (c) 2000-2009 Peter Kovesi
% School of Computer Science & Software Engineering
% The University of Western Australia
% http://www.csse.uwa.edu.au/
%
% Permission is hereby granted, free of charge, to any person obtaining a copy
% of this software and associated documentation files (the "Software"), to deal
% in the Software without restriction, subject to the following conditions:
%
% The above copyright notice and this permission notice shall be included in
% all copies or substantial portions of the Software.
%
% The Software is provided "as is", without warranty of any kind.

% October   2000  Original version
% March     2003  Mods to alow figure name in window bar and allow for subplots.
% April     2007  Proper recording and restoring of MATLAB warning state.
% September 2008  Octave compatible
% May       2009  Reworked argument handling logic for extra flexibility

function h = show(im, param2, param3)

    Octave = exist('OCTAVE_VERSION') ~= 0;  % Are we running under Octave?

    if ~Octave
	s = warning('query','all'); % Record existing warning state.
	warning('off');             % Turn off warnings that might arise if image
                                    % has to be rescaled to fit on screen
    end

    % Check case where im is an image filename rather than image data
    if ~isnumeric(im) & ~islogical(im)
        Title = im;            % Default title is file name
        im = imread(im);
    else
        Title = inputname(1);  % Default title is variable name of image data
    end

    figNo = -1;                % Default value indicating create new figure

    % If two arguments check type of 2nd argument to see if it is the title or
    % figure number that has been supplied
    if nargin == 2
        if strcmp(class(param2),'char')
            Title = param2;
        elseif isnumeric(param2) && length(param2) == 1
            figNo = param2;
        else
            error('2nd argument must be a figure number or title');
        end
    elseif nargin == 3
        figNo = param2;
        Title = param3;

        if ~strcmp(class(Title),'char')
            error('Title must be a string');
        end

        if ~isnumeric(param2) || length(param2) ~= 1
            error('Figure number must be an integer');
        end
    end

    if figNo > 0           % We have a valid figure number
        figure(figNo);     % Reuse or create a figure window with this number
        if ~Octave
            subplot('position',[0 0 1 1]); % Use the whole window
        end
    elseif figNo == -1
        figNo = figure;        % Create new figure window
        if ~Octave
            subplot('position',[0 0 1 1]); % Use the whole window
        end
    end

    if ndims(im) == 2          % Display as greyscale
	imagesc(im);
	colormap('gray');
    else
	imshow(im);            % Display as RGB
    end

    if figNo == 0                 % Assume we are trying to do a subplot
        figNo = gcf;              % Get the current figure number
	axis('image'); axis('off');
	title(Title);             % Use a title rather than rename the figure
    else
	axis('image'); axis('off');
	set(figNo,'name', ['  ' Title])
	if ~Octave
	    truesize(figNo);
	end
    end

    if nargout == 1
       h = figNo;
    end

    if ~Octave
	warning(s);  % Restore warnings
    end

author: jianxiongxiao

repo: ProfXkit


text:
function [resdQ, rhsQ] = volumeAveragedConcentrationFlux(dQ,Q,jflux,T,param)
%   volumeAveragedConcentrationFlux is used to implement the three parameters reduced model for solid phase diffusion.

%   This file is part of the LIONSIMBA Toolbox
%
%	Official web-site: 	http://sisdin.unipv.it/labsisdin/lionsimba.php
% 	Official GitHUB: 	https://github.com/lionsimbatoolbox/LIONSIMBA
%
%   LIONSIMBA: A Matlab framework based on a finite volume model suitable for Li-ion battery design, simulation, and control
%   Copyright (C) 2016-2018 :Marcello Torchio, Lalo Magni, Davide Raimondo,
%                            University of Pavia, 27100, Pavia, Italy
%                            Bhushan Gopaluni, Univ. of British Columbia,
%                            Vancouver, BC V6T 1Z3, Canada
%                            Richard D. Braatz,
%                            Massachusetts Institute of Technology,
%                            Cambridge, Massachusetts 02142, USA
%
%   Main code contributors to LIONSIMBA 2.0:
%                           Ian Campbell, Krishnakumar Gopalakrishnan,
%                           Imperial college London, London, UK
%
%   LIONSIMBA is a free Matlab-based software distributed with an MIT
%   license.

% Diffusion coefficients for the solid phase
[Dps_eff, Dns_eff] = param.SolidDiffusionCoefficientsFunction(T,param);

rhsQ_p      = -30*Dps_eff./param.Rp_p^2.*Q(1:param.Np)    -45./(2*param.Rp_p^2).*jflux(1:param.Np);

resdQ_p     = dQ(1:param.Np)        -rhsQ_p;

rhsQ_n      = -30*Dns_eff./param.Rp_n^2.*Q(param.Np+1:end)-45./(2*param.Rp_n^2).*jflux(param.Np+1:end);

resdQ_n     = dQ(param.Np+1:end) - rhsQ_n   ;

rhsQ        = [rhsQ_p;rhsQ_n];
resdQ       = [resdQ_p;resdQ_n];

% This model has been taken from the paper, "Efficient Macro-Micro Scale Coupled
% Modeling of Batteries" -  Subramanian,Diwakar,Tapriyal - 2005 JES

end

author: lionsimbatoolbox

repo: LIONSIMBA


text:
function e = end(T,i,n)
% overloaded end function

if n==1
  e = prod(size(T)); %#ok<PSIZE>
else
  e = size(T,i);
end

author: mtex-toolbox

repo: mtex

save_path: github-repos/MATLAB/mtex-toolbox-mtex


text:
function [ obj ] = vb_invert( obj )
% Run VB update equations until convergence.
%
% This is a protected method of the tapas_Huge class. It cannot be called
% from outside the class.
%

% Author: Yu Yao (yao@biomed.ee.ethz.ch)
% Copyright (C) 2019 Translational Neuromodeling Unit
%                    Institute for Biomedical Engineering,
%                    University of Zurich and ETH Zurich.
%
% This file is part of TAPAS, which is released under the terms of the GNU
% General Public Licence (GPL), version 3. For further details, see
% <https://www.gnu.org/licenses/>.
%
% This software is provided "as is", without warranty of any kind, express
% or implied, including, but not limited to the warranties of
% merchantability, fitness for a particular purpose and non-infringement.
%
% This software is intended for research only. Do not use for clinical
% purpose. Please note that this toolbox is under active development.
% Considerable changes may occur in future releases. For support please
% refer to:
% https://github.com/translationalneuromodeling/tapas/issues
%

%% initialization
obj = obj.vb_init( );
bUpdateKmeans       = false;
bUpdateAssignments  = false;
bUpdateClusters     = false;
bVerbose            = obj.options.nvp.verbose;
nIt = obj.options.nvp.numberofiterations;
if isempty(nIt)
    nIt = 999;
end
obj.trace.nfe = repmat(obj.trace.nfe, nIt, 1);

%% ======= MAIN LOOP =======
for iIt = 1:nIt
% schedule posterior updates

    if bUpdateKmeans
        % initialize cluster parameters with kmeans
        obj = kmeans_init( obj );
        bUpdateKmeans = false;
    end

    if bUpdateAssignments
        % update alpha, q_nk
        obj = update_assigments( obj );
    end

    if bUpdateClusters
        if obj.options.confVar
            % update m_beta, S_beta
            obj = update_confounds( obj );
        end
        % update m_k, tau_k, nu_k, S_k
        obj = update_clusters( obj ); %%% TODO add support for cluster-specific prior mean and scale
    end

    % update mu_n, Sigma_n, b_nr
    obj = update_dcm( obj, bUpdateClusters );

    % update negative free energy
    obj.trace.nfe(iIt) = obj.vb_nfe( );
    dNfe = obj.trace.nfe(iIt) - obj.posterior.nfe;
    obj.posterior.nfe = obj.trace.nfe(iIt);
    if bVerbose
        if isnan(dNfe)
            fprintf('iteration %3u.\n', iIt);
        else
            fprintf('iteration %3u, dF: %.2E.\n', iIt, dNfe);
        end
    end

    % check stopping conditions
    if ~bUpdateClusters
        % start updating clustering parameters and confound coefficients
        bUpdateClusters = (dNfe < obj.options.convergence.dDcm && dNfe >= 0);% || iIt > 32;
        if bUpdateClusters
            obj.trace.convergence(1) = iIt;
        end

    elseif ~bUpdateAssignments
        % start updating assignment labels
        bUpdateAssignments = dNfe < obj.options.convergence.dClusters;
        if bUpdateAssignments
            bUpdateKmeans = obj.options.convergence.bUseKmeans;
            obj.trace.convergence(2) = iIt;
            % force next iteration to happen
            obj.posterior.nfe = NaN;
        end

    elseif dNfe < obj.options.convergence.dVb
        % convergence of VB scheme
        obj.trace.convergence(3) = iIt;
        obj.trace.nfe(iIt+1:end) = [];
        obj.trace.epsilon = obj.aux.epsilon;
        if bVerbose
            fprintf('VB converged after %u iterations. ', iIt)
        end
        break;

    end
end

if bVerbose
    if iIt == nIt
        fprintf('VB reached the maximum number of iterations. ')
    end
    fprintf('Negative free energy: %.2E\n', obj.posterior.nfe);
end

%% post-processing
% estimate confound-free DCM parameters
if obj.options.confVar
    mu_beta = obj.aux.mu_beta;
    Sigma_beta = obj.aux.Sigma_beta;
    if obj.options.confVar > 1
        mu_beta = sum(bsxfun(@times, mu_beta, permute(obj.posterior.q_nk, ...
            [1 3 2])), 3);
        Sigma_beta = sum(bsxfun(@times, Sigma_beta, permute( ...
            obj.posterior.q_nk, [3 4 1 2])), 4);
    end
    obj.posterior.mu_r = obj.posterior.mu_n(:, 1:obj.idx.P_c) - mu_beta;
    obj.posterior.Sigma_r = obj.posterior.Sigma_n(1:obj.idx.P_c, 1:obj.idx.P_c, :) ...
        + Sigma_beta;
else
    obj.posterior.mu_r = [];
    obj.posterior.Sigma_r = [];
end

% normalized variance of residual (by subject and region)
obj.posterior.nrv = zeros(obj.N, obj.R);
for n = 1:obj.N
    varRes = var(obj.aux.epsilon{n});
    varSig = var(obj.data(n).bold);
    obj.posterior.nrv(n,:) = varRes./varSig;
end

end

%% UPDATE: posterior over cluster parameters
function [ obj ] = update_clusters( obj )

% update parameters of inverse Wishart
obj.posterior.tau = obj.prior.tau_0 + obj.aux.q_k;
obj.posterior.nu  = obj.prior.nu_0  + obj.aux.q_k;

% auxiliary variables
Sigma_n_c = obj.posterior.Sigma_n(1:obj.idx.P_c,1:obj.idx.P_c,:);
mu_n_c = obj.posterior.mu_n(:,1:obj.idx.P_c);
if obj.options.confVar == 1
    Sigma_n_c = Sigma_n_c + obj.aux.Sigma_beta;
    mu_n_c = mu_n_c - obj.aux.mu_beta;
end

for k = 1:obj.K

    if obj.options.confVar > 1
        Sigma_n_c = obj.posterior.Sigma_n(1:obj.idx.P_c,1:obj.idx.P_c,:) + ...
            obj.aux.Sigma_beta(:,:,:,k);
        mu_n_c = obj.posterior.mu_n(:,1:obj.idx.P_c) - ...
            obj.aux.mu_beta(:,:,k);
    end
    % Eq (17)
    mu_k_c = obj.posterior.q_nk(:,k)'*mu_n_c;
    Sigma_k_c = sum(bsxfun(@times, Sigma_n_c, ...
        permute(obj.posterior.q_nk(:,k),[3,2,1])), 3);
    % Eq (16)
    % posterior cluster mean
    obj.posterior.m(k,:) = (mu_k_c + obj.prior.tau_0*obj.prior.m_0)/...
        (obj.aux.q_k(k) + obj.prior.tau_0);
    % posterior cluster covariance
    mu_k_c = mu_k_c./obj.aux.q_k(k);
    dmu = bsxfun(@minus, mu_n_c, mu_k_c);
    obj.posterior.S(:,:,k) = ...
        Sigma_k_c + ...
        dmu'*bsxfun(@times, dmu, obj.posterior.q_nk(:,k)) + ...
        obj.aux.q_k(k)*obj.prior.tau_0/...
            (obj.aux.q_k(k) + obj.prior.tau_0)*...
            (mu_k_c - obj.prior.m_0)'*(mu_k_c - obj.prior.m_0) + ...
        obj.prior.S_0;
    % update auxiliary variables
    obj.aux.ldS(k) = tapas_huge_logdet(obj.posterior.S(:,:,k));
    obj.aux.nu_inv_S(:,:,k) = obj.posterior.nu(k).*inv(obj.posterior.S(:,:,k));

end
end

%% UPDATE: posterior over cluster assignments
function [ obj ] = update_assigments( obj )

Sigma_n_c = obj.posterior.Sigma_n(1:obj.idx.P_c,1:obj.idx.P_c,:);
mu_n_c = obj.posterior.mu_n(:,1:obj.idx.P_c);
if obj.options.confVar == 1
    Sigma_n_c = Sigma_n_c + obj.aux.Sigma_beta;
    mu_n_c = mu_n_c - obj.aux.mu_beta;
end
psiNu = sum(psi(0, bsxfun(@minus, obj.posterior.nu, 0:obj.idx.P_c-1)/2), 2);
psiAlpha = psi(0, obj.posterior.alpha);
log_q_nk = zeros(obj.N, obj.K);
% Eq (18)
for k = 1:obj.K

    if obj.options.confVar > 1
        Sigma_n_c = obj.posterior.Sigma_n(1:obj.idx.P_c,1:obj.idx.P_c,:) + ...
            obj.aux.Sigma_beta(:,:,:,k);
        mu_n_c = obj.posterior.mu_n(:,1:obj.idx.P_c) - ...
            obj.aux.mu_beta(:,:,k);
    end
    dm = bsxfun(@minus,mu_n_c,obj.posterior.m(k,:));
    log_q_nk(:,k) = -.5*sum((dm*obj.aux.nu_inv_S(:,:,k)).*dm, 2) ...
        -.5*squeeze(sum(sum(bsxfun(@times, Sigma_n_c, ...
            obj.aux.nu_inv_S(:,:,k)), 1), 2));

end
log_q_nk = bsxfun(@plus,log_q_nk, (-.5*obj.aux.ldS +.5*psiNu ...
    -.5*obj.idx.P_c./obj.posterior.tau + psiAlpha)');
log_q_nk = exp(bsxfun(@minus, log_q_nk, max(log_q_nk, [], 2)));
obj.posterior.q_nk = bsxfun(@rdivide, log_q_nk, sum(log_q_nk, 2));

% Eq (17)
obj.aux.q_k = sum(obj.posterior.q_nk,1)' + realmin;

% Eq (15)
obj.posterior.alpha = obj.prior.alpha_0 + obj.aux.q_k;

end

%% UPDATE: posterior over confound coefficients
function [ obj ] = update_confounds( obj )

m_k = permute(obj.posterior.m, [3,2,1]);
mu_n_c = obj.posterior.mu_n(:, 1:obj.idx.P_c);

for p = 1:obj.idx.P_c

    % Eq (XX)
    tmp = obj.aux.nu_inv_S(p, p, :);
    tmp = bsxfun(@times, obj.posterior.q_nk, tmp(:)');
    if obj.options.confVar == 1
        tmp = sum(tmp, 2);
    end
    % posterior precision over beta
    obj.posterior.S_beta(:, :, p, :) = bsxfun(@plus, obj.aux.Pi_beta_0, ...
        sum(bsxfun(@times, obj.aux.x_n_2, reshape(tmp, 1, 1, obj.N, [])), 3));

    % Eq (XX+1)
    obj.aux.mu_beta(:,p,:) = 0;
    dmu = bsxfun(@minus, mu_n_c, obj.aux.mu_beta);
    dmu = bsxfun(@minus, dmu, m_k);
    tmp = permute(obj.posterior.q_nk, [1 3 2]).*...
        sum(bsxfun(@times, obj.aux.nu_inv_S(p, :, :), dmu), 2);
    if obj.options.confVar == 1
        tmp = sum(tmp, 3);
    end
    tmp = permute(sum(bsxfun(@times, tmp, obj.aux.x_n), 1), [2 1 3]);
    for k = 1:size(obj.posterior.m_beta, 3)
        obj.posterior.m_beta(:,p,k) = obj.posterior.S_beta(:,:,p,k)\...
            (obj.aux.Pi_m_beta_0 + tmp(:,1,k));
        % posterior covariance over beta
        obj.posterior.S_beta(:,:,p,k) = inv(obj.posterior.S_beta(:,:,p,k));
        obj.aux.ldS_beta(p,k) = tapas_huge_logdet( ...
            obj.posterior.S_beta(:,:,p,k));
    end
    % update auxiliary variables
    obj.aux.mu_beta(:,p,:) = sum(bsxfun(@times, obj.aux.x_n, ...
        permute(obj.posterior.m_beta(:,p,:), [2 1 3])), 2);
    obj.aux.Sigma_beta(p,p,:,:) = sum(sum(bsxfun(@times, obj.aux.x_n_2, ...
        obj.posterior.S_beta(:,:,p,:)), 1), 2);

end

end

%% UPDATE: posterior over DCM and noise parameters
function [ obj ] = update_dcm( obj, bCheckNfe )

% cache DCM parameters
tmpMu       = obj.posterior.mu_n;
tmpSigma    = obj.posterior.Sigma_n;
tmpLdSigma  = obj.aux.ldSigma;
tmpEpsilon  = obj.aux.epsilon;
tmpG        = obj.aux.G;
tmpB        = obj.posterior.b;
tmpBp       = obj.aux.b_prime;
tmpLb       = obj.aux.lambda_bar;
tmpNfe      = obj.vb_nfe( );

% subject-wise update DCM and noise parameters
for n = 1:obj.N
    varEps = sum(var(obj.aux.epsilon{n}));
    % Eq (20)
    m_k = obj.posterior.m;
    if obj.options.confVar
        tmp = permute(obj.aux.mu_beta(n,:,:), [3 2 1]);
        m_k = bsxfun(@plus, m_k, tmp);
    end
    lambda_prime_c = zeros(obj.idx.P_c);
    mu_prime_c = zeros(1,obj.idx.P_c);
    for k = 1:obj.K
        tmps = obj.posterior.q_nk(n,k)*obj.aux.nu_inv_S(:,:,k);
        lambda_prime_c = lambda_prime_c + tmps;
        mu_prime_c = mu_prime_c + m_k(k,:)*tmps;
    end

    % Eq (19)
    % posterior covariance
    Lambda_bar = repmat(obj.aux.lambda_bar(n,:), obj.aux.q_r(n), 1);
    G_Lambda = bsxfun(@times, obj.aux.G{n}, Lambda_bar(:))';
    Pi_n = G_Lambda*obj.aux.G{n};
    Pi_n(1:obj.idx.P_c, 1:obj.idx.P_c) = lambda_prime_c + ...
        Pi_n(1:obj.idx.P_c, 1:obj.idx.P_c);
    Pi_n(obj.idx.P_c+1:end, obj.idx.P_c+1:end) = obj.aux.Pi_h + ...
        Pi_n(obj.idx.P_c+1:end, obj.idx.P_c+1:end);
    obj.posterior.Sigma_n(:,:,n) = inv(Pi_n);
    obj.aux.ldSigma(n) = - tapas_huge_logdet(Pi_n);
    % posterior mean
    tmp = obj.aux.epsilon{n}(:) + obj.aux.G{n}*obj.posterior.mu_n(n,:)';
    % regularization factor
    tmp = tmp*(.9.^mod(obj.trace.nRetract(n),5));
    obj.posterior.mu_n(n,:) = ...
        ((G_Lambda*tmp)' + [mu_prime_c, obj.aux.mu_prime_h])/Pi_n;
    % recalculate jacobian and residual
    obj = obj.options.fncBold(obj, n);

    % Eq (22)
    tmp = sum((obj.aux.G{n}*obj.posterior.Sigma_n(:,:,n)).*obj.aux.G{n},2);
    obj.aux.b_prime(n,:) = sum(obj.aux.epsilon{n}.^2 + ...
        reshape(tmp, obj.aux.q_r(n), obj.R), 1);
    % Eq (21)
    obj.posterior.b(n,:) = obj.prior.b_0 + obj.aux.b_prime(n,:)/2;
    % Eq (23)
    obj.aux.lambda_bar(n,:) = obj.posterior.a(n,:)./obj.posterior.b(n,:);

    newNfe = obj.vb_nfe( );
    % check negative free energy (during first few iterations, accept
    % decrease in F if fit improves)
    bUnstable = any(isnan(obj.aux.G{n}(:))) || any(isinf(obj.aux.G{n}(:)));
    bOverride = bCheckNfe || varEps < sum(var(obj.aux.epsilon{n}));
    if bUnstable || (tmpNfe > newNfe && bOverride)
        % retract update if new parameters unstable or if F decreased
        obj.posterior.mu_n(n,:)      = tmpMu(n,:);
        obj.posterior.Sigma_n(:,:,n) = tmpSigma(:,:,n);
        obj.aux.ldSigma(n)           = tmpLdSigma(n);
        obj.aux.epsilon{n}           = tmpEpsilon{n};
        obj.aux.G{n}                 = tmpG{n};
        obj.posterior.b(n,:)         = tmpB(n,:);
        obj.aux.b_prime(n,:)         = tmpBp(n,:);
        obj.aux.lambda_bar(n,:)      = tmpLb(n,:);
        obj.trace.nRetract(n)        = obj.trace.nRetract(n) + 1;
    else
        % otherwise accept update
        tmpNfe = newNfe;
        obj.trace.nDcmUpdate(n) = obj.trace.nDcmUpdate(n) + 1;
        obj.trace.nRetract(n)   = 0;
    end

end

end

%% initialize cluster parameters with K-means
function [ obj ] = kmeans_init( obj )
% kmeans
mu_n_c = obj.posterior.mu_n(:,1:obj.idx.P_c);
if obj.options.confVar
    mu_n_c = mu_n_c - obj.aux.mu_beta(:,:,1);
end
idx = kmeans(mu_n_c, obj.K, 'Replicates', obj.const.nKmeans);

% initialize additional clusters
obj.posterior.q_nk = zeros(obj.N,obj.K);
obj.posterior.q_nk(sub2ind([obj.N,obj.K], (1:obj.N)', idx)) = 1;
obj.aux.q_k = sum(obj.posterior.q_nk, 1)' + realmin;
obj.posterior.alpha = obj.prior.alpha_0 + obj.aux.q_k;
if obj.options.confVar > 1
    obj.posterior.m_beta = repmat(obj.posterior.m_beta(:,:,1), [1 1 obj.K]);
    obj.posterior.S_beta = repmat(obj.posterior.S_beta(:,:,:,1), ...
        [1 1 1 obj.K]);
    obj.aux.Sigma_beta = repmat(obj.aux.Sigma_beta(:,:,:,1), [1 1 1 obj.K]);
    obj.aux.ldS_beta = repmat(obj.aux.ldS_beta(:,1), 1, obj.K);
    obj.aux.mu_beta = repmat(obj.aux.mu_beta(:,:,1), [1 1 obj.K]);
end
obj = update_clusters( obj );

% debug info
obj.trace.kmeans.labels = idx;
obj.trace.kmeans.inputs = mu_n_c;

end

author: translationalneuromodeling

repo: tapas


text:
function MRSIStruct = setDimension(MRSIStruct, dimLabel, value)
    switch(lower(dimLabel))
        case 'x'
            MRSIStruct.dims.x = value;
        case 'y'
            MRSIStruct.dims.y = value;
        case 'z'
            MRSIStruct.dims.z = value;
        case 'coils'
            MRSIStruct.dims.coils = value;
        case 'averages'
            MRSIStruct.dims.averages = value;
        case 'kx'
            MRSIStruct.dims.kx = value;
        case 'ky'
            MRSIStruct.dims.ky = value;
        case 'kz'
            MRSIStruct.dims.kz = value;
        case 'subspecs'
            MRSIStruct.dims.subspecs = value;
        case 'extras'
            MRSIStruct.dims.extras = value;
        case 't'
            MRSIStruct.dims.t = value;
        otherwise
            error('Dimension not found %s', dimLabel);
    end
end

author: CIC-methods

repo: FID-A

save_path: github-repos/MATLAB/CIC-methods-FID-A


text:
function poly_fprintf ( file_name, node_num, segment, edge_num, ...
  edge_nodes, hole_num, hole_point )

%*****************************************************************************80
%
%% POLY_WRITE writes data to a POLY file.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    12 April 2009
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, string FILE_NAME, the name of the file to create.
%
%    Input, integer NODE_NUM, the number of nodes.
%
%    Input, real SEGMENT(2,NODE_NUM), the nodes.
%
%    Input, integer EDGE_NUM, the number of edges.
%
%    Input, integer EDGE_NODES(2,EDGE_NUM), the nodes that form each edge.
%
%    Input, integer HOLE_NUM, the number of holes in the mesh.
%
%    Input, real HOLE_POINT(2,HOLE_NUM), a point in each hole.
%
  region_num = 0;

  poly_unit = fopen ( file_name, 'w' );

  if ( poly_unit < 0 )
    fprintf ( 1, '\n' );
    fprintf ( 1, 'POLY_WRITE - Fatal error!\n' );
    fprintf ( 1, '  Could not open the output POLY file.\n' );
    error ( 'POLY_WRITE - Fatal error!' );
  end

  fprintf ( poly_unit, '#  %s\n', file_name );
  fprintf ( poly_unit, '#  Created by poly_write.m\n' );
  fprintf ( poly_unit, '#\n' );
  fprintf ( poly_unit, '#  Vertex  Dimension  Attribute  Marker\n' );
  fprintf ( poly_unit, '#  Count              Count      0/1\n' );
  fprintf ( poly_unit, '#\n' );
  fprintf ( poly_unit, '  %d  2  0  0\n', node_num );
  fprintf ( poly_unit, '#\n' );
  fprintf ( poly_unit, '#  Vertex  X  Y  Attributes  Marker\n' );
  fprintf ( poly_unit, '#  Index\n' );
  fprintf ( poly_unit, '#\n' );
  for node = 1 : node_num
    fprintf ( poly_unit, '  %8d  %10f  %10f\n', node, segment(1:2,node) );
  end
  fprintf ( poly_unit, '#\n' );
  fprintf ( poly_unit, '#  Segment  Marker\n' );
  fprintf ( poly_unit, '#  Count    0/1\n' );
  fprintf ( poly_unit, '#\n' );
  fprintf ( poly_unit, '  %d  0\n', edge_num );
  fprintf ( poly_unit, '#\n' );
  fprintf ( poly_unit, '#  Segment_index  Node1  Node2  Marker\n' );
  fprintf ( poly_unit, '#\n' );
  for edge = 1 : edge_num
    fprintf ( poly_unit, '  %8d  %8d  %8d  0\n', edge, edge_nodes(1:2,edge) );
  end
  fprintf ( poly_unit, '#\n' );
  fprintf ( poly_unit, '#  Hole\n' );
  fprintf ( poly_unit, '#  Count\n' );
  fprintf ( poly_unit, '#\n' );
  fprintf ( poly_unit, '  %d\n', hole_num );
  fprintf ( poly_unit, '#\n' );
  fprintf ( poly_unit, '#  Hole_index  X  Y\n' );
  fprintf ( poly_unit, '#\n' );
  for hole = 1 : hole_num
    fprintf ( poly_unit, '  %8d  %10f  %10f\n', hole, hole_point(1:2,hole) );
  end
  fprintf ( poly_unit, '#\n' );
  fprintf ( poly_unit, '#  Region\n' );
  fprintf ( poly_unit, '#  Count\n' );
  fprintf ( poly_unit, '#\n' );
  fprintf ( poly_unit, '  %d\n', region_num );

  fclose ( poly_unit );

  return
end


text:
classdef VademecumLoader < handle

    properties (SetAccess = private, GetAccess = public)
        dataBase
    end

    properties (Access = private)
        index
        path
        iIndex
        jIndex
        tensor
        microName
        tensorCase
        tensName
    end

    methods (Access = public)

        function obj = VademecumLoader(d)
            obj.init(d)
            obj.loadMxMy();
            obj.loadHomogenizedTensor()
            obj.loadAmplificatorTensor()
        end

    end

    methods (Access = private)

        function init(obj,d)
            obj.index = d.index;
            obj.path = d.path;
            obj.microName = d.microName;
        end

        function loadMxMy(obj)
            m1Path = fullfile(obj.path,'m1m2.txt');
            obj.dataBase.mx = load(m1Path);
            obj.dataBase.my = obj.dataBase.mx;
            obj.dataBase.nmx = length(obj.dataBase.mx);
            obj.dataBase.nmy = length(obj.dataBase.my);
        end

        function loadHomogenizedTensor(obj)
            obj.tensorCase = 'C_';
            obj.tensName = 'ConstitutiveTensor';
            obj.loadTensors();
            obj.dataBase.C = obj.tensor;
        end

        function loadAmplificatorTensor(obj)
            obj.tensorCase = 'P_';
            obj.tensName = 'AmplificatorTensor';
            obj.loadTensors();
            obj.dataBase.P = obj.tensor;
        end

        function loadTensors(obj)
            nPlot = length(obj.index);
            for iplot = 1:nPlot
                obj.obtainIJindex(iplot);
                obj.loadTensor();
            end
        end

       function obtainIJindex(obj,iplot)
            obj.iIndex = obj.index(iplot,1);
            obj.jIndex = obj.index(iplot,2);
       end

        function loadTensor(obj)
            mN    = obj.microName;
            tName = obj.tensName;
            i = obj.iIndex;
            j = obj.jIndex;
            name = [obj.tensorCase,num2str(i),num2str(j)];
            fullPath = fullfile(obj.path,mN,tName,[mN,name,'.txt']);
            obj.tensor(i,j,:,:) = load(fullPath);
        end

    end

end

author: SwanLab

repo: Swan

save_path: github-repos/MATLAB/SwanLab-Swan


text:
function varargout = DQ_Spectra_GUI(varargin)
% DQ_SPECTRA_GUI MATLAB code for DQ_Spectra_GUI.fig
%      DQ_SPECTRA_GUI, by itself, creates a new DQ_SPECTRA_GUI or raises the existing
%      singleton*.
%
%      H = DQ_SPECTRA_GUI returns the handle to a new DQ_SPECTRA_GUI or the handle to
%      the existing singleton*.
%
%      DQ_SPECTRA_GUI('CALLBACK',hObject,eventData,handles,...) calls the local
%      function named CALLBACK in DQ_SPECTRA_GUI.M with the given input arguments.
%
%      DQ_SPECTRA_GUI('Property','Value',...) creates a new DQ_SPECTRA_GUI or raises the
%      existing singleton*.  Starting from the left, property value pairs are
%      applied to the GUI before DQ_Spectra_GUI_OpeningFcn gets called.  An
%      unrecognized property name or invalid value makes property application
%      stop.  All inputs are passed to DQ_Spectra_GUI_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 pushbutton_help DQ_Spectra_GUI

% Last Modified by GUIDE v2.5 10-Feb-2023 16:27:47

% Begin initialization code - DO NOT EDIT
gui_Singleton = 0;
gui_State = struct('gui_Name',       mfilename, ...
    'gui_Singleton',  gui_Singleton, ...
    'gui_OpeningFcn', @DQ_Spectra_GUI_OpeningFcn, ...
    'gui_OutputFcn',  @DQ_Spectra_GUI_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 DQ_Spectra_GUI is made visible.
function DQ_Spectra_GUI_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 DQ_Spectra_GUI (see VARARGIN)

% Choose default command line output for DQ_Spectra_GUI
handles.output = [];

% check input EEG
EEG = varargin{1};
avg_fft = varargin{2};
fft_labels = varargin{3};
yout = varargin{4};
yout_lab = varargin{5};
chanArray = varargin{6};
%guiwin_num = varargin{4};
guiwin_num = 1;

%GUI positions (for multiple DQ tables)

if guiwin_num ~= 1
    offset_GUI = 25 * randi([1 10],1); %offset by
    % %set(hObject, 'Units', 'normalized');
    figureposition = get(handles.figure1, 'Position');
    set(hObject, 'Position', ...
        [ figureposition(1)-offset_GUI figureposition(2)-offset_GUI  ...
        figureposition(3) figureposition(4)])
end

%check contents

% Update GUI with ERPSET info
%n_erp = numel(ALLERP);

%for the case of continuous EEG
eeg_names = EEG.setname;
handles.text_EEGSET_title.String = 'Selected EEG Dataset:';
handles.text_EEGSET_title.FontSize = 12;
handles.active_eegset.Style = 'text';
if numel(eeg_names) > 15
   handles.active_eegset.FontSize = 12;
   set(handles.active_eegset,'FontUnits', 'normalized');
end
%update handles
guidata(hObject, handles);

%ERPSET_title_str = ['ERPSET - ' erp_names(current_ERP)];
handles.active_eegset.String = eeg_names;
handles.active_eegset.Value = 1;

%AMS: Multiple EEGset options
handles.newerpwin.String = eeg_names;
handles.newerpwin.Value = 1;

%n_dq = 2; %power or amplitude
typenames = {'Power', 'Amplitude'};
n_dq = length(typenames);

handles.popupmenu_DQ_type.String = typenames;
%Amplitude first by default
handles.popupmenu_DQ_type.Value = n_dq;
[n_elec, n_bands] = size(avg_fft);
selected_DQ_type = handles.popupmenu_DQ_type.Value;

if selected_DQ_type == 1
    table_data = (avg_fft)^2; %Power

else
    table_data = avg_fft;

end

handles.dq_table.Data = table_data;
handles.orig_data = table_data;

nchan = EEG.nbchan;
chanArray = chanArray(chanArray<=nchan);
handles.chanArray = chanArray;

% electrode labels from EEGset(rows)
for i=1:length(chanArray)
    elec_labels{i} = EEG.chanlocs(chanArray(i)).labels;
end
listch = elec_labels;
handles.dq_table.RowName = elec_labels;
handles.orig_RowName = elec_labels;

%selections
handles.sel_row = n_elec; % set to max on load
handles.sel_col = n_bands; % set to max on load
handles.listch = listch;
handles.indxlistch = chanArray;
%set(handles.edit_channels,'String', vect2colon(chanArray, 'Delimiter', 'off'));

% frequency band column labels
handles.dq_table.ColumnName = fft_labels;
handles.orig_ColName = fft_labels;

% Set font size of DQ Table
desired_fontsize = erpworkingmemory('fontsizeGUI');
handles.dq_table.FontSize = desired_fontsize;
handles.heatmap_on = 0;

% Set outliers text window off
set(handles.stdwindow,'Enable','Off');

%set outliers value in handles off initially
handles.outliers_on = 0;

% Update handles structure
handles.EEG = EEG;
handles.avg_fft = avg_fft;
handles.yout = yout;
handles.yout_lab = yout_lab;
%
% if n_erp == 0
%     %case of aSME preavg
%     handles.ALLERP = ERP;
% else
%     handles.ALLERP = ALLERP;
% end

guidata(hObject, handles);

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

% --- Outputs from this function are returned to the command line.
function varargout = DQ_Spectra_GUI_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
varargout{1} = [];

% --- Executes during object creation, after setting all properties.
function text_ERPSET_title_CreateFcn(hObject, eventdata, handles)
% hObject    handle to text_ERPSET_title (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    empty - handles not created until after all CreateFcns called

% --- Executes on selection change in popupmenu_DQ_type.
function popupmenu_DQ_type_Callback(hObject, eventdata, handles)
% hObject    handle to popupmenu_DQ_type (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

% Hints: contents = cellstr(get(hObject,'String')) returns popupmenu_DQ_type contents as cell array
%        contents{get(hObject,'Value')} returns selected item from popupmenu_DQ_type
% handles.popupmenu_DQ_type.Value saves the indx of the selected DQ Measure
selected_DQ_type = handles.popupmenu_DQ_type.Value;
avg_fft = handles.avg_fft;
%selected_bin = handles.popupmenu_bin.Value;

if strcmp(handles.popupmenu_DQ_type.String{selected_DQ_type},'Amplitude')
    table_data = avg_fft;
else
    table_data = (avg_fft).^2; %power
end

handles.dq_table.Data = table_data;
handles.orig_data = table_data;

if handles.heatmap_on
    redraw_heatmap(hObject, eventdata, handles);
end

if handles.outliers_on
   clear_heatmap(hObject, eventdata, handles);
   handles.heatmap_on = 0;

   redraw_outliers(hObject, eventdata, handles);
   set(handles.checkbox_outliers,'Value',1);
end

% Update handles structure
guidata(hObject, handles);

% --- Executes during object creation, after setting all properties.
function popupmenu_DQ_type_CreateFcn(hObject, eventdata, handles)
% hObject    handle to popupmenu_DQ_type (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    empty - handles not created until after all CreateFcns called

% Hint: popupmenu 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

% --- Executes on selection change in popupmenu_bin.
% function popupmenu_bin_Callback(hObject, eventdata, handles)
% hObject    handle to popupmenu_bin (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

% Hints: contents = cellstr(get(hObject,'String')) returns popupmenu_bin contents as cell array
%        contents{get(hObject,'Value')} returns selected item from popupmenu_bin
% handles.popupmenu_DQ_type.Value saves the indx of the selected DQ Measure
% selected_DQ_type = handles.popupmenu_DQ_type.Value;
% selected_bin = handles.popupmenu_bin.Value;
%
% % Check data exists, plot
% if isempty(handles.ERP.dataquality(selected_DQ_type).data)
%
%     % if pointwise SEM, use ERP.binerror
%     if strcmp(handles.ERP.dataquality(selected_DQ_type).type,'Point-wise SEM') && isempty(handles.ERP.binerror) == 0
%         table_data = handles.ERP.binerror(:,:,selected_bin);
%     else
%
%         % Data empty here.
%         table_data = nan(1);
%         disp('DQ data not found for this DQ type and bin. Perhaps it has been cleared?');
%     end
%
% else
%     % data is present and correct
%     table_data = handles.ERP.dataquality(selected_DQ_type).data(:,:,selected_bin);
% end
%
%
% handles.dq_table.Data = table_data;
% handles.orig_data = table_data;
%
%
% %if user switches bin and both heatmap & outliers on,
% %use shortcircuit to clear outliers and redraw heatmap only
% if handles.heatmap_on %|| handles.outliers_on
%     clear_outliers(hObject, eventdata, handles);
%     redraw_heatmap(hObject, eventdata, handles);
%     %turn outliers off
%     handles.outliers = 0;
%     set(handles.checkbox_outliers,'Value',0);
%
% end
%
% %if only outliers on, then clear heatmap and redraw outliers
% if handles.outliers_on
%    clear_heatmap(hObject, eventdata, handles);
%    handles.heatmap_on = 0;
%
%    redraw_outliers(hObject, eventdata, handles);
%    set(handles.checkbox_outliers,'Value',1);
% end
%
%
% % Update handles structure
% guidata(hObject, handles);

% --- Executes during object creation, after setting all properties.
function popupmenu_bin_CreateFcn(hObject, eventdata, handles)
% hObject    handle to popupmenu_bin (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    empty - handles not created until after all CreateFcns called

% Hint: popupmenu 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

% --- Executes on button press in checkbox_text_labels.
% function checkbox_text_labels_Callback(hObject, eventdata, handles)
% % hObject    handle to checkbox_text_labels (see GCBO)
% % eventdata  reserved - to be defined in a future version of MATLAB
% % handles    structure with handles and user data (see GUIDATA)
%
% % Hint: get(hObject,'Value') returns toggle state of checkbox_text_labels
% %selected_DQ_type = handles.popupmenu_DQ_type.Value;
% %selected_bin = handles.popupmenu_bin.Value;
%
% % if strcmp(handles.ERP.dataquality(selected_DQ_type).type,'SD Across Trials') && isempty(handles.ERP.binerror) == 0
% %     table_data = handles.ERP.dataquality(selected_DQ_type).data.SD_bias(:,:,selected_bin);
% % else
% %     table_data = handles.ERP.dataquality(selected_DQ_type).data(:,:,selected_bin);
% %
% % end
% % handles.dq_table.Data = table_data;
% % handles.orig_data = table_data;
%
% % Time-window labels
% clear tw_labels
% ERP = handles.ERP;
% n_dq = numel(ERP.dataquality);
% %[n_elec, n_tw, n_bin] = size(ERP.dataquality(selected_DQ_type).data);
% if strcmp(handles.ERP.dataquality(selected_DQ_type).type,'SD Across Trials') && isempty(handles.ERP.binerror) == 0
%     [n_elec, n_tw, n_bin] = size(ERP.dataquality(selected_DQ_type).data.SD_bias);
% else
%     [n_elec, n_tw, n_bin] = size(ERP.dataquality(selected_DQ_type).data);
% end
% if isfield(ERP.dataquality(selected_DQ_type),'time_window_labels') && isempty(ERP.dataquality(selected_DQ_type).time_window_labels) == 0 && handles.checkbox_text_labels.Value == 1
%     tw_labels = ERP.dataquality(selected_DQ_type).time_window_labels;
% elseif isempty(ERP.dataquality(selected_DQ_type).times)
%     tw_labels = [];
% elseif n_tw == 0  % some data problem, no tw here
%     tw_labels{1} = 'No data here - measure cleared?';
% else
%     for i=1:n_tw
%         tw_labels{i} = [num2str(ERP.dataquality(selected_DQ_type).times(i,1)) ' : ' num2str(ERP.dataquality(selected_DQ_type).times(i,2))];
%     end
% end
% handles.dq_table.ColumnName = tw_labels;
%
% if handles.heatmap_on
%     redraw_heatmap(hObject, eventdata, handles);
% end
%
%
%  % Update handles structure
%  guidata(hObject, handles);

% --- Executes on button press in pushbutton_xls.
function pushbutton_xls_Callback(hObject, eventdata, handles)
% hObject    handle to pushbutton_xls (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)
empty_filename = [];
selected_DQ_type = handles.popupmenu_DQ_type.Value;
avg_fft = handles.avg_fft; %only the average across frequencies within band
chans_per_row = handles.orig_RowName;
band_per_col = handles.orig_ColName;

%sd_correction = handles.sdcorrection.Value;
save_spectral_dq(avg_fft,chans_per_row, band_per_col, selected_DQ_type, empty_filename,'xlsx')

% --- Executes on button press in pushbutton_mat.
function pushbutton_mat_Callback(hObject, eventdata, handles)
% hObject    handle to pushbutton_mat (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)
empty_filename = [];
selected_DQ_type = handles.popupmenu_DQ_type.Value;
avg_fft = handles.avg_fft; %only the average across frequencies within band
chans_per_row = handles.orig_RowName;
band_per_col = handles.orig_ColName;

%sd_correction = handles.sdcorrection.Value;
save_spectral_dq(avg_fft,chans_per_row, band_per_col, selected_DQ_type, empty_filename,'mat')

% --- Executes on button press in pushbutton_help.
function pushbutton_help_Callback(hObject, eventdata, handles)
% hObject    handle to pushbutton_help (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)
web https://github.com/lucklab/erplab/wiki/Spectral-Data-Quality-(continuous-eeg) -browser

% --- Executes on button press in pushbutton_done.
function pushbutton_done_Callback(hObject, eventdata, handles)
% hObject    handle to pushbutton_done (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)
%figure1_CloseRequestFcn(hObject,eventdata,handles)
delete(handles.figure1);

% --- Executes when user attempts to close figure1.
function figure1_CloseRequestFcn(hObject, eventdata, handles)
% hObject    handle to figure1 (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

% Hint: delete(hObject) closes the figure
delete(hObject);

% --- Executes on button press in checkbox_heatmap.
function checkbox_heatmap_Callback(hObject, eventdata, handles)
% hObject    handle to checkbox_heatmap (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

% Hint: get(hObject,'Value') returns toggle state of checkbox_heatmap
heatmap_on = get(hObject,'Value');

if heatmap_on == 1

    %if heatmap_on, then outliers not possible so clear outliers
    clear_outliers(hObject, eventdata, handles);
    handles.outliers_on = 0;
    set(handles.checkbox_outliers, 'Value', 0);

    redraw_heatmap(hObject, eventdata, handles);
    handles.heatmap_on = 1;
else
    clear_heatmap(hObject, eventdata, handles);
    handles.heatmap_on = 0;

    %if outliers is on, then keep outliers active
    if handles.outliers_on == 1
        redraw_outliers(hOBject, eventdata,handles);
    end

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

function redraw_heatmap(hObject, eventdata, handles)

color_map = viridis;
data = handles.dq_table.Data;

if any([isnan(data(1,:))])
    disp('Cannot create heatmap if NaN is present')
    return
end

data_min = min(data(:));
range_val = max(data(:)) - data_min;
range_colormap = size(color_map,1); % as in, 256 shades?
val_increase_per_shade = range_val / range_colormap;

% Use this @ anonymous function to make HTML tag in box in loop below
colergen = @(color,text) ['

',text,'

'];

for cell = 1:numel(data)
    data_here = data(cell);
    shades_up_here = round((data_here - data_min) / val_increase_per_shade);
    if shades_up_here < 1
        shades_up_here = 1;
    end
    if shades_up_here > range_colormap
        shades_up_here = range_colormap;
    end
    RGB_here = color_map(shades_up_here,:);
    hex_color_here = ['#' dec2hex(round(255*RGB_here(1)),2) dec2hex(round(255*RGB_here(2)),2) dec2hex(round(255*RGB_here(3)),2)];
    data2{cell} = colergen(hex_color_here,num2str(data_here));
end

data2 = reshape(data2,size(data));

handles.dq_table.Data = data2;

function clear_heatmap(hObject, eventdata, handles)
handles.dq_table.Data = handles.orig_data;

% --- Executes when entered data in editable cell(s) in dq_table.
function dq_table_CellEditCallback(hObject, eventdata, handles)
% hObject    handle to dq_table (see GCBO)
% eventdata  structure with the following fields (see MATLAB.UI.CONTROL.TABLE)
%	Indices: row and column indices of the cell(s) edited
%	PreviousData: previous data for the cell(s) edited
%	EditData: string(s) entered by the user
%	NewData: EditData or its converted form set on the Data property. Empty if Data was not changed
%	Error: error string when failed to convert EditData to appropriate value for Data
% handles    structure with handles and user data (see GUIDATA)

% --- Executes on selection change in newerpwin.
function newerpwin_Callback(hObject, eventdata, handles)
% hObject    handle to newerpwin (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

% Hints: contents = cellstr(get(hObject,'String')) returns newerpwin contents as cell array
%        contents{get(hObject,'Value')} returns selected item from newerpwin

% --- Executes during object creation, after setting all properties.
function newerpwin_CreateFcn(hObject, eventdata, handles)
% hObject    handle to newerpwin (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    empty - handles not created until after all CreateFcns called

% Hint: popupmenu 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

% --- Executes on button press in pushbutton_selERP.
function pushbutton_selERP_Callback(hObject, eventdata, handles)
% hObject    handle to pushbutton_selERP (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)
selected_erpset = handles.newerpwin.Value;
sel_ERP = handles.ALLERP(selected_erpset);
DQ_Table_GUI(sel_ERP, handles.ALLERP, selected_erpset,2);

% --- Executes on button press in checkbox_outliers.
function checkbox_outliers_Callback(hObject, eventdata, handles)
% hObject    handle to checkbox_outliers (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

% Hint: get(hObject,'Value') returns toggle state of checkbox_outliers

% Hint: get(hObject,'Value') returns toggle state of checkbox_heatmap
outliers_on = get(hObject,'Value'); % this is not the same as handles.outliers_on

if outliers_on == 1
    set(handles.stdwindow,'Enable','On');
%     set(handles.chwindow,'Enable','On');
%     set(handles.chbutton,'Enable','On');
else
    set(handles.stdwindow,'Enable','Off');
%     set(handles.chwindow,'Enable','Off');
%     set(handles.chbutton,'Enable','Off');
end

%updated the handles object for outliers to be in sync with hObject
if outliers_on == 1
    handles.outliers_on = 1;
else
    handles.outliers_on = 0;
end

if handles.outliers_on == 1

    %if heatmap is currently on, stop it and turn its states off
    %before applying outliers fxn

    if handles.heatmap_on == 1
        clear_heatmap(hObject, eventdata, handles);
        handles.heatmap_on = 0;
        set(handles.checkbox_heatmap, 'Value', 0);

        %handles.outliers_on = 1;
        redraw_outliers(hObject, eventdata, handles);

    else
        %handles.outliers_on = 1;
        clear_outliers(hObject, eventdata, handles);
        %guidata(hObject, handles);
        redraw_outliers(hObject, eventdata, handles);

    end

else %but if outliers is off, and heatmap is wanted, you can do it

%     if handles.heatmap_on ==1
%
%         clear_outliers(hObject, eventdata, handles);
%         redraw_heatmap(hObject, eventdata, handles);
%         handles.outliers_on = 0;
%
%     else
%         clear_outliers(hObject, eventdata, handles);
%         handles.outliers_on = 0;
%
%     end

    %also, redraw initial data/channels without outliers
    clear_outliers(hObject, eventdata, handles)

end

% %if heatmap is currently on, stop it and turn its states off
% %before applying outliers fxn
% if handles.heatmap_on == 1
%
%     clear_heatmap(hObject, eventdata, handles);
%     handles.heatmap_on = 0;
%     set(handles.checkbox_heatmap, 'Value', 0);
%
%     if outliers_on == 1
%         handles.outliers_on = 1;
%         redraw_outliers(hObject, eventdata, handles);
%     else
%
%         set(handles.stdwindow,'Enable','Off');
%         set(handles.chwindow,'Enable','Off');
%         set(handles.chbutton,'Enable','Off');
%         clear_outliers(hObject, eventdata, handles);
%         handles.outliers_on = 0;
%
%     end
%
% else
%
%      if outliers_on == 1
%         handles.outliers_on = 1;
%         redraw_outliers(hObject, eventdata, handles);
%      else
%         set(handles.stdwindow,'Enable','Off');
%         set(handles.chwindow,'Enable','Off');
%         set(handles.chbutton,'Enable','Off');
%         clear_outliers(hObject, eventdata, handles);
%         handles.outliers_on = 0;
%
%      end
%
% end
% Update handles structure
guidata(hObject, handles);

function redraw_outliers(hObject, eventdata, handles)
handles = guidata(hObject); %get the most updated handles struct
data = handles.dq_table.Data;

if iscell(data)
    %in case of new std from stdwindow
    data = handles.orig_data;
end
%
% %either default(all) or from chbutton
% chans_to_use = handles.indxlistch;
% data = data(chans_to_use,:);
% n_elec = length(chans_to_use);
% ERPtouse = handles.ERP.chanlocs;
%
% for i=1:n_elec
%     curr_chanInd = chans_to_use(i);
%     elec_labels{i} = ERPtouse(curr_chanInd).labels;
% end

%save new labels to handles
% handles.dq_table.RowName = elec_labels;

%compute outliers as X standard deviations from mean across channels

col_means = mean(data,1);
col_std = std(data,0,1);
chosen_std = str2double((handles.stdwindow.String));

% Use this @ anonymous function to make HTML tag in box in loop below
colergen = @(color,text) ['

',text,'

'];

for cell = 1:numel(data)

    data_here = data(cell);
    [d_row,d_col] = ind2sub(size(data),cell);

    neg_thresh = col_means(d_col) - (chosen_std*col_std(d_col));
    pos_thresh = col_means(d_col) + (chosen_std*col_std(d_col));

    if data_here < neg_thresh || data_here > pos_thresh
        hex_color_here = ['#' dec2hex(255,2) dec2hex(0,2) dec2hex(0,2)];
        data2{cell} = colergen(hex_color_here,num2str(data_here));
    else
        data2{cell} = num2str(data_here);

    end

end

data2= reshape(data2,size(data));
handles.dq_table.Data = data2;

function clear_outliers(hObject, eventdata, handles)
%handles = guidata(hObject);
handles.dq_table.Data = handles.orig_data;
handles.dq_table.RowName = handles.orig_RowName;
%handles.indxlistch = handles.orig_indxlistch;

% Re-Prepare List of current Channels
%chanArray = handles.orig_indxlistch; %default
%handles.indxlistch = chanArray; %default
%set(handles.chwindow, 'String', vect2colon(chanArray, 'Delimiter', 'off'));
 % Update handles structure
guidata(hObject, handles);

function stdwindow_Callback(hObject, eventdata, handles)
% hObject    handle to stdwindow (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

newstd = get(hObject,'String');
handles.stdwindow.String= newstd;
redraw_outliers(hObject, eventdata, handles);

guidata(hObject, handles);

% Hints: get(hObject,'String') returns contents of stdwindow as text
%        str2double(get(hObject,'String')) returns contents of stdwindow as a double

% --- Executes during object creation, after setting all properties.
function stdwindow_CreateFcn(hObject, eventdata, handles)
% hObject    handle to stdwindow (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    empty - handles not created until after all CreateFcns called

% 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 chwindow_Callback(hObject, eventdata, handles)
% % hObject    handle to chwindow (see GCBO)
% % eventdata  reserved - to be defined in a future version of MATLAB
% % handles    structure with handles and user data (see GUIDATA)
%
% % Hints: get(hObject,'String') returns contents of chwindow as text
% %        str2double(get(hObject,'String')) returns contents of chwindow as a double
%
% ch = get(hObject,'string');
% nch = length(handles.listch);
%
% ch_val = eval(ch);
%
% if ~isempty(ch)
%
% %     testch = regexpi(ch,':') ;
% %
% %     if isempty(testch)
% %
% %         try
% %             ch_val = eval(ch);
% %             handles.listch(ch_val); %attempt to index
% %
% %         catch
% %             msgboxText =  'Invalid channel input: please enter channels as integers, with a space between each number, or ch1:chN array syntax';
% %             title = 'ERPLAB: channel GUI input';
% %             errorfound(msgboxText, title);
% %
% %         end
% %     else
% %         ch_val = eval(ch);
% %
% %     end
% %
%
%     tf = checkchannels(ch_val, nch, 1);
%
%     if tf
%         return
%     end
%
%
% %     if length(ch_val) > length(handles.orig_indxlistch)
% %         msgboxText =  'Exceeded Number of Available Channels';
% %         title = 'ERPLAB: channel GUI input';
% %         errorfound(msgboxText, title);
% %         return
% %     end
% %
%     set(handles.chwindow, 'String', vect2colon(ch_val, 'Delimiter', 'off'));
%     handles.indxlistch = ch_val;
%     % Update handles structure
%     guidata(hObject, handles);
%     redraw_outliers(hObject, eventdata, handles);
%     guidata(hObject, handles);
%
% else
%     msgboxText =  'Not valid channel';
%     title = 'ERPLAB: channel GUI input';
%     errorfound(msgboxText, title);
%     return
% end
%

% --- Executes during object creation, after setting all properties.
function chwindow_CreateFcn(hObject, eventdata, handles)
% hObject    handle to chwindow (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    empty - handles not created until after all CreateFcns called

% 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

% --- Executes on button press in chbutton.
function chbutton_Callback(hObject, eventdata, handles)
% hObject    handle to chbutton (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)
listch = handles.listch; %true channels list as labels
indxlistch = handles.indxlistch; %true array of channel indexs
indxlistch = indxlistch(indxlistch<=length(listch)); % true array of channel index if less than expected labels
titlename = 'Select Channel(s)';
if get(hObject, 'Value')
        if ~isempty(listch)
                ch = browsechanbinGUI(listch, indxlistch, titlename);
                if ~isempty(ch)
                        set(handles.chwindow, 'String', vect2colon(ch, 'Delimiter', 'off'));
                        handles.indxlistch = ch;
                        % Update handles structure
                        guidata(hObject, handles);
                        redraw_outliers(hObject, eventdata, handles);
                        guidata(hObject, handles);

                else
                        disp('User selected Cancel')
                        return
                end
        else
                msgboxText =  'No channel information was found';
                title = 'ERPLAB: ploterp GUI input';
                errorfound(msgboxText, title);
                return
        end
end

% --- Executes during object creation, after setting all properties.
function pushbutton_selERP_CreateFcn(hObject, eventdata, handles)
% hObject    handle to pushbutton_selERP (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    empty - handles not created until after all CreateFcns called

% --- Executes on selection change in active_eegset.
function active_eegset_Callback(hObject, eventdata, handles)
% hObject    handle to active_eegset (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

newactive = get(hObject,'Value');
ERPs = handles.ALLERP;
ERP = ERPs(newactive);

n_dq = numel(ERP.dataquality);
for i=1:n_dq
    type_names{i} = ERP.dataquality(i).type;
end

handles.popupmenu_DQ_type.String = type_names;
handles.popupmenu_DQ_type.Value = n_dq;

[n_elec, n_tw, n_bin] = size(ERP.dataquality(n_dq).data);

if n_elec == 0
    % if the dq data is of zero size, then the 3rd dim may be mis-sized.
    % So set to default of ERP.nbin instead.
    n_bin = ERP.nbin;
end

n_bin_names = numel(ERP.bindescr);
if n_bin_names == n_bin

    for i=1:n_bin
        bin_names{i} = ['BIN ' num2str(i) ' - ' ERP.bindescr{i}];
    end
else
    % if not every bin has a bin description, leave then off
    for i=1:n_bin
        bin_names{i} = ['BIN ' num2str(i)];
    end
end

handles.popupmenu_bin.String = bin_names;

% handles.popupmenu_DQ_type.Value saves the indx of the selected DQ Measure
selected_DQ_type = handles.popupmenu_DQ_type.Value;
selected_bin = handles.popupmenu_bin.Value;

table_data = ERP.dataquality(selected_DQ_type).data(:,:,selected_bin);
handles.dq_table.Data = table_data;
handles.orig_data = table_data;

% electrode labels from ERPset, iff present and number matches
if isfield(ERP.chanlocs,'labels') && numel(ERP.chanlocs) == n_elec
    for i=1:n_elec
        elec_labels{i} = ERP.chanlocs(i).labels;
    end
else

    for i=1:ERP.nchan
        elec_labels{i} = i;
    end
end
handles.dq_table.RowName = elec_labels;

% Time-window labels
if isfield(ERP.dataquality(selected_DQ_type),'time_window_labels') && isempty(ERP.dataquality(selected_DQ_type).time_window_labels) == 0  && handles.checkbox_text_labels.Value == 1
    tw_labels = ERP.dataquality(selected_DQ_type).time_window_labels;
elseif isempty(ERP.dataquality(selected_DQ_type).times)
    tw_labels = [];
else
    for i=1:n_tw
        tw_labels{i} = [num2str(ERP.dataquality(selected_DQ_type).times(i,1)) ' : ' num2str(ERP.dataquality(selected_DQ_type).times(i,2))];
    end
end
if n_tw == 0
    tw_labels = 'No data here. Perhaps it was cleared?';
end
handles.dq_table.ColumnName = tw_labels;

% Set font size of DQ Table
desired_fontsize = erpworkingmemory('fontsizeGUI');
handles.dq_table.FontSize = desired_fontsize;
handles.heatmap_on = 0;

%
% Prepare List of current Channels
%
nchan  = ERP.nchan; % Total number of channels
if ~isfield(ERP.chanlocs,'labels')
        for e=1:nchan
                ERP.chanlocs(e).labels = ['Ch' num2str(e)];
        end
end
listch = {''};
for ch =1:nchan
        listch{ch} = [num2str(ch) ' = ' ERP.chanlocs(ch).labels ];
end
handles.listch     = listch;

chanArray = 1:ERP.nchan; %default
handles.indxlistch = chanArray; %default
set(handles.chwindow, 'String', vect2colon(chanArray, 'Delimiter', 'off'));

% Set outliers text window off
set(handles.stdwindow,'Enable','Off');
set(handles.chwindow,'Enable','Off');
set(handles.chbutton,'Enable','Off');

% Update handles structure
handles.ERP = ERP;
handles.ALLERP = ERPs;

%clear all perviously set options
set(handles.checkbox_outliers,'Value',0);
set(handles.checkbox_heatmap,'Value',0);
set(handles.checkbox_text_labels,'Value',0);
set(handles.stdwindow,'Enable','Off');
set(handles.chwindow,'Enable','Off');
set(handles.chbutton,'Enable','Off');

%check to see if outliers is on prior to switching active erpset?
%checkbox_outliers_Callback(hObject, eventdata, handles)

% outliers_on = get(handles.checkbox_outliers,'Value');
%
% if outliers_on == 1
%     handles.outliers_on = 1;
%     set(handles.stdwindow,'Enable','On');
%     set(handles.chwindow,'Enable','On');
%     set(handles.chbutton,'Enable','On');
%     redraw_outliers(hObject, eventdata, handles);
%
% else
%     handles.heatmap_on = 0;
%     set(handles.stdwindow,'Enable','Off');
%     set(handles.chwindow,'Enable','Off');
%     set(handles.chbutton,'Enable','Off');
%     clear_outliers(hObject, eventdata, handles);
%
%
% end

guidata(hObject, handles);

% Hints: contents = cellstr(get(hObject,'String')) returns active_eegset contents as cell array
%        contents{get(hObject,'Value')} returns selected item from active_eegset

% --- Executes during object creation, after setting all properties.
function active_eegset_CreateFcn(hObject, eventdata, handles)
% hObject    handle to active_eegset (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    empty - handles not created until after all CreateFcns called

% Hint: popupmenu 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 tf = checkchannels(chx, nchan, showmsg)

if nargin<3
        showmsg = 1;
end
tf = 0; % no problem by default

if ~mod(chx, 1) == 0
    if showmsg
        msgboxText =  'Invalid channel indexing.';
        title = 'ERPLAB: basicfilterGUI() error:';
        errorfound(msgboxText, title);
    end
    tf = 1; %

end

if isempty(chx)
        if showmsg
                msgboxText =  'Invalid channel indexing.';
                title = 'ERPLAB: basicfilterGUI() error:';
                errorfound(msgboxText, title);
        end
        tf = 1; %
        return
end
if ~isempty(find(chx>nchan))
        if showmsg
                msgboxText =  ['You only have %g channels,\n'...
                        'so you cannot specify indices greater than this.'];
                title = 'ERPLAB: basicfilterGUI() error:';
                errorfound(sprintf(msgboxText, nchan), title);
        end
        tf = 1; %
        return
end
if ~isempty(find(chx<1))
        if showmsg
                msgboxText =  'You cannot use zero or a negative number as a channel indexing';
                title = 'ERPLAB: basicfilterGUI() error:';
                errorfound(msgboxText, title);
        end
        tf = 1; %
        return
end
if length(chx)>length(unique_bc2(chx))
        if showmsg
                msgboxText =  ['Repeated channels are not allowed.\n'...
                        'Therefore, ERPLAB will get rid of them.'];
                title = 'ERPLAB: basicfilterGUI() error:';
                errorfound(sprintf(msgboxText), title, [1 1 0], [0 0 0], 0)
        end
        tf = 0; %
        return
end

% --- Executes on button press in sdcorrection.
function sdcorrection_Callback(hObject, eventdata, handles)
% hObject    handle to sdcorrection (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)
selected_DQ_type = handles.popupmenu_DQ_type.Value;
selected_bin = handles.popupmenu_bin.Value;

if handles.sdcorrection.Value == 1
    table_data = handles.ERP.dataquality(selected_DQ_type).data.SD_unbias(:,:,selected_bin); %SD (Gurland & Tripathi, 1971)
else
    table_data = handles.ERP.dataquality(selected_DQ_type).data.SD_bias(:,:,selected_bin); %SD (N-1)
end
handles.dq_table.Data = table_data;
handles.orig_data = table_data;

%update handles struct
guidata(hObject, handles);

% Hint: get(hObject,'Value') returns toggle state of sdcorrection

% --- Executes when selected cell(s) is changed in dq_table.
function dq_table_CellSelectionCallback(hObject, eventdata, handles)
% hObject    handle to dq_table (see GCBO)
% eventdata  structure with the following fields (see MATLAB.UI.CONTROL.TABLE)
%	Indices: row and column indices of the cell(s) currently selecteds
% handles    structure with handles and user data (see GUIDATA)
if numel(eventdata.Indices)
    row_here = eventdata.Indices(:,1);
    col_here = eventdata.Indices(:,2);
    old_row = handles.sel_row;
    old_col = handles.sel_col;
    handles.sel_row = row_here;
    handles.sel_col = col_here;
%     if isequal(row_here,old_row) == 0
%         row_text = ['Row number ' num2str(row_here) ' now selected'];
%         col_text = ['Column Number' num2str(col_here) ' now selected'];
%         disp(row_text)
%     end
    guidata(hObject, handles);

% catch
%    eventdata = handles.ch;
%    guidata(hObject,handles);
%

end

% --- Executes on button press in plot_button.
function plot_button_Callback(hObject, eventdata, handles)
% hObject    handle to plot_button (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

chans_to_plot = handles.sel_row;
labels_to_plot = handles.sel_col;
yout = handles.yout;
yout_lab = handles.yout_lab;
msgn = 'whole';

if numel(unique_bc2(labels_to_plot)) > 1

    msgboxText =  'Sorry, you may only plot one frequency band at a time!';
    title_header = 'ERPLAB: pop_continuousFFT() error';
    errorfound(msgboxText, title_header);

    return

elseif any(isnan(yout{unique_bc2(labels_to_plot)}))
    msgboxText =  'Sorry, you cannot plot this frequency band due to NaNs';
    title_header = 'ERPLAB: pop_continuousFFT() error';
    errorfound(msgboxText, title_header);

else

    dq_selected = handles.popupmenu_DQ_type.Value;

    if (strcmp(handles.popupmenu_DQ_type.String{dq_selected},'Amplitude'))

        fname = handles.EEG.setname;
        h = figure('Name',['<< ' fname ' >>  ERPLAB Amplitude Spectrum'],...
            'NumberTitle','on', 'Tag','Plotting Spectrum',...
            'Color',[1 1 1]);
        p_yout = yout{unique_bc2(labels_to_plot)};

        if numel(chans_to_plot) == 1
            p_yout_sel = p_yout(:,chans_to_plot); %only use selected channels
            p_yout_lab = yout_lab{unique_bc2(labels_to_plot)}';
            lege = sprintf('EEG Channel: ');
        else
            p_yout_sel = p_yout(:,chans_to_plot); %only use selected channels
            p_yout_sel = mean(p_yout_sel,2); %average across selected channels;
            p_yout_lab = yout_lab{unique_bc2(labels_to_plot)}';
            lege = sprintf('EEG Channel Average: ');
        end

        if numel(p_yout_lab) == 1 %impulse at one frequency/ create impulse
            plot([p_yout_lab p_yout_lab], [0 p_yout_sel]);
        else
            plot(p_yout_lab,p_yout_sel);
        end

        try
            axis([min(p_yout_lab)  max(p_yout_lab)  min(p_yout_sel)*0.9 max(p_yout_sel)*1.1])
        catch % in the case of impulse
            axis([min(p_yout_lab)*0.99  max(p_yout_lab)*1.01  min(p_yout_sel)*0.99 max(p_yout_sel)*1.01])
        end

        if isfield(handles.EEG.chanlocs,'labels')
           % lege = sprintf('EEG Channel Average: ');
            for i=1:length(chans_to_plot)
                lege =   sprintf('%s %s', lege, handles.EEG.chanlocs(chans_to_plot(i)).labels);
            end
            lege = sprintf('%s *%s', lege, msgn);
            legend(lege)
        else
            legend(['EEG Channel: ' vect2colon(chanArray,'Delimiter', 'off') '  *' msgn])
        end

        title('Single-Sided Amplitude Spectrum of y(t)')
        xlabel('Frequency (Hz)')
        ylabel('Amplitude - absolute single-sided (original units)')

        %     if plot_type == 1
        %         set(gca,'XScale','log')
        %         set(gca,'XTick',[1 10 60 100])
        %         if f1 == 0
        %             xstart = 0.1
        %         else
        %             xstart = f1;
        %         end
        %         xlim = [xstart f2];
        %         xlabel('Frequency (Hz) - log scale')
        %
        %     end
    else %choose power

        fname = handles.EEG.setname;
        h = figure('Name',['<< ' fname ' >>  ERPLAB Power Spectrum'],...
            'NumberTitle','on', 'Tag','Plotting Spectrum',...
            'Color',[1 1 1]);
        p_yout = yout{unique_bc2(labels_to_plot)};

        if numel(chans_to_plot) == 1
            p_yout_sel = (p_yout(:,chans_to_plot)).^2; %only use selected channels
            p_yout_lab = yout_lab{unique_bc2(labels_to_plot)}';
            lege = sprintf('EEG Channel: ');
        else
            p_yout_sel = (p_yout(:,chans_to_plot)).^2; %only use selected channels
            p_yout_sel = mean(p_yout_sel,2); %average across selected channels;
            p_yout_lab = yout_lab{unique_bc2(labels_to_plot)}';
            lege = sprintf('EEG Channel Average: ');
        end

%         plot(p_yout_lab,p_yout_sel);
%         axis([min(p_yout_lab)  max(p_yout_lab)  min(p_yout_sel)*0.9 max(p_yout_sel)*1.1])
%
        if numel(p_yout_lab) == 1 %impulse at one frequency/ create impulse
            plot([p_yout_lab p_yout_lab], [0 p_yout_sel]);
        else
            plot(p_yout_lab,p_yout_sel);
        end

        try
            axis([min(p_yout_lab)  max(p_yout_lab)  min(p_yout_sel)*0.9 max(p_yout_sel)*1.1])
        catch % in the case of impulse
            axis([min(p_yout_lab)*0.99  max(p_yout_lab)*1.01  min(p_yout_sel)*0.99 max(p_yout_sel)*1.01])
        end

        if isfield(handles.EEG.chanlocs,'labels')
            %lege = sprintf('EEG Channel Average: ');
            for i=1:length(chans_to_plot)
                lege =   sprintf('%s %s', lege, handles.EEG.chanlocs(chans_to_plot(i)).labels);
            end
            lege = sprintf('%s *%s', lege, msgn);
            legend(lege)
        else
            legend(['EEG Channel: ' vect2colon(chanArray,'Delimiter', 'off') '  *' msgn])
        end

        title('Single-Sided Power Spectrum of y(t)')
        xlabel('Frequency (Hz)')
        ylabel('Power - absolute single-sided (original units)')

        %     if plot_type == 1
        %         set(gca,'XScale','log')
        %         set(gca,'XTick',[1 10 60 100])
        %         if f1 == 0
        %             xstart = 0.1
        %         else
        %             xstart = f1;
        %         end
        %         xlim = [xstart f2];
        %         xlabel('Frequency (Hz) - log scale')
        %
        %     end

    end

end

function edit_channels_Callback(hObject, eventdata, handles)
% hObject    handle to edit_channels (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

% Hints: get(hObject,'String') returns contents of edit_channels as text
%        str2double(get(hObject,'String')) returns contents of edit_channels as a double

ch = get(hObject,'string');
nch = length(handles.listch);

ch_val = eval(ch);

if ~isempty(ch)

    tf = checkchannels(ch_val, nch, 1);

    if tf
        return
    end

%     if length(ch_val) > length(handles.orig_indxlistch)
%         msgboxText =  'Exceeded Number of Available Channels';
%         title = 'ERPLAB: channel GUI input';
%         errorfound(msgboxText, title);
%         return
%     end
%
    set(handles.edit_channels, 'String', vect2colon(ch_val, 'Delimiter', 'off'));
    handles.listch = ch_val;
    %

    % Update handles structure
    guidata(hObject, handles);

else
    msgboxText =  'Not valid channel';
    title = 'ERPLAB: channel GUI input';
    errorfound(msgboxText, title);
    return
end

% --- Executes during object creation, after setting all properties.
function edit_channels_CreateFcn(hObject, eventdata, handles)
% hObject    handle to edit_channels (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    empty - handles not created until after all CreateFcns called

% 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

% --- Executes on button press in browse_button.
function browse_button_Callback(hObject, eventdata, handles)
% hObject    handle to browse_button (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)
listch     = handles.listch;
nchan = length(listch);
for ch =1:nchan
        listch{ch} = [num2str(ch) ' = ' listch{ch}];
end

indxlistch = handles.indxlistch;
indxlistch = indxlistch(indxlistch<=length(listch));
titlename  = 'Select Channel(s)';

if get(hObject, 'Value')
        if ~isempty(listch)
                ch = browsechanbinGUI(listch, indxlistch, titlename);
                if ~isempty(ch)
                        set(handles.edit_channels, 'String', vect2colon(ch, 'Delimiter', 'off'));
                        handles.indxlistch = ch;
                        handles.sel_row = ch;

                        %dq table callback = eventdata.Indices(:,1);
                        test.Indices(:,1) = ch;

                        jUIscrollPane = findobj(handles.dq_table);
                        jUITable = jUIscrollPane.getViewport.getView;

%                         %dq_table_CellSelectionCallback(hObject, eventdata, handles)
%                         dataTB = handles.dq_table;
%                         dataTB.Selection = [ch];

                        dq_table_CellSelectionCallback(hObject, test, handles);

                        % Update handles structure
                        guidata(hObject, handles);
                else
                        disp('User selected Cancel')
                        return
                end
        else
                msgboxText =  'No channel information was found';
                title = 'ERPLAB: fourieeg GUI input';
                errorfound(msgboxText, title);
                return
        end
end

% --- Executes on button press in all_button.
function all_button_Callback(hObject, eventdata, handles)
% hObject    handle to all_button (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)
listch     = handles.listch; %original list
nChan = 1:length(listch);
handles.indxlistch = nChan;
set(handles.edit_channels, 'String', vect2colon(nChan, 'Delimiter', 'off'));
guidata(hObject,handles);

% --------------------------------------------------------------------
function dq_table_ButtonDownFcn(hObject, eventdata, handles)
% hObject    handle to dq_table (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

author: ucdavis

repo: erplab

save_path: github-repos/MATLAB/ucdavis-erplab


text:
%  Checks whether the fitness is too low
%
%  ::
%
%    flag=is_violation(f,penalty)
%
%  Args:
%
%     - **f** [numeric]: fitness
%     - **penalty** [numeric]: lowest possible fitness
%
%  Returns:
%     :
%
%     - **flag** [true|false]:
%
%

author: jmaih


text:
function frequency = Frequency(timestamps,varargin)

%Frequency - Compute instantaneous frequency for a point process (e.g. a spike train).
%
% Compute instantaneous frequency either by smoothing the point process using a fixed
% or adaptive Gaussian kernel, or by simply taking the inverse of the interevent intervals.
% Adaptive kernel smoothing is described in Richmond et al. (1990).
%
%  USAGE
%
%    frequency = Frequency(timestamps,<options>)
%
%    timestamps     list of timestamps
%    <options>      optional list of property-value pairs (see table below)
%
%    =========================================================================
%     Properties    Values
%    -------------------------------------------------------------------------
%     'method'      convolute with a fixed ('fixed', default) or adaptive
%                   ('adaptive') Gaussian kernel, or compute inverse inter-
%                   event intervals ('inverse')
%     'limits'      [start stop] in seconds (default = approx. first and last
%                   timestamps)
%     'binSize'     bin size in seconds (default = 0.050)
%     'smooth'      Gaussian kernel width in number of samples (default = 2)
%     'show'        plot results (default = 'off')
%    =========================================================================
%
%  NOTE
%
%    While the convolution method returns a Mx2 matrix (time bins in column 1,
%    frequencies in column 2), the inverse interspike intervals method returns
%    a Nx2 matrix (original spike timestamps in column 1, frequencies in column 2)

% Copyright (C) 2004-2011 by Michaël Zugaro
%
% 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.

% Defaults
method = 'fixed';
binSize = 0.05;
smooth = 2;
limits = [];
show = 'off';

% Check number of parameters
if nargin < 1,
  error('Incorrect number of parameters (type ''help Frequency'' for details).');
end

% Check parameter sizes
if size(timestamps,2) ~= 1,
	error('Parameter ''timestamps'' is not a vector (type ''help Frequency'' for details).');
end

% Parse parameter list
for i = 1:2:length(varargin),
	if ~ischar(varargin{i}),
		error(['Parameter ' num2str(i+2) ' is not a property (type ''help Frequency'' for details).']);
	end
	switch(lower(varargin{i})),

		case 'method',
			method = lower(varargin{i+1});
			if ~isstring_FMAT(method,'fixed','adaptive','inverse','iisi'),
				error('Incorrect value for property ''method'' (type ''help Frequency'' for details).');
			end

		case 'binsize',
			binSize = varargin{i+1};
			if ~isdscalar(binSize,'>=0'),
				error('Incorrect value for property ''binSize'' (type ''help Frequency'' for details).');
			end

		case 'limits',
			limits = varargin{i+1};
			if ~isdvector(limits,'#2'),
				error('Incorrect value for property ''limits'' (type ''help Frequency'' for details).');
			end

		case 'smooth',
			smooth = varargin{i+1};
			if ~isdscalar(smooth,'>=0'),
				error('Incorrect value for property ''smooth'' (type ''help Frequency'' for details).');
			end

		case 'show',
			show = varargin{i+1};
			if ~isstring_FMAT(show,'on','off'),
				error('Incorrect value for property ''show'' (type ''help Frequency'' for details).');
			end

		otherwise,
			error(['Unknown property ''' num2str(varargin{i}) ''' (type ''help Frequency'' for details).']);

  end
end

if isempty(limits),
	limits = [timestamps(1)-10*binSize timestamps(end)+10*binSize];
end

if strcmp(method,'inverse') | strcmp(method,'iisi'),
	% Inter-spike intervals
	ds = diff(timestamps);
	f1 = 1./ds;
	t = timestamps(1:end-1) + ds/2;
	f2 = Interpolate([t f1],timestamps(2:end-1));
	f2 = [timestamps(1) 0;f2;timestamps(end) 0];
	frequency = f2;
else
	% Smoothing
	% 1) Fixed-kernel
	t = (limits(1):binSize:limits(2))';
	T = Bin(timestamps,t);
	binned = Accumulate(T,1,size(t));
	f = Smooth(binned/binSize,smooth);
	frequency = [t f];
	if strcmp(method,'adaptive'),
		% 2) Variable-kernel (requires the above 'pilot' fixed-kernel estimate)
		% Compute variable-kernel sigma
		N = length(t);
		mu = exp(1/N*sum(log(f(f~=0)))); % Correction: geometric mean is computed using only non-zero values
		lambda = sqrt(f/mu);
		sigma = (smooth*binSize)./lambda;
		sigma = Clip(sigma,0,N*binSize/3);
		% Perform variable-kernel smoothing
		binned = [flipud(binned);binned;flipud(binned)];
		for i = 1:N,
			% Gaussian
			x = (0:binSize:3*sigma(i))';
			x = [flipud(-x);x(2:end)];
			kernel = exp(-x.^2/sigma(i)^2);
			kernel = kernel/sum(kernel);
			% 'Pointwise' convolution
			n = (length(kernel)-1)/2;
			bins = N + i + (-n:n);
			S2(i) = sum(binned(bins).*kernel)/binSize;
		end
		frequency = [t S2'];
	end
end

if strcmp(lower(show),'on'),
	PlotXY(frequency);
	hold on;
	PlotTicks(timestamps,'size',10,'k');
end