text:
function [cfg] = ft_singleplotER(cfg, varargin)

% FT_SINGLEPLOTER plots the event-related fields or potentials of a single
% channel or the average over multiple channels. Multiple datasets can be
% overlayed.
%
% Use as
%   ft_singleplotER(cfg, data)
% or
%   ft_singleplotER(cfg, data1, data2, ..., datan)
%
% The data can be an erp/erf produced by FT_TIMELOCKANALYSIS, a power
% spectrum produced by FT_FREQANALYSIS or connectivity spectrum produced by
% FT_CONNECTIVITYANALYSIS.
%
% The configuration can have the following parameters:
%   cfg.parameter     = field to be plotted on y-axis, for example 'avg', 'powspctrm' or 'cohspctrm' (default is automatic)
%   cfg.maskparameter = field in the first dataset to be used for masking of data; this is not supported when
%                       computing the mean over multiple channels, or when giving multiple input datasets (default = [])
%   cfg.maskstyle     = style used for masking of data, 'box', 'thickness' or 'saturation' (default = 'box')
%   cfg.maskfacealpha = mask transparency value between 0 and 1
%   cfg.xlim          = 'maxmin' or [xmin xmax] (default = 'maxmin')
%   cfg.ylim          = 'maxmin', 'maxabs', 'zeromax', 'minzero', or [ymin ymax] (default = 'maxmin')
%   cfg.channel       = Nx1 cell-array with selection of channels (default = 'all'), see FT_CHANNELSELECTION for details
%   cfg.title         = string, title of plot
%   cfg.showlegend    = 'yes' or 'no', show the legend with the colors (default = 'no')
%   cfg.refchannel    = name of reference channel for visualising connectivity, can be 'gui'
%   cfg.baseline      = 'yes', 'no' or [time1 time2] (default = 'no'), see ft_timelockbaseline
%   cfg.baselinetype  = 'absolute' or 'relative' (default = 'absolute')
%   cfg.trials        = 'all' or a selection given as a 1xn vector (default = 'all')
%   cfg.fontsize      = font size of title (default = 8)
%   cfg.hotkeys       = enables hotkeys (leftarrow/rightarrow/uparrow/downarrow/m) for dynamic zoom and translation (ctrl+) of the axes
%   cfg.interactive   = interactive plot 'yes' or 'no' (default = 'yes')
%                       in a interactive plot you can select areas and produce a new
%                       interactive plot when a selected area is clicked. multiple areas
%                       can be selected by holding down the shift key.
%   cfg.renderer      = 'painters', 'zbuffer', ' opengl' or 'none' (default = [])
%   cfg.linestyle     = linestyle/marker type, see options of the PLOT function (default = '-')
%                       can be a single style for all datasets, or a cell-array containing one style for each dataset
%   cfg.linewidth     = linewidth in points (default = 0.5)
%   cfg.linecolor     = color(s) used for plotting the dataset(s) (default = 'brgkywrgbkywrgbkywrgbkyw')
%                       alternatively, colors can be specified as nx3 matrix of RGB values
%   cfg.directionality = '', 'inflow' or 'outflow' specifies for
%                       connectivity measures whether the inflow into a
%                       node, or the outflow from a node is plotted. The
%                       (default) behavior of this option depends on the dimor
%                       of the input data (see below).
%
% The following options for the scaling of the EEG, EOG, ECG, EMG, MEG and NIRS channels
% is optional and can be used to bring the absolute numbers of the different
% channel types in the same range (e.g. fT and uV). The channel types are determined
% from the input data using FT_CHANNELSELECTION.
%   cfg.eegscale       = number, scaling to apply to the EEG channels prior to display
%   cfg.eogscale       = number, scaling to apply to the EOG channels prior to display
%   cfg.ecgscale       = number, scaling to apply to the ECG channels prior to display
%   cfg.emgscale       = number, scaling to apply to the EMG channels prior to display
%   cfg.megscale       = number, scaling to apply to the MEG channels prior to display
%   cfg.gradscale      = number, scaling to apply to the MEG gradiometer channels prior to display (in addition to the cfg.megscale factor)
%   cfg.magscale       = number, scaling to apply to the MEG magnetometer channels prior to display (in addition to the cfg.megscale factor)
%   cfg.nirsscale      = number, scaling to apply to the NIRS channels prior to display
%   cfg.mychanscale    = number, scaling to apply to the channels specified in cfg.mychan
%   cfg.mychan         = Nx1 cell-array with selection of channels
%   cfg.chanscale      = Nx1 vector with scaling factors, one per channel specified in cfg.channel
%
% For the plotting of directional connectivity data the cfg.directionality
% option determines what is plotted. The default value and the supported
% functionality depend on the dimord of the input data. If the input data
% is of dimord 'chan_chan_XXX', the value of directionality determines
% whether, given the reference channel(s), the columns (inflow), or rows
% (outflow) are selected for plotting. In this situation the default is
% 'inflow'. Note that for undirected measures, inflow and outflow should
% give the same output. If the input data is of dimord 'chancmb_XXX', the
% value of directionality determines whether the rows in data.labelcmb are
% selected. With 'inflow' the rows are selected if the refchannel(s) occur in
% the right column, with 'outflow' the rows are selected if the
% refchannel(s) occur in the left column of the labelcmb-field. Default in
% this case is '', which means that all rows are selected in which the
% refchannel(s) occur. This is to robustly support linearly indexed
% undirected connectivity metrics. In the situation where undirected
% connectivity measures are linearly indexed, specifying 'inflow' or
% 'outflow' can result in unexpected behavior.
%
% to facilitate data-handling and distributed computing you can use
%   cfg.inputfile   =  ...
% if you specify this option the input data will be read from a *.mat
% file on disk. this mat files should contain only a single variable named 'data',
% corresponding to the input structure.
%
% See also FT_SINGLEPLOTTFR, FT_MULTIPLOTER, FT_MULTIPLOTTFR, FT_TOPOPLOTER, FT_TOPOPLOTTFR

% Undocumented local options:
% cfg.zlim/xparam (set to a specific frequency range or time range [zmax zmin] for an average
% over the frequency/time bins for TFR data.  Use in conjunction with e.g. xparam = 'time', and cfg.parameter = 'powspctrm').

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

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% DEVELOPERS NOTE: This code is organized in a similar fashion for multiplot/singleplot/topoplot
% and for ER/TFR and should remain consistent over those 6 functions.
% Section 1: general cfg handling that is independent from the data
% Section 2: data handling, this also includes converting bivariate (chan_chan and chancmb) into univariate data
% Section 3: select the data to be plotted and determine min/max range
% Section 4: do the actual plotting
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%% Section 1: general cfg handling that is independent from the data

% these are used by the ft_preamble/ft_postamble function and scripts
ft_revision = '$Id$';
ft_nargin   = nargin;
ft_nargout  = nargout;

% do the general setup of the function
ft_defaults
ft_preamble init
ft_preamble debug
ft_preamble loadvar varargin
ft_preamble provenance varargin
ft_preamble trackconfig

% the ft_abort variable is set to true or false in ft_preamble_init
if ft_abort
  return
end

Ndata = numel(varargin);
for i=1:Ndata
  % check if the input data is valid for this function
  varargin{i} = ft_checkdata(varargin{i}, 'datatype', {'timelock', 'freq'});
end

% check if the input cfg is valid for this function
cfg = ft_checkconfig(cfg, 'unused',     {'cohtargetchannel'});
cfg = ft_checkconfig(cfg, 'renamedval', {'zlim', 'absmax', 'maxabs'});
cfg = ft_checkconfig(cfg, 'renamedval', {'directionality', 'feedforward', 'outflow'});
cfg = ft_checkconfig(cfg, 'renamedval', {'directionality', 'feedback',    'inflow'});
cfg = ft_checkconfig(cfg, 'renamed',    {'matrixside',     'directionality'});
cfg = ft_checkconfig(cfg, 'renamed',    {'channelindex',   'channel'});
cfg = ft_checkconfig(cfg, 'renamed',    {'channelname',    'channel'});
cfg = ft_checkconfig(cfg, 'renamed',    {'cohrefchannel',  'refchannel'});
cfg = ft_checkconfig(cfg, 'renamed',	  {'zparam',         'parameter'});
cfg = ft_checkconfig(cfg, 'renamed',    {'graphcolor', 'linecolor'});
cfg = ft_checkconfig(cfg, 'deprecated', {'xparam'});

% set the defaults
cfg.viewmode        = ft_getopt(cfg, 'viewmode',      'average'); % average or butterfly
cfg.baseline        = ft_getopt(cfg, 'baseline',      'no');
cfg.trials          = ft_getopt(cfg, 'trials',        'all', 1);
cfg.xlim            = ft_getopt(cfg, 'xlim',          'maxmin');
cfg.ylim            = ft_getopt(cfg, 'ylim',          'maxmin');
cfg.zlim            = ft_getopt(cfg, 'zlim',          'maxmin');
cfg.comment         = ft_getopt(cfg, 'comment',        strcat([date '\n']));
cfg.axes            = ft_getopt(cfg, ' axes',         'yes');
cfg.fontsize        = ft_getopt(cfg, 'fontsize',       8);
cfg.hotkeys         = ft_getopt(cfg, 'hotkeys',       'yes');
cfg.interactive     = ft_getopt(cfg, 'interactive',   'yes');
cfg.renderer        = ft_getopt(cfg, 'renderer',       []); % let MATLAB decide on the default
cfg.maskparameter   = ft_getopt(cfg, 'maskparameter',  []);
cfg.colorgroups     = ft_getopt(cfg, 'colorgroups',   'condition'); % this is the only supported option
cfg.linecolor       = ft_getopt(cfg, 'linecolor',     'brgkywrgbkywrgbkywrgbkyw');
cfg.linestyle       = ft_getopt(cfg, 'linestyle',     '-');
cfg.linewidth       = ft_getopt(cfg, 'linewidth',      0.5);
cfg.maskstyle       = ft_getopt(cfg, 'maskstyle',     'box');
cfg.maskfacealpha   = ft_getopt(cfg, 'maskfacealpha', 1);
cfg.channel         = ft_getopt(cfg, 'channel',       'all');
cfg.title           = ft_getopt(cfg, 'title',          []);
cfg.directionality  = ft_getopt(cfg, 'directionality', []);
cfg.figurename      = ft_getopt(cfg, 'figurename',     []);
cfg.preproc         = ft_getopt(cfg, 'preproc',        []);
cfg.frequency       = ft_getopt(cfg, 'frequency',     'all'); % needed for frequency selection with TFR data
cfg.latency         = ft_getopt(cfg, 'latency',       'all'); % needed for latency selection with TFR data, FIXME, probably not used
cfg.showlegend      = ft_getopt(cfg, 'showlegend',    'no');

% check for linestyle being a cell-array
if ischar(cfg.linestyle)
  cfg.linestyle = repmat({cfg.linestyle}, 1, Ndata);
end
% check it's length, and lengthen it if does not have enough styles in it
if (length(cfg.linestyle) > 1) && (length(cfg.linestyle) < Ndata )
  ft_error('either specify cfg.linestyle as a cell-array with one cell for each dataset, or only specify one linestyle')
elseif (length(cfg.linestyle) == 1)
  cfg.linestyle = repmat(cfg.linestyle, 1, Ndata);
end

% this is needed for the figure title and correct labeling of linecolor later on
if isfield(cfg, 'dataname') && ~isempty(cfg.dataname)
  dataname = cfg.dataname;
elseif isfield(cfg, 'inputfile') && ~isempty(cfg.inputfile)
  dataname = cfg.inputfile;
elseif nargin>1
  dataname = arrayfun(@inputname, 2:nargin, 'UniformOutput', false);
else
  dataname = {};
end

%% Section 2: data handling, this also includes converting bivariate (chan_chan and chancmb) into univariate data

for i=1:Ndata
  dtype{i}    = ft_datatype(varargin{i});
  hastime(i)  = isfield(varargin{i}, 'time');
  hasfreq(i)  = isfield(varargin{i}, 'freq');
end

% check if the input has consistent datatypes
if ~all(strcmp(dtype, dtype{1})) || ~all(hastime==hastime(1)) || ~all(hasfreq==hasfreq(1))
  ft_error('different datatypes are not allowed as input');
else
  dtype   = dtype{1};
  hastime = hastime(1);
  hasfreq = hasfreq(1);
end

% Set x/y/parameter according to datatype and dimord
switch dtype
  case 'timelock'
    xparam = 'time';
    if isfield(varargin{1}, 'trial')
      cfg.parameter = ft_getopt(cfg, 'parameter', 'trial');
    elseif isfield(varargin{1}, 'individual')
      cfg.parameter = ft_getopt(cfg, 'parameter', 'individual');
    elseif isfield(varargin{1}, 'avg')
      cfg.parameter = ft_getopt(cfg, 'parameter', 'avg');
    end
  case 'freq'
    if hastime && hasfreq
      xparam = 'time';   % average over selected frequencies
    else
      xparam = 'freq';
    end
    cfg.parameter = ft_getopt(cfg, 'parameter', 'powspctrm');
  case 'comp'
    % not supported
  otherwise
    % not supported
end

% check whether rpt/subj is present and remove if necessary
dimord = getdimord(varargin{1}, cfg.parameter);
dimtok = tokenize(dimord, '_');
hasrpt = any(ismember(dimtok, {'rpt' 'subj'}));

if ~hasrpt
  assert(isequal(cfg.trials, 'all') || isequal(cfg.trials, 1), 'incorrect specification of cfg.trials for data without repetitions');
else
  assert(~isempty(cfg.trials), 'empty specification of cfg.trials for data with repetitions');
end

% parse cfg.channel
if isfield(cfg, 'channel') && isfield(varargin{1}, 'label')
  cfg.channel = ft_channelselection(cfg.channel, varargin{1}.label);
elseif isfield(cfg, 'channel') && isfield(varargin{1}, 'labelcmb')
  cfg.channel = ft_channelselection(cfg.channel, unique(varargin{1}.labelcmb(:)));
end

% apply baseline correction
if ~strcmp(cfg.baseline, 'no')
  tmpcfg = keepfields(cfg, {'baseline', 'baselinetype', 'baselinewindow', 'demean', 'parameter', 'channel'});
  for i=1:Ndata
    % keep mask-parameter if it is set
    if ~isempty(cfg.maskparameter)
      tempmask = varargin{i}.(cfg.maskparameter);
    end
    if strcmp(dtype, 'timelock') && strcmp(xparam, 'time')
      varargin{i} = ft_timelockbaseline(tmpcfg, varargin{i});
    elseif strcmp(dtype, 'freq') && strcmp(xparam, 'time')
      varargin{i} = ft_freqbaseline(tmpcfg, varargin{i});
    elseif strcmp(dtype, 'freq') && strcmp(xparam, 'freq')
      ft_error('baseline correction is not supported for spectra without a time dimension');
    else
      ft_warning('baseline correction not applied, please set xparam');
    end
    % put mask-parameter back if it is set
    if ~isempty(cfg.maskparameter)
      varargin{i}.(cfg.maskparameter) = tempmask;
    end
  end
end

% channels should NOT be selected and averaged here, since a topoplot might follow in interactive mode
tmpcfg = keepfields(cfg, {'showcallinfo', 'trials'});
if hasrpt
  tmpcfg.avgoverrpt = 'yes';
else
  tmpcfg.avgoverrpt = 'no';
end
if hastime && hasfreq
  tmpcfg.avgoverfreq = 'yes';       % average over selected frequencies
  tmpcfg.frequency = cfg.frequency; % not to be confused with cfg.xlim or cfg.ylim
  tmpcfg.keepfreqdim = 'no';
else
  tmpcfg.avgoverfreq = 'no';
end
tmpvar = varargin{1};
[varargin{:}] = ft_selectdata(tmpcfg, varargin{:});
% restore the provenance information and put back cfg.channel
tmpchannel  = cfg.channel;
[cfg, varargin{:}] = rollback_provenance(cfg, varargin{:});
cfg.channel = tmpchannel;

if isfield(tmpvar, cfg.maskparameter) && ~isfield(varargin{1}, cfg.maskparameter)
  % the mask parameter is not present after ft_selectdata, because it is
  % not included in all input arguments. Make the same selection and copy
  % it over
  tmpvar = ft_selectdata(tmpcfg, tmpvar);
  varargin{1}.(cfg.maskparameter) = tmpvar.(cfg.maskparameter);
end

clear tmpvar tmpcfg dimord dimtok hastime hasfreq hasrpt

% ensure that the preproc specific options are located in the cfg.preproc
% substructure, but also ensure that the field 'refchannel' remains at the
% highest level in the structure. This is a little hack by JM because the field
% refchannel can relate to connectivity or to an EEg reference.

if isfield(cfg, 'refchannel'), refchannelincfg = cfg.refchannel; cfg = rmfield(cfg, 'refchannel'); end
cfg = ft_checkconfig(cfg, 'createsubcfg',  {'preproc'});
if exist('refchannelincfg', 'var'), cfg.refchannel  = refchannelincfg; end

if ~isempty(cfg.preproc)
  % preprocess the data, i.e. apply filtering, baselinecorrection, etc.
  fprintf('applying preprocessing options\n');
  if ~isfield(cfg.preproc, 'feedback')
    cfg.preproc.feedback = cfg.interactive;
  end
  for i=1:Ndata
    varargin{i} = ft_preprocessing(cfg.preproc, varargin{i});
  end
end

% Handle the bivariate case
dimord = getdimord(varargin{1}, cfg.parameter);
if startsWith(dimord, 'chan_chan_') || startsWith(dimord, 'chancmb_')
  % convert the bivariate data to univariate and call this plotting function again
  cfg.originalfunction = 'ft_singleplotER';
  cfg.trials = 'all'; % trial selection has been taken care off
  bivariate_common(cfg, varargin{:});
  return
end

% Apply channel-type specific scaling
fn = fieldnames(cfg);
tmpcfg = keepfields(cfg, fn(endsWith(fn, 'scale') | startsWith(fn, 'mychan') | strcmp(fn, 'channel') | strcmp(fn, 'parameter')));
for i=1:Ndata
  varargin{i}= chanscale_common(tmpcfg, varargin{i});
end

%% Section 3: select the data to be plotted and determine min/max range

% Take the desided subselection of channels, this is the same in all datasets
[selchan] = match_str(varargin{1}.label, cfg.channel);

% Get physical min/max range of x, i.e. time or frequency
if strcmp(cfg.xlim, 'maxmin')
  % Find maxmin throughout all varargins:
  xmin = [];
  xmax = [];
  for i=1:Ndata
    xmin = min([xmin varargin{i}.(xparam)]);
    xmax = max([xmax varargin{i}.(xparam)]);
  end
else
  xmin = cfg.xlim(1);
  xmax = cfg.xlim(2);
end

% Get the index of the nearest bin, this is the same in all datasets
xminindx = nearest(varargin{1}.(xparam), xmin);
xmaxindx = nearest(varargin{1}.(xparam), xmax);
xmin = varargin{1}.(xparam)(xminindx);
xmax = varargin{1}.(xparam)(xmaxindx);
selx = xminindx:xmaxindx;
xval = varargin{1}.(xparam)(selx);

% get physical y-axis range, i.e. parameter to be plotted
if ~isnumeric(cfg.ylim)
  % find maxmin throughout all varargins
  ymin = [];
  ymax = [];
  for i=1:Ndata
    % select the channels in the data that match with the layout and that are selected for plotting
    switch cfg.viewmode
      case 'average'
        dat = nanmean(varargin{i}.(cfg.parameter)(selchan, selx), 1); % mean over channels
      case 'butterfly'
        dat = varargin{i}.(cfg.parameter)(selchan,selx);
    end
    ymin = nanmin([ymin nanmin(nanmin(nanmin(dat(:))))]);
    ymax = nanmax([ymax nanmax(nanmax(nanmax(dat(:))))]);
  end
  if strcmp(cfg.ylim, 'maxabs') % handle maxabs, make y-axis center on 0
    ymax = max([abs(ymax) abs(ymin)]);
    ymin = -ymax;
  elseif strcmp(cfg.ylim, 'zeromax')
    ymin = 0;
  elseif strcmp(cfg.ylim, 'minzero')
    ymax = 0;
  end
else
  ymin = cfg.ylim(1);
  ymax = cfg.ylim(2);
end

% gather the data from all input data structures
datamatrix = zeros(Ndata, length(selx));
for i=1:Ndata
  datamatrix(i,:) = mean(varargin{i}.(cfg.parameter)(selchan, selx), 1); % mean over channels
end

% gather the mask from the first input data structure
if ~isempty(cfg.maskparameter)
  maskmatrix = mean(varargin{1}.(cfg.maskparameter)(selchan, selx), 1); % mean over channels
else
  % create an Nx0 matrix
  maskmatrix = zeros(length(selchan), 0);
end

%% Section 4: do the actual plotting

% determine the coloring of channels/conditions
linecolor = linecolor_common(cfg, varargin{:});

cla
hold on

yval = datamatrix;
mask = maskmatrix;

if strcmp(cfg.maskstyle, 'difference')
  % combine the conditions in a single plot, highlight the difference
  ft_plot_vector(xval, yval, 'color', linecolor, 'linewidth', cfg.linewidth, 'style', cfg.linestyle{1}, 'linewidth', cfg.linewidth, 'highlight', mask, 'highlightstyle', cfg.maskstyle, 'hlim', [xmin xmax], 'vlim', [ymin ymax], 'facealpha', cfg.maskfacealpha);
else
  % loop over the conditions, plot them on top of each other
  for i=1:Ndata
    ft_plot_vector(xval, yval(i,:), 'color', linecolor(i,:), 'linewidth', cfg.linewidth, 'style', cfg.linestyle{i}, 'highlight', mask, 'highlightstyle', cfg.maskstyle, 'linewidth', cfg.linewidth, 'hlim', [xmin xmax], 'vlim', [ymin ymax], 'facealpha', cfg.maskfacealpha);
  end
end

if ischar(linecolor)
  set(gca, 'ColorOrder', char2rgb(linecolor))
elseif isnumeric(linecolor)
  set(gca, 'ColorOrder', linecolor)
end

% show the legend with the colors of the conditions
if istrue(cfg.showlegend) && Ndata>1
  if strcmp(cfg.maskstyle, 'difference')
    colorLabels = {'difference'};
  else
    colorLabels = {};
  end
  for i=1:Ndata
    if ischar(linecolor)
      colorLabels{end+1} = [dataname{i} '='         linecolor(i)     ];
    elseif isnumeric(linecolor)
      colorLabels{end+1} = [dataname{i} '=' num2str(linecolor(i, :)) ];
    end
  end
  legend(colorLabels)
end

% set xlim and ylim
if xmin~=xmax
  xlim([xmin xmax]);
end
if ymin~=ymax
  ylim([ymin ymax]);
end

% adjust mask box extents to ymin/ymax
if ~isempty(cfg.maskparameter)
  ptchs = findobj(gcf, 'type', 'patch');
  for i = 1:length(ptchs)
    YData = get(ptchs(i), 'YData');
    YData(YData == min(YData)) = ymin;
    YData(YData == max(YData)) = ymax;
    set(ptchs(i), 'YData',YData);
  end
end

% Set callback to adjust axes
if strcmp('yes', cfg.hotkeys)
  %  attach data and cfg to figure and attach a key listener to the figure
  set(gcf, 'KeyPressFcn', {@key_sub, xmin, xmax, ymin, ymax})
end

% Create axis title containing channel name(s) and channel number(s):
if ~isempty(cfg.title)
  t = cfg.title;
else
  if length(cfg.channel) == 1
    t = [char(cfg.channel) ' / ' num2str(selchan) ];
  else
    t = sprintf('mean(%0s)', join_str(', ', cfg.channel));
  end
end
title(t, 'fontsize', cfg.fontsize);

% set the figure window title, add channel labels if number is small
if isempty(get(gcf, 'Name'))
  if length(selchan) < 5
    chans = join_str(', ', cfg.channel);
  else
    chans = '<multiple channels>';
  end
  if ~isempty(cfg.figurename)
    set(gcf, 'name', cfg.figurename);
    set(gcf, 'NumberTitle', 'off');
  elseif ~isempty(dataname)
    set(gcf, 'Name', sprintf('%d: %s: %s (%s)', double(gcf), mfilename, join_str(', ', dataname), chans));
    set(gcf, 'NumberTitle', 'off');
  else
    set(gcf, 'Name', sprintf('%d: %s (%s)', double(gcf), mfilename, chans));
    set(gcf, 'NumberTitle', 'off');
  end
end

% make the figure interactive
if strcmp(cfg.interactive, 'yes')
  % add the cfg/data/channel information to the figure under identifier linked to this axis
  ident                  = ['axh' num2str(round(sum(clock.*1e6)))]; % unique identifier for this axis
  set(gca, 'tag', ident);
  info                   = guidata(gcf);
  info.(ident).cfg       = cfg;
  info.(ident).varargin  = varargin;
  info.(ident).dataname  = dataname;
  guidata(gcf, info);
  set(gcf, 'windowbuttonupfcn',     {@ft_select_range, 'multiple', false, 'yrange', false, 'callback', {@select_topoplotER}, 'event', 'windowbuttonupfcn'});
  set(gcf, 'windowbuttondownfcn',   {@ft_select_range, 'multiple', false, 'yrange', false, 'callback', {@select_topoplotER}, 'event', 'windowbuttondownfcn'});
  set(gcf, 'windowbuttonmotionfcn', {@ft_select_range, 'multiple', false, 'yrange', false, 'callback', {@select_topoplotER}, 'event', 'windowbuttonmotionfcn'});
end

% set renderer if specified
if ~isempty(cfg.renderer)
  set(gcf, 'renderer', cfg.renderer)
end

hold off

% do the general cleanup and bookkeeping at the end of the function
ft_postamble debug
ft_postamble trackconfig
ft_postamble previous varargin
ft_postamble provenance
ft_postamble savefig

% add a menu to the figure, but only if the current figure does not have subplots
menu_fieldtrip(gcf, cfg, false);

if ~ft_nargout
  % don't return anything
  clear cfg
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% SUBFUNCTION which is called after selecting a time range
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function select_topoplotER(range, varargin)

% fetch cfg/data based on axis indentifier given as tag
ident    = get(gca, 'tag');
info     = guidata(gcf);
cfg      = info.(ident).cfg;
varargin = info.(ident).varargin;
if ~isempty(range)
  cfg = removefields(cfg, 'inputfile');   % the reading has already been done and varargin contains the data
  cfg = removefields(cfg, 'showlabels');  % this is not allowed in topoplotER
  cfg.baseline = 'no';                    % make sure the next function does not apply a baseline correction again
  cfg.dataname = info.(ident).dataname;   % put data name in here, this cannot be resolved by other means
  cfg.channel = 'all';                    % make sure the topo displays all channels, not just the ones in this singleplot
  cfg.comment = 'auto';
  cfg.trials = 'all';                     % trial selection has already been taken care of
  cfg.xlim = range(1:2);
  % if user specified a ylim, copy it over to the zlim of topoplot
  if isfield(cfg, 'ylim')
    cfg.zlim = cfg.ylim;
    cfg = rmfield(cfg, 'ylim');
  end
  fprintf('selected cfg.xlim = [%f %f]\n', cfg.xlim(1), cfg.xlim(2));
  % ensure that the new figure appears at the same position
  f = figure('position', get(gcf, 'Position'));
  ft_topoplotER(cfg, varargin{:});
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% SUBFUNCTION which handles hot keys in the current plot
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function key_sub(handle, eventdata, varargin)
xlimits = xlim;
ylimits = ylim;
incr_x = abs(xlimits(2) - xlimits(1)) /10;
incr_y = abs(ylimits(2) - ylimits(1)) /10;

if length(eventdata.Modifier) == 1 && strcmp(eventdata.Modifier{:}, 'control')
  % TRANSLATE by 10%
  switch eventdata.Key
    case 'leftarrow'
      xlim([xlimits(1)+incr_x xlimits(2)+incr_x])
    case 'rightarrow'
      xlim([xlimits(1)-incr_x xlimits(2)-incr_x])
    case 'uparrow'
      ylim([ylimits(1)-incr_y ylimits(2)-incr_y])
    case 'downarrow'
      ylim([ylimits(1)+incr_y ylimits(2)+incr_y])
  end % switch
else
  % ZOOM by 10%
  switch eventdata.Key
    case 'leftarrow'
      xlim([xlimits(1)-incr_x xlimits(2)+incr_x])
    case 'rightarrow'
      xlim([xlimits(1)+incr_x xlimits(2)-incr_x])
    case 'uparrow'
      ylim([ylimits(1)-incr_y ylimits(2)+incr_y])
    case 'downarrow'
      ylim([ylimits(1)+incr_y ylimits(2)-incr_y])
    case 'm'
      xlim([varargin{1} varargin{2}])
      ylim([varargin{3} varargin{4}])
  end % switch
end % if

author: spm


text:
function T = convert_CPD_to_table_hidden_ps(CPD, self_val)
% CONVERT_CPD_TO_TABLE_HIDDEN_PS Convert a Gaussian CPD to a table
% function T = convert_CPD_to_table_hidden_ps(CPD, self_val)
%
% self_val must be a non-empty vector.
% All the parents are hidden.
%
% This is used by misc/convert_dbn_CPDs_to_tables

m = CPD.mean;
C = CPD.cov;
W = CPD.weights;

[ssz dpsize] = size(m);

T = zeros(dpsize, 1);
for i=1:dpsize
  T(i) = gaussian_prob(self_val, m(:,i), C(:,:,i));
end

author: bayesnet

repo: bnt

save_path: github-repos/MATLAB/bayesnet-bnt


text:
%% febio
% Below is a demonstration of the features of the |febio| function

%%
clear; close all; clc;

%% Syntax
% |febio(varargin);|

%% Description
% UNDOCUMENTED
%% Examples
%
%%
%
% <<gibbVerySmall.gif>>
%
% _*GIBBON*_
% <www.gibboncode.org>
%
% _Kevin Mattheus Moerman_, <gibbon.toolbox@gmail.com>
%%
% _*GIBBON footer text*_
%
% License: <https://github.com/gibbonCode/GIBBON/blob/master/LICENSE>
%
% GIBBON: The Geometry and Image-based Bioengineering add-On. A toolbox for
% image segmentation, image-based modeling, meshing, and finite element
% analysis.
%
% Copyright (C) 2006-2022 Kevin Mattheus Moerman and the GIBBON contributors
%
% This program is free software: you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation, either version 3 of the License, or
% (at your option) any later version.
%
% This program is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
% GNU General Public License for more details.
%
% You should have received a copy of the GNU General Public License
% along with this program.  If not, see <http://www.gnu.org/licenses/>.

author: gibbonCode

repo: GIBBON

save_path: github-repos/MATLAB/gibbonCode-GIBBON


text:
data = caffe_solver.net.blob_vec(caffe_solver.net.name2blob_index('data')).get_data();
rois = squeeze(caffe_solver.net.blob_vec(caffe_solver.net.name2blob_index('rois')).get_data());
labels = squeeze(caffe_solver.net.blob_vec(caffe_solver.net.name2blob_index('labels')).get_data());
bbox_targets = squeeze(caffe_solver.net.blob_vec(caffe_solver.net.name2blob_index('bbox_targets')).get_data());
bbox_loss_weights = squeeze(caffe_solver.net.blob_vec(caffe_solver.net.name2blob_index('bbox_loss_weights')).get_data());
cls_score = squeeze(caffe_solver.net.blob_vec(caffe_solver.net.name2blob_index('cls_score')).get_data());

scores = squeeze(caffe_solver.net.blobs('cls_score').get_data());
classes = [{'bg'} imdb_train{1}.classes];
figure(1);
for k = 1:size(data,4)

  img = permute((data(:,:,:,k)+128)/256, [2 1 3 4 5 6]);
  img_ = img(:,:, [3 2 1]);
inds = rois(1,:)==k-1;
  boxes = rois(2:end,inds)';

  bbox_reg_boxes = find( bbox_loss_weights(5,inds) > 0 ) ;
  bbox_reg_weights = bbox_targets(:,bbox_reg_boxes) ;

  boxes_cell = {boxes(bbox_reg_boxes,:)};

  subplot(1,size(data,4),k);
  try
    track_deltas = caffe_solver.net.blob_vec(caffe_solver.net.name2blob_index('bbox_disp')).get_data();
    track_boxes = caffe_solver.net.blob_vec(caffe_solver.net.name2blob_index('rois_disp')).get_data();
    track_boxes = squeeze(track_boxes)' +1;
    pred_boxes = rfcn_bbox_transform_inv(track_boxes(:,2:end), squeeze(track_deltas)');
    pred_boxes = pred_boxes(:, 5:end);
    if k==track_boxes(1),   boxes_cell = {track_boxes(:,2:end)};
    else
       boxes_cell = {pred_boxes};
    end
  catch
  end
  this_scrs = scores(:,inds);
  [votes,predictions] = sort(this_scrs(:,bbox_reg_boxes), 1, 'descend') ;
        top5pred =  squeeze(predictions(1:5,:));
        highest_class_names = strrep(classes(top5pred), '_', '\_');

  if size(img,3) > 3
    subplot(2,2,1);showboxes(img_, boxes_cell);
    subplot(2,2,2); image(img(:,:,[4 : 6],:)); axis image;
    if size(img,3) > 10
      subplot(2,2,3); image(img(:,:,[7 : 9],:)); axis image;
      subplot(2,2,4); image(img(:,:,[9 : 11],:)); axis image;
    end
  else
    showboxes(img_, boxes_cell);
  end
  title(['top5 classes=' highest_class_names{1:5} ]);

end

author: feichtenhofer

repo: Detect-Track


text:
function calpak_test56 ( )

%*****************************************************************************80
%
%% TEST56 tests YEAR_IS_LEAP_JULIAN.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    13 April 2013
%
%  Author:
%
%    John Burkardt
%
  fprintf ( 1, '\n' );
  fprintf ( 1, 'CALPAK_TEST56\n' );
  fprintf ( 1, '  For the Julian calendar:\n' );
  fprintf ( 1, '  YEAR_IS_LEAP_JULIAN reports leap years.\n' );
  fprintf ( 1, '\n' );
  fprintf ( 1, '  Year  Leap?\n' );
  fprintf ( 1, '\n' );

  for y = 1990 : 2000
    sy = y_to_s_julian ( y );
    fprintf ( 1, '  %10s  %d\n', sy, year_is_leap_julian ( y ) );
  end

  return
end

author: johannesgerer

repo: jburkardt-m


text:
function exact = p01_exact ( )

%*****************************************************************************80
%
%% P01_EXACT returns the exact integral for problem 1.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    04 November 2009
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Output, real EXACT, the value of the integral.
%
  exact = exp ( 1.0 ) - 1.0;

  return
end


text:
% plot real and imaginary parts of various finite difference derivative operators
% input: c = scalar multiplier
% input: M = number of points on interval (excluding one end point)
% input: dx = uniform grid point separation

function plotmultiplies(c, M, dx)

m=0:(M-1);

thetam = 2*pi*m/M;

deltaR = c*(exp(i*thetam)-1)/dx;
deltaL = c*(1-exp(-i*thetam))/dx;
deltaC = i*c*sin(thetam)/dx;
deltaC4 = i*c*(8*sin(thetam)-sin(2*thetam))/(6*dx);
deltaC6 = i*c*sin(thetam).*(1+(2/3)*sin(thetam/2).^2 +(8/15)*sin(thetam/2).^4)/dx;

% plot real part of eigenvalues
subplot(1,2,1);
plot(thetam, real(deltaR), 'r-d');
hold on;
plot(thetam, real(deltaL), 'b-*');
plot(thetam, real(deltaC), 'g-s');
plot(thetam, real(deltaC4), 'k-o');
hold off;
xlabel('\theta_m = 2\pim/M', 'FontSize', 18);
title('Real part of ODE multipliers', 'FontSize', 18);
legend('\delta_+', '\delta_-', '\delta_0', '\delta_0(1-(1/6)dx^2\delta^2)');

% plot imaginary part of eigenvalues
subplot(1,2,2);
plot(thetam, imag(deltaR), 'r-d');
hold on;
plot(thetam, imag(deltaL), 'b-*');
plot(thetam, imag(deltaC), 'g-s');
plot(thetam, imag(deltaC4), 'k-o');
plot(thetam, imag(deltaC6), 'm-h');
plot([0 pi], [0 pi], 'k-', 'LineWidth', 2);      % first part of correct phase relationship
plot([pi 2*pi], [-pi 0], 'k-', 'LineWidth', 2); % second part

hold off;
xlabel('\theta_m = 2\pim/M', 'FontSize', 18);
title('Imaginary part of ODE multipliers', 'FontSize', 18);
legend('\delta_+', '\delta_-', '\delta_0', '\delta_0(1-(1/6)dx^2\delta^2)', '6th order central difference', 'Exact');

repo: Aero-matlab


text:
function chebyshev_polynomial_test05 ( )

%*****************************************************************************80
%
%% CHEBYSHEV_POLYNOMIAL_TEST05 tests T_QUADRATURE_RULE.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    29 March 2012
%
%  Author:
%
%    John Burkardt
%
  fprintf ( 1, '\n' );
  fprintf ( 1, 'CHEBYSHEV_POLYNOMIAL_TEST05:\n' );
  fprintf ( 1, '  T_QUADRATURE_RULE computes the quadrature rule\n' );
  fprintf ( 1, '  associated with T(n,x);\n' );

  n = 7;
  [ x, w ] = t_quadrature_rule ( n );

  r8vec2_print ( n, x, w, '      X            W' );

  fprintf ( 1, '\n' );
  fprintf ( 1, '  Use the quadrature rule to estimate:\n' );
  fprintf ( 1, '\n' );
  fprintf ( 1, '    Q = Integral ( -1 <= X <= +1 ) X^E / sqrt ( 1-x^2) dx\n' );
  fprintf ( 1, '\n' );
  fprintf ( 1, '   E       Q_Estimate      Q_Exact\n' );
  fprintf ( 1, '\n' );

  for e = 0 : 2 * n - 1
    if ( e == 0 )
      f = ones ( n, 1 );
    else
      f(1:n) = x(1:n).^e;
    end
    q = w' * f;
    q_exact = t_integral ( e );
    fprintf ( 1, '  %2d  %14g  %14g\n', e, q, q_exact );
  end

  return
end


text:
function evaluate(resPath, wwwPath, class, width, height)

if strcmp(class, 'chair') || strcmp(class, 'bed')
  param_length = [5, 6, 12, 14];
elseif strcmp(class, 'swivelchair') || strcmp(class, 'sofa')
  param_length = [7, 8, 14, 16];
end

load(resPath, 'outputs');
numInst = size(outputs, 1);
alphaPred = outputs(:, 1 : param_length(1));
finvPred = outputs(:, (param_length(1) + 1) : param_length(2));
sincosThetaPred = outputs(:, (param_length(2) + 1) : param_length(3));
sincosThetaPred(:, 3) = 0;
sincosThetaPred(:, 6) = 1;
tranPred = outputs(:, (param_length(3) + 1) : param_length(4));

addpath(fullfile('3D', 'genSynData'));
addpath(fullfile('3D', 'tools'));
stickStruct = getStickFigure('class', class);

for i = 1 : numInst
  im = visualize3Dpara(alphaPred(i, :)', sincosThetaPred(i, :)', ...
    tranPred(i, :)', 1 ./ finvPred(i, :), stickStruct.baseShape{1}, ...
    stickStruct.edgeAdj{1}, 'h', height, 'w', width, 'lineWidth', 6);
  imwrite(im, fullfile(wwwPath, sprintf('%08d.jpg', i)));
end

end

author: jiajunwu


text:
classdef ConfusionMatrix < handle
    properties (GetAccess = public)
        confusionMatrixData;
        nClasses;
        classNames;
        containsNullClass;
    end

    methods (Access = public)
        function obj = ConfusionMatrix(truthClasses,predictedClasses,classes,classNames)
            if nargin > 1
                if nargin == 3
                    classNames = [];
                end

                if(any(predictedClasses == Labeling.kNullClass))
                    obj.confusionMatrixData = confusionmat(truthClasses,predictedClasses,'Order',[classes 0]);
                    obj.classNames = [classNames, Labeling.kNullClassStr];
                    obj.containsNullClass = true;
                else
                    obj.confusionMatrixData = confusionmat(truthClasses,predictedClasses,'Order',classes);
                    obj.classNames = classNames;
                    obj.containsNullClass = false;
                end

                obj.nClasses = size(obj.confusionMatrixData,1);
            end
        end
    end

    methods (Access = public, Static)
        function confusionMatrix = CreateConfusionMatrixWithData(confusionMatrixData)
            confusionMatrix = ConfusionMatrix();
            confusionMatrix.confusionMatrixData = confusionMatrixData;
        end
    end

end

author: avenix

repo: WDK

save_path: github-repos/MATLAB/avenix-WDK


text:
function sevec = figresize(bw2)
% This function take the cropped binary image and resize it to 5 x 7
% char representation as single vector.

bw_7050=imresize(bw2,[70,50]);
for cnt=1:7
    for cnt2=1:5
        Atemp=sum(bw_7050((cnt*10-9:cnt*10),(cnt2*10-9:cnt2*10)));
        sevec((cnt-1)*5+cnt2)=sum(Atemp);
    end
end

sevec=((100-sevec)/100);
sevec=sevec';

author: Sable

repo: mcbench-benchmarks


text:
function lpp_to_polynomial_test ( )

%*****************************************************************************80
%
%% LPP_TO_POLYNOMIAL_TEST tests LPP_TO_POLYNOMIAL.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    09 September 2014
%
%  Author:
%
%    John Burkardt
%
  fprintf ( 1, '\n' );
  fprintf ( 1, 'LPP_TO_POLYNOMIAL_TEST:\n' );
  fprintf ( 1, '  LPP_TO_POLYNOMIAL is given a Legendre product polynomial\n' );
  fprintf ( 1, '  and determines its polynomial representation.\n' );

  m = 2;
  fprintf ( 1, '\n' );
  fprintf ( 1, '  Using spatial dimension M = %d:\n', m );

  for rank = 1 : 11

    l = comp_unrank_grlex ( m, rank );

    o_max = prod ( floor ( ( l(1:m) + 2 ) / 2 ) );
    [ o, c, e ] = lpp_to_polynomial ( m, l, o_max );

    label = sprintf ( '  LPP #%d = L(%d,X)*L(%d,Y) =', rank, l(1), l(2) );

    fprintf ( 1, '\n' );
    polynomial_print ( m, o, c, e, label );

  end

  return
end


text:
% ir_denoise_test_sep1
% Examine denoising with ordinary vs separable finite differences
% i.e. |Ch X|_1 + |Cv X|_1 vs |Ch X Cv|_1
% Uses ADMM algorithm to optimize denoising cost function.

if ~isvar('yi'), printm 'generate noisy image'
	xt = ellipse_im(2^7);
	xt = zeros(2^6); xt(20:60,20:48) = 9;
	sig = 2^0; clim = [0 9];

	snr = @(x) -20*log(norm(x(:)-xt(:)) / norm(xt(:)));

	rng(7)
	yi = xt + sig * randn(size(xt));

	im plc 2 4
	im(1, xt, clim, 'True')
	im(2, yi, clim, 'Noisy')
	xlabelf('SNR = %4.1f dB', snr(yi))
end

if ~isvar('Cs'), printm 'C2 and Cs'
	arg.Ch = Cdiffs(size(xt,1), 'type_diff', 'circshift', 'offsets', 1);
	arg.Cv = Cdiffs(size(xt,2), 'type_diff', 'circshift', 'offsets', 1);
	Cs = fatrix2('arg', arg, 'idim', size(xt), 'odim', size(xt), ...
		'forw', @(arg, x) arg.Ch * x * arg.Cv, ... % separable!
		'back', @(arg, y) arg.Ch' * y * arg.Cv');
	C2 = Cdiffs(size(xt), 'type_diff', 'circshift', 'offsets', '2d:hv');

	im(3, 'col', 1, C2 * xt, 'Usual'), cbar
	im(4, Cs * xt, 'Separable'), cbar
	drawnow
end

if ~isvar('x2'), printm 'run admm with usual 2d finite differences'
	[x2s snr2] = ir_denoise_admm1(yi, 'beta', 2^+1, ...
		'stop_diff_tol', 1e-4, 'chat', 0, ...
		'isave', 'all', ...
		'userfun', @(x, iter) snr(x), ...
		'niter', 2^5, 'shrink', [], 'rho', 1);
	x2 = x2s(:,:,end);
end

if ~isvar('x1'), printm 'run admm with separable finite differences'
	[x1s snr1] = ir_denoise_admm1(yi, 'beta', 2^+1, ...
		'stop_diff_tol', 1e-4, 'chat', 0, ...
		'C', 'sep1', ...
		'isave', 'all', ...
		'userfun', @(x, iter) snr(x), ...
		'niter', 2^5, 'shrink', [], 'rho', 1);
	x1 = x1s(:,:,end);
end

if 1 % figures
	im(7, x2, clim, 'Denoised usual')
	xlabelf('SNR = %4.1f dB', snr(x2))
	im(8, x1, clim, 'Denoised separ.')
	xlabelf('SNR = %4.1f dB', snr(x1))

%	im subplot 5
	subplot(223)
	iters = 0:numel(snr1)-1;
	plot(iters, snr2, '-o', iters, snr1, '-x')
	snr_lim = [floor(min([snr1; snr2])) ceil(max([snr1; snr2]))];
	axis([0 numel(iters) snr_lim])
	ytick(snr_lim)
	xlabelf 'ADMM iteration'
	ylabelf 'SNR (dB)'
	legend('usual', 'separable', 'location', 'southeast')
end

author: JeffFessler

repo: mirt

save_path: github-repos/MATLAB/JeffFessler-mirt


text:
function boolean  = isobj(handle,currobj)
    % ISOBJ	True for objects
    %       ISOBJ(X) returns 1's where the elements of X are
    % 		valid handles and 0's where they are not.
    %
    % Keith Rogers 11/30/93
    if ~isempty(handle) ; currobj = 0; boolean = 0;return; end
    if nargin == 1
        currobj = 0;
    end
    boolean = (handle == currobj);
    if(strcmp(get(currobj,'Type'),'axes'));
        boolean = boolean | (handle == get(currobj,'xlabel'));
        boolean = boolean | (handle == get(currobj,'ylabel'));
        boolean = boolean | (handle == get(currobj,'zlabel'));
        boolean = boolean | (handle == get(currobj,'title'));
    end
    children = get(currobj,'Children');
    if(children)
        for i = 1:length(children)
            boolean = boolean | isobj(handle,children(i));
        end
    end

author: CelsoReyes


text:
%  generate the Companion matrix of the polynomials f
%
%  Syntax:   >> C = CompanionMatrix(f)
%
%  Input:      f --- (string or coefficient vector) polynomial
%
%  Output:     C --- (matrix) The companion matrix
%
%  Example:  >> CompanionMatrix('x^4 + 2*x^3 + 3*x^2 + 4*x+5')
%
%            ans =
%
%                -2    -3    -4    -5
%                 1     0     0     0
%                 0     1     0     0
%                 0     0     1     0

author: zarathustr

repo: LibQPEP

save_path: github-repos/MATLAB/zarathustr-LibQPEP


text:
function y = Ackley(x)
  n = 2;
  sum1 = 0;
  sum2 = 0;

  for i=1:n
   sum1 = s1+x(i)^2;
   sum2 = s2+cos((2*pi)*x(i));
  end

  y = 20 + exp(1)-20*exp(-0.2*sqrt(1/n*sum1))-exp(1/n*sum2);

end


text:
function varargout = mapc(values, fcns, varargin)

% varargout = mapc(value, fcns, varargin)
%
% See functional_programming_examples.m for a discussion of this function.
%
% Pass all of the inputs in |values| (a cell array of arguments) to each
% function handle in the cell array |fcns|, allowing output size to vary.
% The output is therefore a cell array instead of an array. Any additional
% arguments are passed to the |cellfun| function, such as the
% 'ErrorHandler' property.
%
% >> mapc({1}, {@(x) 2*x, @asin, @(x) char(x+[115 110 108 96 115 110])})
% ans =
%     [2]    [1.5708]    'tomato'
%
% Tucker McClure
% Copyright 2013 The MathWorks, Inc.

    [varargout{1:nargout}] = cellfun(@(f) f(values{:}), fcns, ...
                                     'UniformOutput', false, varargin{:});

end


text:
function modelConstrained = addCouplingConstraint(model, rxnList, c, d, ineqSense)
% Constrains one (weighted) sum of fluxes to be above a lower bound.
% Appends to existing inequality constraints if they are present
%
% USAGE:
%
%    modelConstrained = constrainRxnListAboveBound(model, rxnList, c, d, ineqSense)
%
% INPUTS:
%    model:               model structure
%    rxnList:             cell array of reaction names
%
% OPTIONAL INPUTS:
%    c:                   `k x 1` vector :math:`c*v \geq d`
%    d:                   `n x 1` vector :math:`c*v \geq d`
%    ineqSense:           `k x 1` inequality sense {'L','G'}
%
% OUTPUT:
%    modelConstrained:    constrained model:
%
%                           * S - Stoichiometric matrix
%                           * b - Right hand side = dx/dt
%                           * C - Inequality constraint matrix
%                           * d - Inequality constraint right hand side
%                             :math:`[S; C] * v {=, \leq, \geq } [dxdt, d]`

% EXAMPLE:
%
%    rxnList = {'PCHOLP_hs_f', 'PLA2_2_f', 'SMS_f','PCHOLP_hs_b', 'PLA2_2_b', 'SMS_b'};
%    c = [1, 1, 1, 1, 1, 1];
%    d = 10;
%    ineqSense = 'G';
%    modelConstrained = constrainRxnListAboveBound(modelIrrev, rxnList, C, d, ineqSense);

modelConstrained = addCOBRAConstraints(model, rxnList ,d ,'c',c,'dsense',ineqSense);

author: opencobra

repo: cobratoolbox


text:
function net = rbfunpak(net, w)
%RBFUNPAK Separates a vector of RBF weights into its components.
%
%	Description
%	NET = RBFUNPAK(NET, W) takes an RBF network data structure NET and  a
%	weight vector W, and returns a network data structure identical to
%	the input network, except that the centres C, the widths WI, the
%	second-layer weight matrix W2 and the second-layer bias vector B2
%	have all been set to the corresponding elements of W.
%
%	See also
%	RBFPAK, RBF
%

%	Copyright (c) Ian T Nabney (1996-2001)

% Check arguments for consistency
errstring = consist(net, 'rbf');
if ~errstring
  error(errstring);
end

if net.nwts ~= length(w)
  error('Invalid length of weight vector')
end

nin 	= net.nin;
nhidden = net.nhidden;
nout 	= net.nout;

mark1 = nin*nhidden;
net.c = reshape(w(1:mark1), nhidden, nin);
if strcmp(net.actfn, 'gaussian')
  mark2 = mark1 + nhidden;
  net.wi = reshape(w(mark1+1:mark2), 1, nhidden);
else
  mark2 = mark1;
  net.wi = [];
end
mark3 = mark2 + nhidden*nout;
net.w2 = reshape(w(mark2+1:mark3), nhidden, nout);
mark4 = mark3 + nout;
net.b2 = reshape(w(mark3+1:mark4), 1, nout);


text:
function classType = mrGrayCheckClassType(fileName)
%
% classType = mrGrayCheckClassType(fileName)
%
% returns 'g' for gray file, 'c' for mrGray class file, 'n' for new
% ITKGray, 'u' for unknown.
% NIFTI class file.
%

[~,~,e] = fileparts(fileName);

switch(lower(e))
    case {'.gray','.grey'}
        classType = 'mrGray';
    case '.class'
        classType = 'mrGray';
    case {'.gz','.nii'}
        classType = 'n';
    otherwise
        classType = 'u';
end

if strcmpi(classType, 'mrGray')
    error('Vistasoft no longer supports mrGray Gray and Class files. Please convert to nifti.')
end

return;

author: vistalab

repo: vistasoft

save_path: github-repos/MATLAB/vistalab-vistasoft


text:
function x = forloop(x, n, f, cleanup_fcn)

% x = forloop(x, n, f, cleanup_fcn)
%
% Loop a given number of times. f should be a function handle taking in
% x and k, the iteration number (which starts at 1), and returning the
% updated x. If x is a cell array, each cell is handed to f as in
% f(x{:}), which allows f to have the form @(x1, x2,..., k).
%
% forloop(1,                        % Start with x = 1.
%         9,                        % Loop 9 times.
%         @(x, k) last(fprintf('Factorial(%d): %d\n', k, k*x), ...
%                      k*x))        % Print the result and return k*x.
%
% Tucker McClure
% Copyright 2013 The MathWorks, Inc.

    for k = 1:n
        if iscell(x)
            x = f(x{:}, k);
        else
            x = f(x, k);
        end
    end

    if nargin == 4
        if iscell(x)
            x = cleanup_fcn(x{:});
        else
            x = cleanup_fcn(x);
        end
    end

end


text:
cutoff = 50;
aspectr = 1;
original = 0;

Ds = smooth3(D);

if original
hiso = patch(isosurface(Ds,cutoff),...
    'FaceColor',[1,.75,.65],...
    'EdgeColor','none');

hcap = patch(isocaps(D,cutoff),...
    'FaceColor','interp',...
    'EdgeColor','none');
else
[hiso,hcap] = imagePatch(X,Y,Z,D,[nan coords(1) nan nan nan 40]);
end

view(45,30)
axis tight
daspect([1,1,aspectr])

myLight = lightangle(45,30);
set(gcf,'Renderer','zbuffer'); lighting phong
isonormals(Ds,hiso)
set(hcap,'AmbientStrength',.6)
set(hiso,'SpecularColorReflectance',0,'SpecularExponent',50)

% this is the key to avoiding dark surface head.
[az,el] = view
lightangle(az,el)

author: canlab

repo: CanlabCore

save_path: github-repos/MATLAB/canlab-CanlabCore


text:
function [c, s, b, g, smin, active_set] = thresholded_oasisAR2(y, g, sn, optimize_b,...
    optimize_g, decimate, maxIter, thresh_factor)
%% Infer the most likely discretized spike train underlying an AR(1) fluorescence trace
% Solves the sparse non-negative deconvolution problem
%  min 1/2|c-y|^2 + lam |s|_1 subject to s_t = c_t-g c_{t-1} >=s_min or =0

%% inputs:
%   y:  T*1 vector, One dimensional array containing the fluorescence intensities
%withone entry per time-bin.
%   g:  scalar, Parameter of the AR(1) process that models the fluorescence ...
%impulse response.
%   sn:  scalar, standard deviation of the noise distribution
%   optimize_b: bool, optimize baseline if True
%   optimize_g: integer, number of large, isolated events to consider for
%       optimizing g
%   decimate: int, decimation factor for estimating hyper-parameters faster
%       on decimated data
%   maxIter:  int, maximum number of iterations
%   active_set: npool x 4 matrix, warm stared active sets
%  thresh_factor: scalar, set the maximum thresh as thresh_factor*sn^2*T

%% outputs
%   c: T*1 vector, the inferred denoised fluorescence signal at each time-bin.
%   s: T*1 vector, discetized deconvolved neural activity (spikes)
%   b: scalar, fluorescence baseline
%   g: scalar, parameter of the AR(1) process
%   smin: scalar, minimum nonzero spike count
%   active_set: npool x 4 matrix, active sets

%% Authors: Pengcheng Zhou, Carnegie Mellon University, 2016
% ported from the Python implementation from Johannes Friedrich

%% References
% Friedrich J et.al., NIPS 2016, Fast Active Set Method for Online Spike Inference from Calcium Imaging

%% input arguments
y = reshape(y, [], 1);
T = length(y);

if ~exist('g', 'var') || isempty(g)
    g = estimate_time_constant(y, 2);
end
if ~exist('sn', 'var') || isempty(sn)
    sn = GetSn(y);
end
if ~exist('optimize_b', 'var') || isempty(optimize_b)
    optimize_b = false;
end
if ~exist('optimize_g', 'var') || isempty(optimize_g)
    optimize_g = 0;
end
if ~exist('decimate', 'var') || isempty(decimate)
    decimate = 1;
else
    decimate = max(1, round(decimate));
end
if ~exist('maxIter', 'var') || isempty(maxIter)
    maxIter = 10;
end
if ~exist('thresh_factor', 'var') || isempty(thresh_factor)
    thresh_factor = 1.0;
end

thresh = thresh_factor* sn * sn * T;
smin = 0;

% change parameters due to downsampling
if decimate>1
    decimate = 1;  %#ok<NASGU>
    disp('to be done');
    %     fluo = y;
    %     y = resample(y, 1, decimate);
    %     g = g^decimate;
    %     thresh = thresh / decimate / decimate;
    %     T = length(y);
end
g_converged = false;

%% optimize parameters
tol = 1e-4;
if ~optimize_b   %% don't optimize the baseline b
    %% initialization
    b = 0;
    [solution, spks, active_set] = oasisAR2(y, g, [], smin);

    %% iteratively update parameters lambda & g
    for miter=1:maxIter
        len_active_set = size(active_set, 1);

        res = y - solution;
        RSS = res' * res;

        % update g
        if and(optimize_g, ~g_converged);
            g0 = g;
            [solution, active_set, g, spks] = update_g(y, g, spks,smin);
            if abs(g-g0)/g0 < 1e-3 % g is converged
                g_converged = true;
            end
        end
%         g_converged = true;
        res = y - solution;
        RSS = res' * res;
        if or(RSS>thresh, sum(solution)<1e-9)  % constrained form has been found, stop
            break;
        else
            % update lam
            [smin, solution, spks, active_set] = update_smin(y, g, smin,...
                solution, spks, active_set, sqrt(thresh), max(spks));

            % no more change of the active set
            if size(active_set,1)==len_active_set
                break;
            end
        end
    end
else
    %% initialization
    b = quantile(y, 0.15);
    [solution, spks, active_set] = oasisAR2(y-b, g, [], smin);

    %% optimize the baseline b and dependends on the optimized g too
    g_converged = false;
    for miter=1:maxIter
           % update b and g

        res = y - solution - b;
        RSS = res' * res;
        len_active_set = size(active_set,1);

        if or(abs(RSS-thresh) < tol, sum(solution)<1e-9)
            break;
        else
            % update smin
            [smin, solution, spks, active_set] = update_smin(y-b, g, smin,...
                solution, spks, active_set, sqrt(thresh), max(spks));
            b = mean(y-solution);

            % update b and g
            if and(optimize_g, ~g_converged);
                g0 = g;
                [solution, active_set, g, spks] = update_g(y-b, g0, spks, smin);
                if abs(g-g0)/g0 < 1e-4;
                    g_converged = true;
                end
            end

        end
    end

end
c = solution;
s = spks;

%% nested functions
    function [smin, solution, spks, active_set] = update_smin(y, g, smin, solution, ...
            spks, active_set, thresh, s_max)
        %%estimate smin to match the thresholded RSS
        len_active_set = size(active_set, 1);
        sv = linspace(smin, s_max, min(9, len_active_set));
        ind_start = 1;
        ind_end = length(sv);

        while (ind_end-ind_start)>1
            ind = floor((ind_start+ind_end)/2);
            tmp_smin = sv(ind);
            [tmp_solution, tmp_spks, tmp_active_set] = oasisAR2(y, g, [], ...
                tmp_smin, [], [], active_set);
            sqRSS = norm(y-tmp_solution,2);
            if sqRSS<thresh % increase smin
                solution = tmp_solution;
                spks = tmp_spks;
                active_set = tmp_active_set;
                smin = tmp_smin;
                ind_start = ind;
            elseif sqRSS>thresh % decrease smin
                ind_end = ind;
            else
                break;
            end
        end
    end

end

%update the AR coefficient: g
function [c, active_set, g, s] = update_g(y, g, spks, smin)
%% inputs:
%   y:  T*1 vector, One dimensional array containing the fluorescence intensities
%withone entry per time-bin.
%   g: 2x1 vector, AR2 parameters
%   active_set: npools*4 matrix, previous active sets.
% smin: scalr, minimize size of nonzero spikes

%% outputs
%   c: T*1 vector
%   s: T*1 vector, spike train
%   active_set: npool x 4 matrix, active sets
%   g: scalar

%% Authors: Pengcheng Zhou, Carnegie Mellon University, 2016

%% initialization
s_th = quantile(spks(spks>1e-3), 0.25);
tsp = find(spks>=s_th);
tsp = reshape(tsp, 1, []);
time_p = find(conv2(double(spks<=s_th), ones(30,1), 'same')>0);
time_p = reshape(time_p,[],1);
y = reshape(y,[],1);    % fluorescence data
yp = y(time_p);
T = length(y);
tau_dr = ar2exp(g);
tau_d = tau_dr(1);
tau_r = tau_dr(2);

%% find the optimal g and get the warm started active_set
tau_d0 = tau_d;
tau_r0 = tau_r;
bnd_d = tau_d0 * [1/4, 4];
bnd_r = tau_r0 * [1/4, 4];
for m=1:10
    tau_r = fminbnd(@rss_taur, bnd_r(1), bnd_r(2));
    tau_d = fminbnd(@rss_taud, bnd_d(1), bnd_d(2));
    if and(abs(tau_d-tau_d0)/tau_d0 < 1e-4, abs(tau_r-tau_r0)/tau_r0 < 1e-4)
        break;
    else
        tau_d0 = tau_d;
        tau_r0 = tau_r;
    end
end

%% copute the optimal solution
g = exp2ar([tau_d, tau_r]);
[c,s,active_set] = oasisAR2(y, g, 0, smin);

%% nested functions

    function rss = rss_taur(tau_r)
        ht = (exp(-(1:T)/tau_d) - exp(-(1:T)/tau_r))/(tau_d-tau_r);
        ht(T) = 0;
        ind = bsxfun(@minus, time_p, tsp);
        ind(ind<=0) = T;
        V = ht(ind);

        % find the best value of st
        s = (V'*V)\(V'*yp);
        res = yp - V*s;
        rss = res' * res;
    end

    function rss = rss_taud(tau_d)
         ht = (exp(-(1:T)/tau_d) - exp(-(1:T)/tau_r))/(tau_d-tau_r);
        ht(T) = 0;
        ind = bsxfun(@minus, time_p, tsp);
        ind(ind<=0) = T;
        V = ht(ind);

        % find the best value of st
        s = (V'*V)\(V'*yp);
        res = yp - V*s;
        rss = res' * res;
    end
end

author: flatironinstitute

repo: CaImAn-MATLAB


text:
function g = p27_g ( n, x )

%*****************************************************************************80
%
%% P27_G evaluates the gradient for problem 27.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    18 January 2001
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, integer N, the number of variables.
%
%    Input, real X(N), the values of the variables.
%
%    Output, real G(N), the gradient of the objective function.
%
  g = zeros ( n, 1 );

  r = sqrt ( x(1)^2 + x(2)^2 );

  if ( r == 0.0 )
    g(1) = 0.0;
    g(2) = 0.0;
    return
  end

  rx1 = x(1) / r;
  rx2 = x(2) / r;

  a = ( 1.0 + 0.001 * r^2 )^( -2 );
  ar = - 0.004 * r * ( 1.0 + 0.001 * r^2 )^( -3 );

  b = ( sin ( r ) )^2 - 0.5;
  br = sin ( 2.0 * r );

  g(1) = ( ar * b + a * br ) * rx1;
  g(2) = ( ar * b + a * br ) * rx2;

  return
end


text:
function mixgauss = mixgauss_classifier_train(trainFeatures, trainLabels, nc, varargin)
% function mixgauss = mixgauss_classifier_train(trainFeatures, trainLabels, nclusters, varargin)
% trainFeatures(:,i) for i'th example
% trainLabels should be 0,1
% To evaluate performance on a tets set, use
% mixgauss = mixgauss_classifier_train(trainFeatures, trainLabels, nc, 'testFeatures', tf, 'testLabels', tl)

[testFeatures, testLabels, max_iter, thresh, cov_type, mu, Sigma, priorC, method, ...
 cov_prior, verbose, prune_thresh] = process_options(...
    varargin, 'testFeatures', [], 'testLabels', [], ...
     'max_iter', 10, 'thresh', 0.01, 'cov_type', 'diag', ...
    'mu', [], 'Sigma', [], 'priorC', [], 'method', 'kmeans', ...
    'cov_prior', [], 'verbose', 0, 'prune_thresh', 0);

Nclasses = 2; % max([trainLabels testLabels]) + 1;

pos = find(trainLabels == 1);
neg = find(trainLabels == 0);

if verbose, fprintf('fitting pos\n'); end
[mixgauss.pos.mu, mixgauss.pos.Sigma, mixgauss.pos.prior] = ...
    mixgauss_em(trainFeatures(:, pos), nc, varargin{:});

if verbose, fprintf('fitting neg\n'); end
[mixgauss.neg.mu, mixgauss.neg.Sigma, mixgauss.neg.prior] = ...
    mixgauss_em(trainFeatures(:, neg), nc, varargin{:});

if ~isempty(priorC)
  mixgauss.priorC = priorC;
else
  mixgauss.priorC = normalize([length(pos) length(neg)]);
end

author: fieldtrip

repo: fieldtrip

save_path: github-repos/MATLAB/fieldtrip-fieldtrip


text:
function [fit, pairs] = calcGBFit(job,varargin)
% fit between child to child misorientations and job.p2c
%
% Syntax
%
%   [fit, pairs] = job.calcGBFit
%
%   % visualize the result
%   [gB, pairId] = job.grains.boundary.selectByGrainId(pairs);
%   plot(gB, fit(pairId) ./ degree, 'linewidth',2);
%   setColorRange([2,8])
%   mtexColorMap white2black
%
%   % consider only P2C neighbors
%   [fit, p2cPairs] = job.calcGBFit('p2c')
%
%   % consider only C2C neighbors
%   [fit, c2cPairs] = job.calcGBFit('c2c')
%s
% Input
%  job - @parentGrainReconstructor
%
% Output
%  fit   - fit between the grain boundary misorientations and job.p2c
%  pairs - list of grainId of the neighboring grains
%

noOpt = ~check_option(varargin,{'p2c','c2c'});

% parent to child neighbors
if noOpt || check_option(varargin,'p2c')

  % extract all parent to child neighbors
  [pairs, oriParent, oriChild] = getP2CPairs(job, varargin{:});

  % compute the misorientation
  mori = inv(oriChild).*oriParent;

  % compute the fit
  fit = angle(mori, job.p2c);

else

  fit = [];
  pairs = [];

end

if noOpt || check_option(varargin,'c2c')

  % extract all child to child neighbors
  [c2cPairs, oriChild] = getC2CPairs(job, varargin{:});
  pairs = [pairs;c2cPairs];

  % compute the corresponding misorientations
  mori = inv(oriChild(:,1)).*oriChild(:,2);

  % misorientation to c2c variants
  fit = [fit;min(angle_outer(mori, job.c2c),[],2)];

end

author: mtex-toolbox

repo: mtex

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


text:
function X=RiccatiSolveC(A,B,Q,R,S,E)
%%RICCATISOLVEC Solve the general continuous-time Riccati equation.
%
%INPUTS: A An nXn matrix.
%        B An nXm matrix, where m<=n.
%        Q An nXn matrix such that Q=Q' and all eigenvalues are non-
%          negative.
%        R An optional mXm matrix such that R=R' and all eigenvalues are
%          non-negative. If omitted, eye(m,m) is used.
%        S An optional nXm matrix. If omitted, zeros(n,n) is used.
%        E An optional nXn matrix. If omitted, eye(n) is used.
%
%OUTPUTS: X The nXn nonnegative definite solution to the discrete-time
%           algebraic Ricatti equation.
%
%This function finds the nonnegative definite solution to the
%continuous-time Ricatti equation having the form
%A'*X*E+E'*X*A-(E'*X*B+S)*inv(R)*(B'*X*E+S')+Q=0
%or, if R, S, and E are omitted, the equation under consideration becomes
%A'*X+X*A-X*B*B'*X+Q=0
%The algorithm of  [1] is used. Note that this function does not work with
%problems of the form -X*B*B'*X+Q=0 due to numerical issues.
%
%The continuous-time Ricatti equation arises when solving for the  steady-
%state covariance of a continuous-time linear Kalman filter, as described
%in Chapter 9.2.3 of [2].
%
%REFERENCES:
%[1] W. F. Arnold III and A. J. Laub, "Generalized eigenproblem algorithms
%    and software for algebraic Riccati equations," Proceedings of the
%    IEEE, vol. 72, no. 12, pp. 1746-1754, Dec. 1984.
%[2] Y. Bar-Shalom, X. R. Li, and T. Kirubarajan, Estimation with
%    Applications to Tracking and Navigation. New York: John Wiley and
%    Sons, Inc, 2001.
%
%October 2013 David F. Crouse, Naval Research Laboratory, Washington D.C.
%(UNCLASSIFIED) DISTRIBUTION STATEMENT A. Approved for public release.

n=size(A,1);
m=size(B,2);

if(nargin<4)
    R=eye(m,m);
end

if(nargin<5)
    S=zeros(n,m);
end

if(nargin<6)
   E=eye(n);
end

L=[E,            zeros(n,n),zeros(n,m);
   zeros(n,n),   E',        zeros(n,m);
   zeros(m,n),   zeros(m,n),zeros(m,m)];
M=[A,   zeros(n,n),  B;
  -Q,  -A',         -S;
   S',  B',          R];

[MHat,LHat,V,U]=qz(M,L,'real');
[~,~,~,U]=ordqz(MHat,LHat,V,U,'lhp');

W=[E,           zeros(n,n),zeros(n,m);
   zeros(n,n),  eye(n),zeros(n,m)]*U;

W11=W(1:n,1:n);
W21=W((n+1):(2*n),1:n);

X=W21/W11;
end

%LICENSE:
%
%The source code is in the public domain and not licensed or under
%copyright. The information and software may be used freely by the public.
%As required by 17 U.S.C. 403, third parties producing copyrighted works
%consisting predominantly of the material produced by U.S. government
%agencies must provide notice with such work(s) identifying the U.S.
%Government material incorporated and stating that such material is not
%subject to copyright protection.
%
%Derived works shall not identify themselves in a manner that implies an
%endorsement by or an affiliation with the Naval Research Laboratory.
%
%RECIPIENT BEARS ALL RISK RELATING TO QUALITY AND PERFORMANCE OF THE
%SOFTWARE AND ANY RELATED MATERIALS, AND AGREES TO INDEMNIFY THE NAVAL
%RESEARCH LABORATORY FOR ALL THIRD-PARTY CLAIMS RESULTING FROM THE ACTIONS
%OF RECIPIENT IN THE USE OF THE SOFTWARE.

author: USNavalResearchLaboratory

repo: TrackerComponentLibrary


text:
function group = BF_ToGroup(groupIndices,maxLength)
% BF_ToGroup    convert from cell form to vector form of group labels
%
% Converts to a vector of group labels from a cell, where each element is a
% vector of indices for that group (or vice-versa).

% ------------------------------------------------------------------------------
% Copyright (C) 2020, Ben D. Fulcher <ben.d.fulcher@gmail.com>,
% <http://www.benfulcher.com>
%
% If you use this code for your research, please cite the following two papers:
%
% (1) B.D. Fulcher and N.S. Jones, "hctsa: A Computational Framework for Automated
% Time-Series Phenotyping Using Massive Feature Extraction, Cell Systems 5: 527 (2017).
% DOI: 10.1016/j.cels.2017.10.001
%
% (2) B.D. Fulcher, M.A. Little, N.S. Jones, "Highly comparative time-series
% analysis: the empirical structure of time series and their methods",
% J. Roy. Soc. Interface 10(83) 20130048 (2013).
% DOI: 10.1098/rsif.2013.0048
%
% This function is free software: you can redistribute it and/or modify it under
% the terms of the GNU General Public License as published by the Free Software
% Foundation, either version 3 of the License, or (at your option) any later
% version.
%
% This program is distributed in the hope that it will be useful, but WITHOUT
% ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
% FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
% details.
%
% You should have received a copy of the GNU General Public License along with
% this program. If not, see <http://www.gnu.org/licenses/>.
% ------------------------------------------------------------------------------

if iscell(groupIndices) % Convert to vector of group indices

    if all(cellfun(@isempty,groupIndices)) % all empty -- rubbish
        error('Nothing in the group indices')
    end

    % Second input only important for cell -> vector transformation
    if nargin < 2 || isempty(maxLength)
        maxLength = max(cellfun(@max,groupIndices));
        % Make the length of output equal to the maximum index
    end

    group = zeros(maxLength,1);

    for i = 1:length(groupIndices)
        group(groupIndices{i}) = i;
    end

    if sum(cellfun(@length,groupIndices)) < length(group)
        warning('Group is missing some labels')
    end

else
    % Convert from vector of categorical group labels to cell of group indices:
    if iscategorical(groupIndices)
        classLabels = categories(groupIndices);
        numGroups = length(classLabels);
        group = cell(numGroups,1);
        for i = 1:numGroups
            group{i} = find(groupIndices==classLabels{i});
        end
    else
        numGroups = max(groupIndices);
        group = cell(numGroups,1);
        for i = 1:numGroups
            group{i} = find(groupIndices==i);
        end
    end
end

end

author: benfulcher

repo: hctsa

save_path: github-repos/MATLAB/benfulcher-hctsa


text:
function p = predictOneVsAll(all_theta, X)
%PREDICT Predict the label for a trained one-vs-all classifier. The labels
%are in the range 1..K, where K = size(all_theta, 1).
%  p = PREDICTONEVSALL(all_theta, X) will return a vector of predictions
%  for each example in the matrix X. Note that X contains the examples in
%  rows. all_theta is a matrix where the i-th row is a trained logistic
%  regression theta vector for the i-th class. You should set p to a vector
%  of values from 1..K (e.g., p = [1; 3; 1; 2] predicts classes 1, 3, 1, 2
%  for 4 examples)

m = size(X, 1);
num_labels = size(all_theta, 1);

% You need to return the following variables correctly
p = zeros(size(X, 1), 1);

% Add ones to the X data matrix
X = [ones(m, 1) X];

% ====================== YOUR CODE HERE ======================
% Instructions: Complete the following code to make predictions using
%               your learned logistic regression parameters (one-vs-all).
%               You should set p to a vector of predictions (from 1 to
%               num_labels).
%
% Hint: This code can be done all vectorized using the max function.
%       In particular, the max function can also return the index of the
%       max element, for more information see 'help max'. If your examples
%       are in rows, then, you can use max(A, [], 2) to obtain the max
%       for each row.
%

% x * all_theta' : 5000 * 10，每一行的每一个列值代表图像是1,2...10的可能性
% use max(),M = max(A,[],dim) 沿着维度 dim 返回最大元素。例如，如果 A 为矩阵，则 max(A,[],2) 是包含每一行的最大值的列向量。
% probaitilty代表每一行最大的值，p代表索引
[probability, p] = max(sigmoid(X * all_theta'), [], 2);

% =========================================================================

end

author: xjwhhh

repo: AndrewNgMachineLearning


text:
function [y,n] = ADEM_plaid(x,n)
% creates a Gaussian modulated n x n visual plaid stimulus
% FORMAT [y,n] = ADEM_plaid(x,[n])
% x(1) - horizontal displacement
% x(2) - vertical displacement
%__________________________________________________________________________
% Copyright (C) 2008 Wellcome Trust Centre for Neuroimaging

% Karl Friston
% $Id: ADEM_plaid.m 3655 2009-12-23 20:15:34Z karl $

% default stimulus
%--------------------------------------------------------------------------
try
    n;
catch
    n = 4;
end

% stimulus
%--------------------------------------------------------------------------
sx  = [1:n]';
sx  = sx - n/2 - n*x(1)/8;
sy  = [1:n] ;
sy  = sy - n/2 - n*x(2)/8;
sx  = exp(-sx.^2/(2*(n/6)^2)).*cos(2*pi*sx*2/n);
sy  = exp(-sy.^2/(2*(n/6)^2)).*cos(2*pi*sy*2/n);

% vectorise under defaults
%--------------------------------------------------------------------------
if nargin > 1
    y = sx*sy;
else
    y = spm_vec(sx*sy);
end

author: spm


text:
function [Z, active_Z] = generate_reference(PopObj, original_active_number, Global, action, mode)

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

global LEARNING_DISABLE_FLAG;
if strcmp(mode, 'strict')
    active_number = 0;
    PopObj = normalization(PopObj);
    active_Z = [];
    while active_number < original_active_number
        instructed_action = action;
        [Z, received_action] = generate_Z(Global, instructed_action, 'commit');
        Z = unique([Z; active_Z], 'rows');
        [~, Allocation] = pair(PopObj, Z, 'sin');
        active_index = unique(Allocation);
        active_Z = Z(active_index, :);
        active_number = size(active_Z, 1);
        fprintf('%d(+%d)\n', active_number, numel(active_index)); % TODO
        if received_action ~= instructed_action
            LEARNING_DISABLE_FLAG = true;
            warning('reference generation malfunctions');
            break;
        end
    end
elseif strcmp(mode, 'normal')
    Z = generate_Z(Global, action, 'commit');
    active_Z = [];
end
end

function varargout = generate_Z(Global, action, commit_str)
global MAX_REFERENCE_NUM_FLAG current_density reference_generation_intial_flag;
persistent incremental_sequence pointer;
if reference_generation_intial_flag == true
    MAX_REFERENCE_NUM_FLAG = false;
    incremental_sequence = density_sequence(Global.M);
    pointer = 1;
    reference_generation_intial_flag = false;
    % TRIM
    TRIM_INDEX = find(incremental_sequence < Global.N);
    incremental_sequence(TRIM_INDEX(1: end)) = [];
end
original_pointer = pointer;
pointer = pointer + action;
if pointer >= numel(incremental_sequence)
    MAX_REFERENCE_NUM_FLAG = true;
end
pointer = min(numel(incremental_sequence), pointer);
pointer = max(1, pointer);
real_action = pointer - original_pointer;

Z = UniformPoint(incremental_sequence(pointer), Global.M);
current_density = size(Z, 1);

if ~strcmp(commit_str, 'commit')
    pointer = original_pointer;
end

varargout{1} = Z;
if nargout >= 2
    varargout{2} = real_action;
end
end

function incremental_sequence = density_sequence(M)
incremental_sequence = [];
R_MAX = 1e6;
step = 1;
R = nchoosek(M + step - 1, step);
while R < R_MAX
    incremental_sequence = [incremental_sequence, R];
    step = step + 1;
    R = nchoosek(M + step - 1, step);
end
end

author: BIMK

repo: PlatEMO

save_path: github-repos/MATLAB/BIMK-PlatEMO


text:
function a = double(s)
%DOUBLE Converts a sparse tensor to a dense multidimensional array.
%
%  See also SPTENSOR, SPTENSOR/FULL.
%
%MATLAB Tensor Toolbox.
%Copyright 2015, Sandia Corporation.

% This is the MATLAB Tensor Toolbox by T. Kolda, B. Bader, and others.
% http://www.sandia.gov/~tgkolda/TensorToolbox.
% Copyright (2015) Sandia Corporation. Under the terms of Contract
% DE-AC04-94AL85000, there is a non-exclusive license for use of this
% work by or on behalf of the U.S. Government. Export of this data may
% require a license from the United States Government.
% The full license terms can be found in the file LICENSE.txt

a = zeros([size(s) 1 1]);
if nnz(s) > 0
    a(tt_sub2ind(size(s),s.subs)) = s.vals;
end

author: zhangqianqianQQ

repo: MachineVisionAlgorithm


text:
function sinh_values_test ( )

%*****************************************************************************80
%
%% SINH_VALUES_TEST demonstrates the use of SINH_VALUES.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    23 June 2007
%
%  Author:
%
%    John Burkardt
%
  fprintf ( 1, '\n' );
  fprintf ( 1, 'SINH_VALUES_TEST:\n' );
  fprintf ( 1, '  SINH_VALUES stores values of \n' );
  fprintf ( 1, '  the hyperbolic sine function.\n' );
  fprintf ( 1, '\n' );
  fprintf ( 1, '      X           FX\n' );
  fprintf ( 1, '\n' );

  n_data = 0;

  while ( 1 )

    [ n_data, x, fx ] = sinh_values ( n_data );

    if ( n_data == 0 )
      break
    end

    fprintf ( 1, '  %12f  %24.16f\n', x, fx );

  end

  return
end

author: johannesgerer

repo: jburkardt-m


text:
% startup.m
%
%  Discussion:
%
%    MATLAB always has a "current directory".
%
%    This current directory starts out with a default directory often located
%    somewhere in the file system associated with the MATLAB application.
%
%    It's convenient to have this initial directory changed to somewhere
%    associated with YOUR files, so that you can easily draw MATLAB's attention
%    to various files you want to run.
%
%    By placing this file, "startup.m", somewhere in MATLAB's search path,
%    you can move MATLAB's current directory to a more suitable location.
%
%    For that matter, you can also add directories to the search path,
%    define constants, and take care of other matters that produce a pleasing
%    environment for you.
%
%    For instance, under the Macintosh OSX operating system, the main MATLAB
%    directory might be /Applications//MATLAB_R2007b or
%    /Applications/MATLAB_R2009b.app.  It's enough to put the startup.m file
%    directly in this folder.  If you're neat, you might prefer to place it in
%    the toolbox/local subdirectory.  MATLAB will still be able to find and
%    execute it.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    24 August 2011
%
%  Author:
%
%    John Burkardt
%
  beep off

  fprintf ( 1, '\n' );
  fprintf ( 1, 'MATLAB''s default startup directory is: "%s".\n', pwd );
  cd /Users/jburkardt/public_html/m_src
  fprintf ( 1, 'Moving to preferred directory: "%s".\n', pwd );


text:
function cicros(opt)
    %   This subroutine "circle"  selects the Ni closest earthquakes
    %   around a interactively selected point.  Resets ZG.newcat and ZG.newt2
    %   Operates on "primeCatalog".
    %
    % axis: h1
    % plots to: plos1 as xk
    % inCatalog: a
    % outCatalog: newt2, newcat, newa, newa2
    % mouse controlled
    % closest events OR radius
    % calls: bdiff
    %
    %  Input Ni:
    %
    persistent ic

    report_this_filefun();
    ZG=ZmapGlobal.Data;

    if opt==0 && ~isempty(ic)
        opt=ic;
    end

    axes(h1)

    titStr ='Selecting EQ in Circles                         ';
    messtext= ...
        ['                                                '
        '  Please use the LEFT mouse button              '
        ' to select the center point.                    '
        ' The "ni" events nearest to this point          '
        ' will be selected and displayed in the map.     '];

    msg.dbdisp(messtext, titStr);

    % Input center of circle with mouse
    %
    [xa0,ya0]  = ginput(1);

    stri1 = [ 'Circle: lon = ' num2str(xa0) '; lat= ' num2str(ya0)];
    stri = stri1;
    pause(0.1)
    %  calculate distance for each earthquake from center point
    %  and sort by distance
    %
    l = sqrt(((xsecx' - xa0)).^2 + ((xsecy + ya0)).^2) ;
    [s,is] = sort(l);
    ZG.newt2 = newa(is(:,1),:) ;

    switch(ic)
        case 1 % select  N clostest events

            l =  sort(l);
            messtext = ['Radius of selected Circle:' num2str(l(ni))  ' km' ];
            disp(messtext)
            msg.dbdisp(messtext, 'Message')
            %
            % take first ni and sort by time
            %
            ZG.newt2 = ZG.newt2(1:ni,:);
            ZG.newt2.sort('Date');
            %
            % plot Ni clostest events on map as 'x':

            set(gca,'NextPlot','add')
            [na,ma] = size(ZG.newt2);
            plot(ZG.newt2(:,ma),-ZG.newt2.Depth,'xk','Tag','plos1');
            set(gcf,'Pointer','arrow')
            %
            % plot circle containing events as circle
            x = -pi-0.1:0.1:pi;
            plot(xa0+sin(x)*l(ni), ya0+cos(x)*l(ni),'w')
            l(ni)

            %
            ZG.newcat = ZG.newt2;                   % resets ZG.newcat and ZG.newt2

            clear l s is
            bdiff(ZG.newt2)

        case 2 % select  events within ra

            l =  sort(l);
            ll = l <=ra;
            messtext = ['Number of events in Circle :' num2str(sum(ll)) ];
            disp(messtext)
            msg.dbdisp(messtext, 'Message')
            %
            % take first ni and sort by time
            %
            ZG.newt2 = ZG.newt2.subset(ll);
            ZG.newt2.sort('Date');
            %
            % plot Ni clostest events on map as 'x':

            set(gca,'NextPlot','add')
            [na,ma] = size(ZG.newt2);
            plot(ZG.newt2(:,ma),-ZG.newt2.Depth,'xk','Tag','plos1');
            set(gcf,'Pointer','arrow')
            %
            % plot circle containing events as circle
            x = -pi-0.1:0.1:pi;
            plot(xa0+sin(x)*ra, ya0+cos(x)*ra,'w')
            l(ni)

            %
            ZG.newcat = ZG.newt2;                   % resets ZG.newcat and ZG.newt2

            clear l s is
            bdiff(ZG.newt2)

        case 3 % select  events within ra

            ax=findobj(gcf,'Tag','mainmap_ax');
            [x,y, mouse_points_overlay] = select_polygon(ax);
            plot(x,y,'b-');
            YI = -newa(:,7);          % this substitution just to make equation below simple
            XI = newa(:,end);
            ll = polygon_filter(x,y, XI, YI, 'inside');

            %plot the selected eqs and mag freq curve
            newa2 = newa.subset(ll);
            ZG.newt2 = newa2;
            ZG.newcat = newa.subset(ll);
            pl = plot(newa2(:,end),-newa2(:,7),'xk');
            set(pl,'MarkerSize',5,'LineWidth',1)
            bdiff(newa2)
        otherwise
            error('no option specified. 1: N closest events within radius, 2: radius')
    end
    ic=opt;
end


text:
function out = isempty(F)
%ISEMPTY   Empty Boolean check for a CHEBFUN3V object.
%   ISEMPTY(F) returns 1 if every component of F is an empty CHEBFUN3V, and
%   return 0 otherwise.

% Copyright 2017 by The University of Oxford and The Chebfun Developers.
% See http://www.chebfun.org/ for Chebfun information.

if ( isempty(F.components) )
    out = 1;
    return
end

% Take isempty of each component:
out = cellfun(@isempty, F.components, 'UniformOutput', false);
out = all(cell2mat(out));

end

author: chebfun

repo: chebfun

save_path: github-repos/MATLAB/chebfun-chebfun


text:
%SETFEATSIZE Set the feature size of a dataset
%
%   A = SETFEATSIZE(A,FEATSIZE)
%
% INPUT
%   A         Dataset
%   FEATSIZE  Feature size vector
%
% OUTPUT
%   A         Dataset
%
% DESCRIPTION
% By default the feature size of a dataset is its number of features, i.e.
% the number of columns in the DATA field of A. If the features are samples
% of a multi-dimensional data item, e.g. the pixels of an image, the
% original size of this data item may be stored in FEATSIZE. The product of
% all elements in FEATSIZE has to be equal to the number of columns in the
% DATA field of A.

author: marianux

repo: ecg-kit

save_path: github-repos/MATLAB/marianux-ecg-kit


text:
function [transformed_sequence] = sequence_pixelchange(sequence, operation, name)
% sequence_pixelchange Returns sequence with arbitrary pixel transformation
%
% This sequence converter returns sequence with arbitrary pixel transformation that preserves the
% size of the image and the groundtruth.
%
% Cache notice: The results of this function are cached in the workspace cache directory.
%
% Input:
% - sequence (structure): A valid sequence structure.
% - operation (function): A handle of the pixel transformation function.
% - name (string, optional): Name of the operation for caching purposes.
%
% Output:
% - tranform_sequence (structure): A sequence descriptor of a converted sequence.

if ischar(operation)
    operation_name = operation;
    operation = str2func(operation_name);
else
    operation_name = func2str(operation);
end;

if ~isempty(name) && ischar(name)
    operation_name = name;
end;

cache_directory = fullfile(get_global_variable('directory'), 'cache', ...
    sprintf('pixelchange_%s', operation_name), sequence.name);

mkpath(cache_directory);

cache_groundtruth = fullfile(cache_directory, 'groundtruth.txt');

sequence_groundtruth = fullfile(sequence.directory, 'groundtruth.txt');

if file_newer_than(cache_groundtruth, sequence_groundtruth)
    transformed_sequence = sequence_create(cache_directory, 'name', sequence.name);
    transformed_sequence.values.names = sequence.values.names;
    transformed_sequence.values.data = sequence.values.data;
    return;
end;

print_debug('Generating cached sequence ''%s'' for operation ''%s''...', sequence.name, operation_name);

tags_cache_names = cell(0, 0);
tags_cache_data = false(sequence.length, 0);

for i = 1:sequence.length

    original_image = imread(sequence_get_image(sequence, i));
    original_tags = get_tags(sequence, i);

    [transformed_image, transformed_tags] = operation(original_image, ...
        original_tags, i, sequence.length);

    for l = 1:length(transformed_tags)
        if ~ismember(transformed_tags{l}, tags_cache_names)
            tags_cache_data = cat(2, tags_cache_data, false(sequence.length, 0));
            tags_cache_names{end+1} = transformed_tags{l};  %#ok<AGROW>
        end;
        tags_cache_data(i, strcmp(tags_cache_names, transformed_tags{l})) = 1;
    end;

    if size(transformed_image, 3) == 1
        transformed_image = repmat(rgb2gray(transformed_image), [1 1 3]);
    end;

    if size(original_image, 1) ~= size(transformed_image, 1) || ...
            size(original_image, 2) ~= size(transformed_image, 2)
        error('The operation should return image of same width and height.');
    end;

    imwrite(transformed_image, fullfile(cache_directory, sprintf('%08d.jpg', i)));

end;

for l = 1:length(tags_cache_names)
    csvwrite(fullfile(cache_directory, sprintf('%s.tag', ...
        tags_cache_names{l})), tags_cache_data(:, l));
end;

write_trajectory(cache_groundtruth, sequence.groundtruth);

transformed_sequence = sequence_create(cache_directory, 'name', sequence.name);

transformed_sequence.values.names = sequence.values.names;
transformed_sequence.values.data = sequence.values.data;

author: votchallenge

repo: toolkit-legacy


text:
function box = mergeBoxes(box1, box2)
%MERGEBOXES Merge two boxes, by computing their greatest extent.
%
%   BOX = mergeBoxes(BOX1, BOX2);
%
%   Example
%   box1 = [5 20 5 30];
%   box2 = [0 15 0 15];
%   mergeBoxes(box1, box2)
%   ans =
%       0 20 0 30
%
%
%   See also
%   boxes2d, drawBox, intersectBoxes
%

% ------
% Author: David Legland
% E-mail: david.legland@grignon.inra.fr
% Created: 2010-07-26, using Matlab 7.9.0.529 (R2009b)
% Copyright 2010-2022 INRA - Cepia Software Platform

% unify sizes of data
if size(box1,1) == 1
    box1 = repmat(box1, size(box2,1), 1);
elseif size(box2, 1) == 1
    box2 = repmat(box2, size(box1,1), 1);
elseif size(box1,1) ~= size(box2,1)
    error('Bad size for inputs');
end

% compute extreme coords
mini = min(box1(:,[1 3]), box2(:,[1 3]));
maxi = max(box1(:,[2 4]), box2(:,[2 4]));

% concatenate result into a new box structure
box = [mini(:,1) maxi(:,1) mini(:,2) maxi(:,2)];

author: mattools

repo: matGeom

save_path: github-repos/MATLAB/mattools-matGeom


text:
function [ as, bs, cs ] = sphere01_triangle_vertices_to_sides ( v1, v2, v3 )

%*****************************************************************************80
%
%% SPHERE01_TRIANGLE_VERTICES_TO_SIDES computes spherical triangle sides on unit sphere.
%
%  Discussion:
%
%    We can use the ACOS system call here, but the ARC_COSINE routine
%    will automatically take care of cases where the input argument is
%    (usually slightly) out of bounds.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    22 September 2010
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, real V1(3), V2(3), V3(3), the vertices of the spherical
%    triangle.
%
%    Output, real AS, BS, CS, the (geodesic) length of the sides
%    of the triangle.
%
  as = r8_acos ( v2(1:3)' * v3(1:3) );
  bs = r8_acos ( v3(1:3)' * v1(1:3) );
  cs = r8_acos ( v1(1:3)' * v2(1:3) );

  return
end


text:
classdef asRoiClass < impoly
    properties (Access = private)
        parentAxesHandle   = 0;     % parent axes handle
        parentImageHandle  = 0;     % parent image handle
        parentFigureHandle = 0;

        textHandle            = [];
        textContextMenuHandle = [];
        precision = '%g';      % standard notation is compact

        filterStr = ''      % WARNING: test

        guiHandle    = struct('figure','text')
        cmenuHandle  = struct('sendPosition','ignoreZeros');

        sendPositionCallback = [];
        sendPositionCallbackId = [];

        fullStrToggle = true;  % show full string or just compact form
    end

    properties (Access = public)

    end

    methods
        function obj = asRoiClass(parentAxesHandle, roiPos, sendPositionCallback)
            if nargin < 2
                roiPos = [];
                if nargin < 1
                    parentAxesHandle = gca;
                end
            end

            obj     = obj@impoly(parentAxesHandle, roiPos);
            obj.setColor('green')
            obj.parentAxesHandle   = parentAxesHandle;
            obj.parentImageHandle  = findall(get(parentAxesHandle,'Children'),'Type','image');
            obj.parentFigureHandle = get(get(obj.parentAxesHandle,'Parent'),'Parent');

            if nargin == 3
                obj.sendPositionCallback = sendPositionCallback;
            end

            obj.updateImpolyContextMenu;
            obj.createTextContextMenu;

%             addNewPositionCallback(obj,@(pos)obj.showRoiGui);
%             obj.showRoiGui;
            addNewPositionCallback(obj,@(pos)obj.updateRoiString);
            obj.updateRoiString;

        end

        function roi = getRoiData(obj)
            ud = get(obj.parentAxesHandle,'UserData');
            if isempty(ud)
                refImg = get(obj.parentImageHandle,'CData');
            else
                if isfield(ud,'isComplex') && ud.isComplex
                    refImg = ud.cplxImg;
                else
                    refImg = get(obj.parentImageHandle,'CData');
                end
            end
            mask   = obj.createMask;
            roi    = refImg(mask == 1);
            if obj.getIgnoreZerosToggle
                roi = roi(roi~=0);
            end
            if ~isempty(obj.filterStr)
                eval(['roi = roi(roi',obj.filterStr,');']);
            end

        end

        function [m, s] = getMeanAndStd(obj)
           roi = obj.getRoiData;
           if ~isempty(roi)
               m = mean(roi);
               s = std(roi);
           else
               m = 0;
               s = 0;
           end
        end

        function N = getN(obj)
           roi = obj.getRoiData;
           N = numel(roi);
        end
        function mini = getMin(obj)
            mini = min(obj.getRoiData);
        end
        function maxi = getMax(obj)
            maxi = max(obj.getRoiData);
        end

        function updateRoiString(obj)
            if obj.fullStrToggle
                obj.drawFullString;
            else
                obj.drawMeanAndStdString;
            end
        end

        function str = getMeanAndStdString(obj)
            [m, s] = obj.getMeanAndStd;
            str = [num2str(m,obj.precision), ' +- ', num2str(s,obj.precision)];
        end

        function drawMeanAndStdString(obj)
            % writes roi statistics as a text in the image window
            if ~isempty(obj.textHandle)
                delete(obj.textHandle);
            end
            obj.textHandle = text(.015,.99,obj.getMeanAndStdString,'Units','normalized','parent',obj.parentAxesHandle,...
                'Color','green','BackgroundColor','Black',...
                'VerticalAlignment','top',...
                'UIContextMenu',obj.textContextMenuHandle);
        end

        function drawFullString(obj)
            % writes roi statistics as a text in the image window
            if ~isempty(obj.textHandle)
                delete(obj.textHandle);
            end
            str = strvcat(obj.getMeanAndStdString,num2str(obj.getN),num2str(obj.getMin,obj.precision),num2str(obj.getMax,obj.precision));
            obj.textHandle = text(.015,.98,str,'Units','normalized','parent',obj.parentAxesHandle,...
                'Color','green','BackgroundColor','Black',...
                'VerticalAlignment','top',...
                'UIContextMenu',obj.textContextMenuHandle);
        end

        function h = getTextHandle(obj)
            h = obj.textHandle;
        end

        function bool = guiIsPresent(obj)
            bool = false;
            if ishandle(obj.guiHandle.figure)
                if strcmp(get(obj.guiHandle.figure,'Tag') ,'roiGui')
                    bool = true;
                end
            end
        end

        function showRoiGui(obj)
            if obj.guiIsPresent
                set(obj.guiHandle.text,'String',obj.getMeanAndStdString);
            else
                pos = [800, 800, 40, 40];
                str = obj.getMeanAndStdString;

                obj.guiHandle.figure = figure('MenuBar','none','Toolbar','none',...
                    'Position',pos, 'Tag', 'roiGui');
                obj.guiHandle.text = uicontrol('Style','Text','String',str,'HorizontalAlignment','left',...
                    'Units','normalized','pos',[0 0 1 1],...
                    'parent',obj.guiHandle.figure,'HandleVisibility','on',...
                    'FontUnits','normalized','FontSize',.35);
            end
        end

        function delete(obj)
            if obj.guiIsPresent
                close(obj.guiHandle.figure);
            end
            if ~isempty(obj.textHandle) && ishandle(obj.textHandle)
                delete(obj.textHandle);
            end
        end

        function toggle = getSendPositionToggle(obj)
            switch get(obj.cmenuHandle.sendPosition,'Checked')
                case 'on'
                    toggle = true;
                case 'off'
                    toggle = false;
            end
        end

        function toggle = getIgnoreZerosToggle(obj)
            switch get(obj.cmenuHandle.ignoreZeros,'Checked')
                case 'on'
                    toggle = true;
                case 'off'
                    toggle = false;
            end
        end

        function setSendPositionToggle(obj, toggle)
            if nargin < 2
                toggle = ~obj.getSendPositionToggle;
            end

            switch toggle
                case 1
                    set(obj.cmenuHandle.sendPosition,'Checked','on');
                    obj.sendPositionCallbackId = addNewPositionCallback(obj,obj.sendPositionCallback);
                    obj.callSendPositionCallback(); % execute callback once
                case 0
                    set(obj.cmenuHandle.sendPosition,'Checked','off');
                    removeNewPositionCallback(obj,obj.sendPositionCallbackId);
            end
        end

        function addFilterString(obj,str)
            % WARNING, INCOMPLETE IMPLEMENTATION
            obj.filterStr = str;
            obj.updateRoiString;
        end

        function callSendPositionCallback(obj)
            obj.sendPositionCallback(obj.getPosition);
        end

        function setIgnoreZerosToggle(obj, toggle)
            if nargin < 2
                toggle = ~obj.getIgnoreZerosToggle;
            end

            switch toggle
                case 1
                    set(obj.cmenuHandle.ignoreZeros,'Checked','on');
                case 0
                    set(obj.cmenuHandle.ignoreZeros,'Checked','off');
            end
            obj.updateRoiString;
        end

    end

    methods (Access = private)
        function updateImpolyContextMenu(obj)
            % add some features to the impoly context menu

            cmh = obj.getContextMenu;

            uimenu(cmh,'Label','Delete ROI'   ,...
                'callback',@(src,evnt)obj.delete);

            uimenu(cmh,'Label','Delete all ROIs'   ,...
                'callback',@(src,evnt)evalin('base','asDeleteAllRois'));

            if ~isempty(obj.sendPositionCallback)
                obj.cmenuHandle.sendPosition = uimenu(cmh,'Label','Send Position'   ,...
                    'callback',@(src,evnt)obj.setSendPositionToggle);
            end

            obj.cmenuHandle.ignoreZeros = uimenu(cmh,'Label','Ignore Zeros'   ,...
                'callback',@(src,evnt)obj.setIgnoreZerosToggle,...
                'Checked','on');

        end

        function createTextContextMenu(obj)
            % create a context menu to choose between different notations
            % in the image text

            if isempty(obj.textContextMenuHandle)
                obj.textContextMenuHandle = uicontextmenu('Parent',obj.parentFigureHandle);

                uimenu(obj.textContextMenuHandle,'Label','Decimal notation'   ,...
                    'callback',@(src,evnt)obj.setNotation('%d'));
                uimenu(obj.textContextMenuHandle,'Label','Fixed-point notation'   ,...
                    'callback',@(src,evnt)obj.setNotation('%2.2f'));
                uimenu(obj.textContextMenuHandle,'Label','Exponential notation'   ,...
                    'callback',@(src,evnt)obj.setNotation('%2.2e'));
                uimenu(obj.textContextMenuHandle,'Label','Compact exponential notation'   ,...
                    'callback',@(src,evnt)obj.setNotation('%g'));

            end
        end

        function setNotation(obj, precision)
            obj.precision = precision;
            obj.updateRoiString % update text
        end

    end
end

author: leoliuf

repo: MRiLab

save_path: github-repos/MATLAB/leoliuf-MRiLab