text:
function [success,msg] = ConvertLarvaeRiveraAlba2JAABA(varargin)

success = false;
msg = '';

[inmoviefile,intrxfile,...
  expdir,moviefilestr,trxfilestr,perframedirstr,...
  dosoftlink] = myparse(varargin,...
  'inmoviefile','','intrxfile','',...
  'expdir','','moviefilestr','movie.ufmf','trxfilestr','trx.mat','perframedirstr','perframe',...
  'dosoftlink',false);

% check that required inputs are given
if isempty(inmoviefile),
  msg = 'Input movie file is empty';
  return;
end
if isempty(intrxfile),
  msg = 'Input trx mat file is empty';
  return;
end

% output file locations
moviefile = fullfile(expdir,moviefilestr);
trxfile = fullfile(expdir,trxfilestr);
perframedir = fullfile(expdir,perframedirstr);

% load in the trx
try
  td = load(intrxfile,'-regexp','^(trx|timestamps)$');
  if ~isfield(td,'trx'),
    msg = sprintf('No trx variable in file %s',intrxfile);
    return;
  end
  trx = td.trx;
  if isfield(td,'timestamps'),
    timestamps = td.timestamps;
  else
    timestamps = [];
  end
  clear td;
catch ME,
  msg = sprintf('Could not load trx from file %s: %s',intrxfile,getReport(ME));
  return;
end

% compute timestamps
if isempty(timestamps),
  nframes = max([trx.endframe]);
  dt = zeros(1,nframes-1);
  n = zeros(1,nframes-1);
  for i = 1:numel(trx),
    dt(trx(i).firstframe:trx(i).endframe-1) = ...
      dt(trx(i).firstframe:trx(i).endframe-1) + trx(i).dt;
    n(trx(i).firstframe:trx(i).endframe-1) = ...
      n(trx(i).firstframe:trx(i).endframe-1) + 1;
  end
  dt(n>0) = dt(n>0) / n(n>0);
  mediandt = median(dt(dt>0));
  dt(n==0) = mediandt;
  timestamps = [0,cumsum(dt)]; %#ok<NASGU>
end

% create the experiment directory
if ~exist(expdir,'dir'),
  [success1,msg1] = mkdir(expdir);
  if ~success1,
    msg = msg1;
    return;
  end
end

% save the trx file
try
  save(trxfile,'trx','timestamps');
catch ME,
  msg = sprintf('Could not save to file %s: %s',trxfile,getReport(ME));
  return;
end
if ~exist(trxfile,'file'),
  msg = sprintf('Failed to save trx to file %s',trxfile);
  return;
end

% copy/soft-link movie
if strcmp(fullfile(inmoviefile),fullfile(moviefile)),
  fprintf('Input and out movie files are the same, not copying/linking.\n');
else

  if dosoftlink,
    if isunix,
      cmd = sprintf('ln -s %s %s',inmoviefile,moviefile);
      unix(cmd);
      % test to make sure it worked
      [status,result] = unix(sprintf('readlink %s',moviefile));
      result = strtrim(result);
      if status ~= 0 || ~strcmp(result,inmoviefile),
        warndlg(sprintf('Failed to make soft link, copying %s to %s instead',inmoviefile,moviefile));
        dosoftlink = false;
      end
    elseif ispc,
      cmd = sprintf('mkshortcut.vbs /target:"%s" /shortcut:"%s"',inmoviefile,moviefile);
      fprintf('Making a Windows shortcut file at "%s" with target "%s"\n',inmoviefile,moviefile);
      system(cmd);
      % test to make sure that worked
      [equalmoviefile,didfind] = GetPCShortcutFileActualPath(moviefile);
      if ~didfind || ~strcmp(equalmoviefile,inmoviefile),
        warndlg(sprintf('Failed to make shortcut, copying %s to %s instead',inmoviefile,moviefile));
        dosoftlink = false;
      end
    else
      warndlg(sprintf('Unknown OS, not soft-linking movie file %s',inmoviefile));
      dosoftlink = false;
    end
  end

  if ~dosoftlink,
    [success1,msg1] = copyfile(inmoviefile,moviefile);
    if ~success1,
      msg = msg1;
      success = false;
      return;
    end
  end
end

% make per-frame directory
if ~exist(perframedir,'dir'),
  [success1,msg1] = mkdir(perframedir);
  if ~success1,
    msg = msg1;
    return;
  end
end

success = true;

author: kristinbranson

repo: JAABA


text:
function hf_out = lhs_operation_joint(hf, samplesf, reg_filter, init_samplef, XH, init_hf, proj_reg)

% This is the left-hand-side operation in Conjugate Gradient

hf_out = cell(size(hf));

% Extract projection matrix and filter separately
P = cellfun(@real, hf(2,1,:), 'uniformoutput',false);
hf = hf(1,1,:);

% Get sizes
num_features = length(hf);
filter_sz = zeros(num_features,2);
for k = 1:num_features
    filter_sz(k,:) = [size(hf{k},1), size(hf{k},2)];
end
[~, k1] = max(filter_sz(:,1));  % Index for the feature block with the largest spatial size
block_inds = 1:num_features;
block_inds(k1) = [];
output_sz = [size(hf{k1},1), 2*size(hf{k1},2)-1];

% Compute the operation corresponding to the data term in the optimization
% (blockwise matrix multiplications)
%implements: A' diag(sample_weights) A f

% sum over all features and feature blocks
sh = sum(bsxfun(@times, samplesf{k1}, hf{k1}), 3);    % assumes the feature with the highest resolution is first
pad_sz = cell(1,1,num_features);
for k = block_inds
    pad_sz{k} = (output_sz - [size(hf{k},1), 2*size(hf{k},2)-1]) / 2;

    sh(1+pad_sz{k}(1):end-pad_sz{k}(1), 1+pad_sz{k}(2):end,1,:) = ...
        sh(1+pad_sz{k}(1):end-pad_sz{k}(1), 1+pad_sz{k}(2):end,1,:) + sum(bsxfun(@times, samplesf{k}, hf{k}), 3);
end

% weight all the samples and take conjugate
% sh = bsxfun(@times,sample_weights,sh);
sh = conj(sh);

% multiply with the transpose
hf_out1 = cell(1,1,num_features);
hf_out1{k1} = conj(sum(bsxfun(@times, sh, samplesf{k1}), 4));
for k = block_inds
    hf_out1{k} = conj(sum(bsxfun(@times, sh(1+pad_sz{k}(1):end-pad_sz{k}(1), 1+pad_sz{k}(2):end,1,1,:), samplesf{k}), 4));
end

% compute the operation corresponding to the regularization term (convolve
% each feature dimension with the DFT of w, and the tramsposed operation)
% add the regularization part
% hf_conv = cell(1,1,num_features);
for k = 1:num_features
    reg_pad = min(size(reg_filter{k},2)-1, size(hf{k},2)-1);

    % add part needed for convolution
    hf_conv = cat(2, hf{k}, conj(rot90(hf{k}(:, end-reg_pad:end-1, :), 2)));

    % do first convolution
    hf_conv = convn(hf_conv, reg_filter{k});

    % do final convolution and put toghether result
    hf_out1{k} = hf_out1{k} + convn(hf_conv(:,1:end-reg_pad,:), reg_filter{k}, 'valid');
end

% Stuff related to the projection matrix

% B * P
BP_cell = cell(1,1,num_features);
for k = 1:num_features
    BP_cell{k} = sum(bsxfun(@times, reshape(reshape(init_samplef{k}, [], size(init_samplef{k},3)) * P{k}, size(init_samplef{k},1), size(init_samplef{k},2), []), init_hf{k}), 3);
end

BP = BP_cell{k1};
for k = block_inds
    BP(1+pad_sz{k}(1):end-pad_sz{k}(1), 1+pad_sz{k}(2):end) = ...
        BP(1+pad_sz{k}(1):end-pad_sz{k}(1), 1+pad_sz{k}(2):end) + BP_cell{k};
end

% multiply with the transpose: A^H * BP
hf_out{1,1,k1} = hf_out1{k1} +  bsxfun(@times, BP, conj(samplesf{k1}));

% B^H * BP
fBP = cell(1,1,num_features);
fBP{k1} = reshape(bsxfun(@times, conj(init_hf{k1}), BP), [], size(init_hf{k1},3));

% Compute proj matrix part: B^H * A_m * f
shBP = cell(1,1,num_features);
shBP{k1} = reshape(bsxfun(@times, conj(init_hf{k1}), sh), [], size(init_hf{k1},3));

for k = block_inds
    % multiply with the transpose: A^H * BP
    hf_out{1,1,k} = hf_out1{k} +  bsxfun(@times, BP(1+pad_sz{k}(1):end-pad_sz{k}(1), 1+pad_sz{k}(2):end), conj(samplesf{k}));

    % B^H * BP
    fBP{k} = reshape(bsxfun(@times, conj(init_hf{k}), BP(1+pad_sz{k}(1):end-pad_sz{k}(1), 1+pad_sz{k}(2):end)), [], size(init_hf{k},3));

    % Compute proj matrix part: B^H * A_m * f
    shBP{k} = reshape(bsxfun(@times, conj(init_hf{k}), sh(1+pad_sz{k}(1):end-pad_sz{k}(1), 1+pad_sz{k}(2):end)), [], size(init_hf{k},3));
end

% hf_out2 = cell(1,1,num_features);
for k = 1:num_features
    fi = size(hf{k},1) * (size(hf{k},2)-1) + 1;  % index where the last frequency column starts

    % B^H * BP
    hf_out2 = 2*real(XH{k} * fBP{k} - XH{k}(:,fi:end) * fBP{k}(fi:end,:)) + proj_reg * P{k};

    % Compute proj matrix part: B^H * A_m * f
    hf_out{2,1,k} = hf_out2 + (2*real(XH{k} * shBP{k} - XH{k}(:,fi:end) * shBP{k}(fi:end,:)));
end
end

author: martin-danelljan

repo: ECO


text:
function	W = get_estimated_weight2(Model,parm,wmode);
% Return Full Weight matrix from 'Model'
%  W = get_estimated_weight2(Model,parm);
%  W = get_estimated_weight2(Model,parm,wmode);
% --- input
% Model.W : Weight matrix for active input [Ydim x Nactive]
% Model.ix_act : active input index
% parm.Dtau
% parm.M_all
% parm.xnorm
% parm.ynorm
%  In the default mode, normalization constants, xnorm and ynorm
%  are scaled back into the weight matrix 'W' : weight for original input
%  If wmode is given and wmode = 0,
%  no scale normalization is done: weight for normalized input
%
% --- output
% W : Weight matrix : 3D-array [Ydim x Xdim x Dtau)]
% W(n,m,:) : temporal weight for n-th output & m-th input data
% Dtau : time embedding dim
% Ydim : Output space dim
% Xdim : Input space dim
%
% 2008-5-20 Masa-aki Sato

if isfield(parm,'Dtau')
	Dtau  = parm.Dtau;
else
	Dtau  = 1;
end

M_all = Model.M_all;
Ydim  = size(Model.W,1);
Xdim  = M_all/Dtau;

if isfield(Model,'ix_act')
	% Active index
	ix_act = Model.ix_act;

	W = zeros(Ydim ,M_all);
	W(:,ix_act) = Model.W;
else
	W =  Model.W;
end

if ~exist('wmode','var') || wmode~=0,
	if length(parm.xmean) == Xdim,
		parm.xmean = repmat(parm.xmean ,[Dtau 1]);
		parm.xnorm = repmat(parm.xnorm ,[Dtau 1]);
	end

	% Scale back by normalization factor
	if isfield(parm,'xnorm') & isfield(parm,'ynorm')
		W = (parm.ynorm(:)*(1./parm.xnorm(:)')) .* W;
	end
end

W = reshape( W, [Ydim, Xdim, Dtau]);

author: KamitaniLab

repo: GenericObjectDecoding


text:
function [w,run] = train_bfgs(x,w,lambda)
% TRAIN_BFGS    Train a logistic regression model by BFGS.
%
% W = TRAIN_BFGS(X,W) returns maximum-likelihood weights given data and a
% starting guess.
% Data is columns of X, each column already scaled by the output (+1 or -1).
% W is the starting guess for the parameters (a column).

% Written by Thomas P Minka

if nargin < 3
  lambda = 0;
end
flops(0);
[d,n] = size(x);
old_g = zeros(size(w));
ih = eye(d);
for iter = 1:1000
  old_w = w;
  % s1 = 1-sigma
  s1 = 1./(1+exp(w'*x));
  g = x*s1' - lambda*w;
  flops(flops + flops_mul(w',x) + n*(flops_exp+2) + flops_mul(x,s1') + 2*d);
  if iter > 1
    dw = w - prev_w;
    dg = g - old_g;
    dwdg = dw'*dg;
    ihdg = ih*dg;
    b = 1 + (dg'*ihdg)/dwdg;
    ihdgdw = ihdg*dw';
    ih = ih + (b*dw*dw' - ihdgdw' - ihdgdw)/dwdg;
    flops(flops + d + d + flops_mul(dw',dg) + flops_mul(ih,dg) + ...
	flops_mul(dg',ihdg)+2 + flops_mul(ihdg,dw') + ...
	flops_mul(b,dw) + flops_mul(dw,dw') + 4*d*d);
  end
  u = -ih*g;
  flops(flops + flops_mul(ih,g));
  prev_w = w;

  % line search along u
  ug = u'*g;
  ux = u'*x;
  a = s1.*(1-s1);
  uhu = (ux.^2)*a' + lambda*(u'*u);
  w = w + (ug/uhu)*u;
  old_g = g;
  flops(flops + flops_mul(u',g) + flops_mul(u',x) + 2*n + ...
      n+flops_mul(1,n,1) + 2*d+1);
  if lambda > 0
    flops(flops + 1+flops_mul(u',u));
  end

  run.w(:,iter) = w;
  run.flops(iter) = flops;
  run.e(iter) = logProb(x,w) -0.5*lambda*w'*w;

  if max(abs(w - old_w)) < 1e-5
    break
  end
end
figure(2)
plot(run.e)
if iter == 1000
  warning('not enough iters')
end

author: FuzhenZhuang

repo: Transfer-Learning-Toolkit


text:
function moveSensor(speed, rot_punkt, tot_rot, theta, patch)

    n_patch = length(patch);

    % Loop through all sensors
    for i = 1:n_patch

        Points = get(patch(i),'Vertices');
        r = ones(size(Points));
        r = [speed*cos(tot_rot)*r(:,1), speed*sin(tot_rot)*r(:,2)] + Points(:,1:2);
        set(patch(i), 'XData', r(:,1), 'YData', r(:,2));
        rotate(patch(i),[0,0,1], theta*180/pi, [rot_punkt,0]);
    end

end

author: kennydl

repo: Reinforcment-Learning-With-Q-Learning


text:
function [sensor,num_Lst,InputFile,main_meta] = LoadSensorType(path_data)
%LOADSENSORTYPE Basic metadata should be loaded first to see which sensor
%here.
% Add Landsat 9 , Feb., 17, 2022

%% Search metadate file for Landsat 4-8
    main_meta=dir(fullfile(path_data,'L*MTL.txt'));
    existMTL=size(main_meta);
    InputFile=[];
    if existMTL(1)==0
        main_meta=dir(fullfile(path_data, 'S*MTD*TL*.xml'));
        if isempty(main_meta)
            main_meta=dir(fullfile(path_data, 'MTD*TL*.xml'));
        end
        if ~isempty(main_meta)
            txtstart = strfind(path_data,'S2A') ;
            num_Lst='2A';
            if isempty(txtstart)
                txtstart = strfind(path_data,'S2B') ;    % S2A or S2B  by Shi 10/18/2017
                num_Lst='2B';
            end
            txtend   = strfind(path_data,'.SAFE')-1 ;
            InputFile.DataStrip = path_data(txtstart:txtend) ;

            InputFile.Dmain = path_data(1:txtstart-1) ;

            txtstart = strfind(path_data,'GRANULE')+8 ;

            InputFile.InspireXML=path_data(1:txtstart-10); % extra add to read inspire xml.

            txtend   =  length(path_data);
            InputFile.Granule = path_data(txtstart:txtend) ;
            InputFile.pathh = path_data ;
            main_meta.name=InputFile.Granule;% return file name fro Fmask results.
        else
%             fprintf('No available data in the current folder!\n');
            sensor=[]; % no supported image.
            num_Lst=[];
            InputFile=[];
            main_meta=[];
            return;
        end
    else
        % determine sensor type
        % open and read hdr file
        fid_in=fopen(fullfile(path_data,main_meta.name),'r');
        geo_char=fscanf(fid_in,'%c',inf);
        fclose(fid_in);
        geo_char=geo_char';
        geo_str=strread(geo_char,'%s');

        % Identify Landsat Number (Lnum = 4, 5, 7, or 8)
        LID=char(geo_str(strmatch('SPACECRAFT_ID',geo_str)+2));
%         num_Lst=str2double(LID(end-1));
        num_Lst=(LID(end-1));
    end
    % define Landsat sensor.
    sensor='';
    if strcmp(num_Lst,'8') | strcmp(num_Lst,'9')
        sensor='L_OLI_TIRS';
    else
        if strcmp(num_Lst,'4')||strcmp(num_Lst,'5')||strcmp(num_Lst,'6')
            sensor='L_TM';
        else
            if strcmp(num_Lst,'7')
                sensor='L_ETM_PLUS';
            end
        end
    end
    if strcmp(num_Lst,'2A')||strcmp(num_Lst,'2B')
        sensor='S_MSI';
    end
end

author: GERSL

repo: Fmask

save_path: github-repos/MATLAB/GERSL-Fmask


text:
% pop_timtopo() - call the timtopo() function for epoched EEG datasets.
%                 Plots the epoch mean for each channel on a single axis,
%                 plus scalp maps of the data at specified latencies.
% Usage:
%   >> pop_timtopo( EEG, timerange, topotimes, title, 'key', 'val', ...);
%
% Inputs:
%   EEG         - input dataset
%   timerange   - [min max] epoch time range (in ms) to plot
%   topotimes   - array of times to plot scalp maps {Default: NaN
%                 = display scalp map at frame of max var()}
%
% Optional inputs:
%   title       - optional plot title
%   'key','val' - optional topoplot() arguments (see >> help topoplot)
%
% Author: Arnaud Delorme, CNL / Salk Institute, 2001
%
% See also: timtopo()

% Copyright (C) 2001 Arnaud Delorme, Salk Institute, arno@salk.edu
%
% This program is free software; you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation; either version 2 of the License, or
% (at your option) any later version.
%
% This program is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
% GNU General Public License for more details.
%
% You should have received a copy of the GNU General Public License
% along with this program; if not, write to the Free Software
% Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA

% 01-25-02 reformated help & license -ad
% 02-16-02 text interface editing -sm & ad
% 03-15-02 add all topoplot options -ad
% 03-18-02 added title -ad & sm

function com = pop_timtopo( EEG, timerange, topotime, plottitle, varargin);

com = '';
if nargin < 1
	help pop_timtopo;
	return;
end;

if nargin < 3
	promptstr    = { 'Plotting time range (ms):', ...
			         ['Scalp map latencies (ms, NaN -> max-RMS)'], ...
					 'Plot title:' ...
			         'Scalp map options (see >> help topoplot):' };
	inistr       = { [num2str( EEG.xmin*1000) ' ' num2str(EEG.xmax*1000)], ...
			         'NaN', ...
	                 ['ERP data and scalp maps' fastif(~isempty(EEG.setname), [' of ' EEG.setname ], '') ], ...
			         ''  };
	result       = inputdlg2( promptstr, 'ERP data and scalp maps -- pop_timtopo()', 1, inistr, 'pop_timtopo');
	if size(result,1) == 0 return; end;
	timerange    = eval( [ '[' result{1} ']' ] );
	topotime     = eval( [ '[' result{2} ']' ] );
	plottitle    = result{3};
	options      = [ ',' result{4} ];
	figure;
else
	options = [];
	for i=1:length( varargin )
		if isstr( varargin{ i } )
			options = [ options ', ''' varargin{i} '''' ];
		else
			options = [ options ', [' num2str(varargin{i}) ']' ];
		end;
	end;
end;
try, icadefs; set(gcf, 'color', BACKCOLOR, 'Name', ' timtopo()'); catch, end;

if exist('plottitle') ~= 1
    plottitle = ['ERP data and scalp maps' fastif(~isempty(EEG.setname), [' of ' EEG.setname ], '') ];
end;

if ~isempty(EEG.chanlocs)
    if ~isfield(EEG, 'chaninfo'), EEG.chaninfo = []; end;
	SIGTMP = reshape(EEG.data, size(EEG.data,1), EEG.pnts, EEG.trials);
	posi = round( (timerange(1)/1000-EEG.xmin)/(EEG.xmax-EEG.xmin) * (EEG.pnts-1))+1;
	posf = round( (timerange(2)/1000-EEG.xmin)/(EEG.xmax-EEG.xmin) * (EEG.pnts-1))+1;
	if length( options ) < 2
    	timtopo( mean(SIGTMP(:,posi:posf,:),3), EEG.chanlocs, [timerange(1) timerange(2) 0 0], topotime, '', 0, 0, 'chaninfo', EEG.chaninfo);
        com = sprintf('figure; pop_timtopo(%s, [%s], [%s], ''%s'');', inputname(1), num2str(timerange), num2str(topotime), plottitle);
	else
		com = sprintf('timtopo( mean(SIGTMP(:,posi:posf,:),3), EEG.chanlocs, [timerange(1) timerange(2) 0 0], topotime, '''', 0, 0, ''chaninfo'', EEG.chaninfo %s);', options);
		eval(com)
	    com = sprintf('figure; pop_timtopo(%s, [%s], [%s], ''%s'' %s);', inputname(1), num2str(timerange), num2str(topotime), plottitle, options);
	end;
else
	fprintf('Cannot make plot without channel locations\n');
	return;
end;
return;

author: PatternRecognition

repo: OpenBMI


text:
function o = boxoverlap(a, b)
% Compute the symmetric intersection over union overlap between a set of
% bounding boxes in a and a single bounding box in b.
%
% a  a matrix where each row specifies a bounding box
% b  a single bounding box

% AUTORIGHTS
% -------------------------------------------------------
% Copyright (C) 2011-2012 Ross Girshick
% Copyright (C) 2008, 2009, 2010 Pedro Felzenszwalb, Ross Girshick
%
% This file is part of the voc-releaseX code
% (http://people.cs.uchicago.edu/~rbg/latent/)
% and is available under the terms of an MIT-like license
% provided in COPYING. Please retain this notice and
% COPYING if you use this file (or a portion of it) in
% your project.
% -------------------------------------------------------

x1 = max(a(:,1), b(1));
y1 = max(a(:,2), b(2));
x2 = min(a(:,3), b(3));
y2 = min(a(:,4), b(4));

w = x2-x1+1;
h = y2-y1+1;
inter = w.*h;
aarea = (a(:,3)-a(:,1)+1) .* (a(:,4)-a(:,2)+1);
barea = (b(3)-b(1)+1) * (b(4)-b(2)+1);
% intersection over union overlap
o = inter ./ (aarea+barea-inter);
% set invalid entries to 0 overlap
o(w <= 0) = 0;
o(h <= 0) = 0;

author: rbgirshick

repo: rcnn

save_path: github-repos/MATLAB/rbgirshick-rcnn


text:
function graphToDot(adj, varargin)
% GRAPHTODOT Makes a GraphViz (AT&T) ile representing  an adjacency matrix
% function graphToDot(adj, ...)
% Optional arguments should be passed as name/value pairs [default]
%
% 'filename' - if omitted, writes to 'tmp.dot'
% 'arc_label' - arc_label{i,j} is a string attached to the i-j arc [""]
% 'node_label' - node_label{i} is a string attached to the node i ["i"]
% 'width'     - width in inches [10]
% 'height'    - height in inches [10]
% 'leftright' - 1 means layout left-to-right, 0 means top-to-bottom [0]
% 'directed'  - 1 means use directed arcs, 0 means undirected [1]
%
% For details on graphviz, See http://www.research.att.com/sw/tools/graphviz
%
% See also dot_to_graph and draw_dot
%
% First version written by Kevin Murphy 2002.
% Modified by Leon Peshkin, Jan 2004.

node_label = [];   arc_label = [];   % set default args
width = 10;        height = 10;
leftright = 0;     directed = 1;     filename = 'tmp.dot';

for i = 1:2:nargin-1                    % get optional args
    switch varargin{i}
        case 'filename', filename = varargin{i+1};
        case 'node_label', node_label = varargin{i+1};
        case 'arc_label', arc_label = varargin{i+1};
        case 'width', width = varargin{i+1};
        case 'height', height = varargin{i+1};
        case 'leftright', leftright = varargin{i+1};
        case 'directed', directed = varargin{i+1};
    end
end

fid = fopen(filename, 'w');
if directed
    fprintf(fid, 'digraph G {\n');
    arctxt = '->';
    if isempty(arc_label)
        labeltxt = '';
    else
        labeltxt = '[label="%s"]';
    end
else
    fprintf(fid, 'graph G {\n');
    arctxt = '--';
    if isempty(arc_label)
        labeltxt = '[dir=none]';
    else
        labeltext = '[label="%s",dir=none]';
    end
end
edgeformat = strcat(['%d ',arctxt,' %d ',labeltxt,';\n']);
fprintf(fid, 'center = 1;\n');
fprintf(fid, 'size=\"%d,%d\";\n', width, height);
if leftright
    fprintf(fid, 'rankdir=LR;\n');
end
Nnds = length(adj);
for node = 1:Nnds               %  process nodes
    if isempty(node_label)
        fprintf(fid, '%d;\n', node);
    else
        fprintf(fid, '%d [ label = "%s" ];\n', node,
node_label{node});
    end
end
for node1 = 1:Nnds   % process edges
    if directed
        arcs = find(adj(node1,:));         % children(adj, node);
    else
        arcs = find(adj(node1,node1+1:Nnds)); % remove duplicate arcs
    end
    for node2 = arcs
        fprintf(fid, edgeformat, node1, node2);
    end
end
fprintf(fid, '}');
fclose(fid);

author: bayesnet

repo: bnt

save_path: github-repos/MATLAB/bayesnet-bnt


text:
function check = semicircular_check ( a, b )

%*****************************************************************************80
%
%% SEMICIRCULAR_CHECK checks the parameters of the Semicircular CDF.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    20 September 2004
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, real A, B, the parameter of the PDF.
%    0.0 < B.
%
  if ( b <= 0.0 )
    fprintf ( 1, '\n' );
    fprintf ( 1, 'SEMICIRCULAR_CHECK - Fatal error!\n' );
    fprintf ( 1, '  B <= 0.0\n' );
    check = 0;
    return
  end

  check = 1;

  return
end

author: johannesgerer

repo: jburkardt-m


text:
function svgAnimation(s, type, key, value, duration)
for i = 1:length(s)
    userdata = get(s(i),'UserData');
    if isfield(userdata, 'svg') && isfield(userdata.svg, 'Animation')
        next = length(userdata.svg.Filter) + 1;
    else
        next = 1;
    end
    userdata.svg.Animation(next).SubAnimation.Type = type;
    userdata.svg.Animation(next).SubAnimation.Key = key;
    userdata.svg.Animation(next).SubAnimation.Value = value;
    userdata.svg.Animation(next).SubAnimation.Duration = duration;
    set(s(i),'UserData', userdata);
end

author: marianux

repo: ecg-kit

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


text:
function [ readerobj ] = open_rs_reader( filename )
%OPEN_RS_READER Intiates a reader object for RADARSAT-2 or RCM product.xml file
%
% Written by: Wade Schwartzkopf, NGA/IDT
%
% //////////////////////////////////////////
% /// CLASSIFICATION: UNCLASSIFIED       ///
% //////////////////////////////////////////

%% Read and process XML metadata
xp=javax.xml.xpath.XPathFactory.newInstance.newXPath();
domnode=xmlread(filename);
symmetry=[0 0 0];
% Because RS2/RCM XML uses a default namespace, we have to use a local-name
% query, rather than the much simpler traditional XPath query that assumes
% no namespace.  See isrs2.m for more details.
collector_name=char(xp.evaluate(['/*[local-name()=''product'']'...
    '/*[local-name()=''sourceAttributes'']'...
    '/*[local-name()=''satellite'']'], domnode));
if strcmpi(collector_name, 'RADARSAT-2')
    ia_str = 'imageAttributes';
    ipdf_str=['/*[local-name()=''product'']'...
        '/*[local-name()=''imageAttributes'']'...
        '/*[local-name()=''fullResolutionImageData'']'];
    beta_lut_str=['/*[local-name()=''product'']'...
        '/*[local-name()=''imageAttributes'']'...
        '/*[local-name()=''lookupTable''][@incidenceAngleCorrection="Beta Nought"]'];
    beta_lut_str = fullfile(fileparts(filename), char(xp.evaluate(beta_lut_str, domnode)));
    noise_str = ''; % No separate noise file for RS2
elseif strncmpi(collector_name, 'RCM', 3)
    ia_str = 'imageReferenceAttributes';
    ipdf_str=['/*[local-name()=''product'']'...
        '/*[local-name()=''sceneAttributes'']'...
        '/*[local-name()=''imageAttributes'']'...
        '/*[local-name()=''ipdf'']'];
    beta_lut_str=['/*[local-name()=''product'']'...
        '/*[local-name()=''imageReferenceAttributes'']'...
        '/*[local-name()=''lookupTableFileName''][@sarCalibrationType="Beta Nought"]'];
    beta_lut_str = fullfile(fileparts(filename), 'calibration', char(xp.evaluate(beta_lut_str, domnode)));
    noise_str=['/*[local-name()=''product'']'...
        '/*[local-name()=''imageReferenceAttributes'']'...
        '/*[local-name()=''noiseLevelFileName'']'];
    noise_str = fullfile(fileparts(filename), 'calibration', char(xp.evaluate(noise_str, domnode)));
end
lineOrder=char(xp.evaluate(['/*[local-name()=''product'']'...
    '/*[local-name()=''' ia_str ''']'...
    '/*[local-name()=''rasterAttributes'']'...
    '/*[local-name()=''lineTimeOrdering'']'],...
    domnode));
lookdir=upper(char(xp.evaluate(['/*[local-name()=''product'']'...
    '/*[local-name()=''sourceAttributes'']'...
    '/*[local-name()=''radarParameters'']'...
    '/*[local-name()=''antennaPointing'']'],...
    domnode)));
sampleOrder=char(xp.evaluate(['/*[local-name()=''product'']'...
    '/*[local-name()=''' ia_str ''']'...
    '/*[local-name()=''rasterAttributes'']'...
    '/*[local-name()=''pixelTimeOrdering'']'],...
    domnode));
symmetry(1)=xor(strcmp(lineOrder,'Decreasing'),lookdir(1)=='L');
symmetry(2)=strcmp(sampleOrder,'Decreasing');
if exist(beta_lut_str, 'file')
    betanode={xmlread(beta_lut_str)};
else
    betanode = {};
end
if exist(noise_str, 'file')
    try % Simulated RCM datasets have misformed noise XML
        betanode = {betanode{:} xmlread(noise_str)};
    end
end
meta=meta2sicd_rs_xml(domnode,betanode{:});

%% Get all files from a polarimetric collect
num_files=str2double(xp.evaluate(['count(' ipdf_str ')'],domnode));
readerobj={};
for i=1:num_files
    basepathname=fileparts(filename);
    datafilename=char(xp.evaluate([ipdf_str '[' num2str(i) ']'],domnode));
    % If format is needed, do this:
    % format=char(xp.evaluate(['/*[local-name()=''product'']'...
    %                          '/*[local-name()=''' ia_str ''']'...
    %                          '/*[local-name()=''productFormat'']']),...
    %                         domnode);
    % Return individual reader objects for each polarimetric channel
    if exist(fullfile(basepathname,datafilename),'file')
        readerobj{end+1}=open_tiff_reader_noxml(fullfile(basepathname,datafilename),symmetry);
        meta{i}.native.tiff=readerobj{end}.get_meta();
        readerobj{end}.get_meta = @() meta{i};
    end
end
end

% //////////////////////////////////////////
% /// CLASSIFICATION: UNCLASSIFIED       ///
% //////////////////////////////////////////

author: ngageoint


text:
function XRecovered = recoverData(reducedData, U, K)

% This function is used to recover the original data from the data set with
% reduced features. reducedData = (m x K), where m is the number of examples
% and K is the reduced number of features. U = (n x n). UReduced = (n x K),
% where K is the number of reduced features that we have. Therefore
% reducedData*(UReduced') becomes (m x n) , m examples made up of n features.

UReduced = U(:,1:K);
XRecovered = reducedData*(UReduced');

author: Sable

repo: mcbench-benchmarks


text:
% DEMCMU35GPLVMFGPLVM1 Learn a GPLVM on CMU 35 data set.

% FGPLVM

% Fix seeds
randn('seed', 1e5);
rand('seed', 1e5);

% Get the sequence numbers.
[Y, lbls] = lvmLoadData('cmu35WalkJog');
seq = cumsum(sum(lbls)) - [1:31];

dataSetName = 'cmu35gplvm';
experimentNo = 1;

% load data
[Y, lbls, Ytest, lblstest] = lvmLoadData(dataSetName);

% Set up model
options = fgplvmOptions('fitc');
options.optimiser = 'conjgrad';
options.back = 'mlp';
options.backOptions = mlpOptions(10);
options.numActive = 200;
options.fixInducing = 1;
options.fixIndices = round(linspace(1, size(Y, 1), options.numActive));
latentDim = 4;

d = size(Y, 2);
model = fgplvmCreate(latentDim, d, Y, options);
model.bias = mean(model.y);
model.scale = std(model.y);

% Add dynamics model.
optionsDyn = gpOptions('fitc');
optionsDyn.numActive = 200;
optionsDyn.fixInducing = 1;
optionsDyn.kern = kernCreate(model.X, {'rbf', 'white'});
optionsDyn.kern.comp{1}.inverseWidth = 0.2;
% This gives signal to noise of 0.1:5e-3 or 20:1.
optionsDyn.kern.comp{1}.variance = 0.01;
optionsDyn.kern.comp{2}.variance = 0.95;
diff = 1;
learn = 1;
optionsDyn.fixIndices = round(linspace(1, size(Y, 1)-length(seq), options.numActive));
model = fgplvmAddDynamics(model, 'gp', optionsDyn, diff, learn, seq);

% Optimise the model.
iters = 1000;
display = 1;

model = fgplvmOptimise(model, display, iters);

% Save the results.
modelWriteResult(model, dataSetName, experimentNo);

author: SheffieldML

repo: GPmat

save_path: github-repos/MATLAB/SheffieldML-GPmat


text:
function y = vl_nnloss_modified(x, c, varargin)
%VL_NNLOSS CNN categorical or attribute loss.
%   Y = VL_NNLOSS(X, C) computes the loss incurred by the prediction
%   scores X given the categorical labels C.
%
%   The prediction scores X are organised as a field of prediction
%   vectors, represented by a H x W x D x N array. The first two
%   dimensions, H and W, are spatial and correspond to the height and
%   width of the field; the third dimension D is the number of
%   categories or classes; finally, the dimension N is the number of
%   data items (images) packed in the array.
%
%   While often one has H = W = 1, the case W, H > 1 is useful in
%   dense labelling problems such as image segmentation. In the latter
%   case, the loss is summed across pixels (contributions can be
%   weighed using the `InstanceWeights` option described below).
%
%   The array C contains the categorical labels. In the simplest case,
%   C is an array of integers in the range [1, D] with N elements
%   specifying one label for each of the N images. If H, W > 1, the
%   same label is implicitly applied to all spatial locations.
%
%   In the second form, C has dimension H x W x 1 x N and specifies a
%   categorical label for each spatial location.
%
%   In the third form, C has dimension H x W x D x N and specifies
%   attributes rather than categories. Here elements in C are either
%   +1 or -1 and C, where +1 denotes that an attribute is present and
%   -1 that it is not. The key difference is that multiple attributes
%   can be active at the same time, while categories are mutually
%   exclusive. By default, the loss is *summed* across attributes
%   (unless otherwise specified using the `InstanceWeights` option
%   described below).
%
%   DZDX = VL_NNLOSS(X, C, DZDY) computes the derivative of the block
%   projected onto the output derivative DZDY. DZDX and DZDY have the
%   same dimensions as X and Y respectively.
%
%   VL_NNLOSS() supports several loss functions, which can be selected
%   by using the option `type` described below. When each scalar c in
%   C is interpreted as a categorical label (first two forms above),
%   the following losses can be used:
%
%   Classification error:: `classerror`
%     L(X,c) = (argmax_q X(q) ~= c). Note that the classification
%     error derivative is flat; therefore this loss is useful for
%     assessment, but not for training a model.
%
%   Top-K classification error:: `topkerror`
%     L(X,c) = (rank X(c) in X <= K). The top rank is the one with
%     highest score. For K=1, this is the same as the
%     classification error. K is controlled by the `topK` option.
%
%   Log loss:: `log`
%     L(X,c) = - log(X(c)). This function assumes that X(c) is the
%     predicted probability of class c (hence the vector X must be non
%     negative and sum to one).
%
%   Softmax log loss (multinomial logistic loss):: `softmaxlog`
%     L(X,c) = - log(P(c)) where P(c) = exp(X(c)) / sum_q exp(X(q)).
%     This is the same as the `log` loss, but renormalizes the
%     predictions using the softmax function.
%
%   Multiclass hinge loss:: `mhinge`
%     L(X,c) = max{0, 1 - X(c)}. This function assumes that X(c) is
%     the score margin for class c against the other classes.  See
%     also the `mmhinge` loss below.
%
%   Multiclass structured hinge loss:: `mshinge`
%     L(X,c) = max{0, 1 - M(c)} where M(c) = X(c) - max_{q ~= c}
%     X(q). This is the same as the `mhinge` loss, but computes the
%     margin between the prediction scores first. This is also known
%     the Crammer-Singer loss, an example of a structured prediction
%     loss.
%
%   When C is a vector of binary attribures c in (+1,-1), each scalar
%   prediction score x is interpreted as voting for the presence or
%   absence of a particular attribute. The following losses can be
%   used:
%
%   Binary classification error:: `binaryerror`
%     L(x,c) = (sign(x - t) ~= c). t is a threshold that can be
%     specified using the `threshold` option and defaults to zero. If
%     x is a probability, it should be set to 0.5.
%
%   Binary log loss:: `binarylog`
%     L(x,c) = - log(c(x-0.5) + 0.5). x is assumed to be the
%     probability that the attribute is active (c=+1). Hence x must be
%     a number in the range [0,1]. This is the binary version of the
%     `log` loss.
%
%   Logistic log loss:: `logistic`
%     L(x,c) = log(1 + exp(- cx)). This is the same as the `binarylog`
%     loss, but implicitly normalizes the score x into a probability
%     using the logistic (sigmoid) function: p = sigmoid(x) = 1 / (1 +
%     exp(-x)). This is also equivalent to `softmaxlog` loss where
%     class c=+1 is assigned score x and class c=-1 is assigned score
%     0.
%
%   Hinge loss:: `hinge`
%     L(x,c) = max{0, 1 - cx}. This is the standard hinge loss for
%     binary classification. This is equivalent to the `mshinge` loss
%     if class c=+1 is assigned score x and class c=-1 is assigned
%     score 0.
%
%   VL_NNLOSS(...,'OPT', VALUE, ...) supports these additionals
%   options:
%
%   InstanceWeights:: []
%     Allows to weight the loss as L'(x,c) = WGT L(x,c), where WGT is
%     a per-instance weight extracted from the array
%     `InstanceWeights`. For categorical losses, this is either a H x
%     W x 1 or a H x W x 1 x N array. For attribute losses, this is
%     either a H x W x D or a H x W x D x N array.
%
%   TopK:: 5
%     Top-K value for the top-K error. Note that K should not
%     exceed the number of labels.
%
%   See also: VL_NNSOFTMAX().

% Copyright (C) 2014-15 Andrea Vedaldi.
% Copyright (C) 2016 Karel Lenc.
% All rights reserved.
%
% This file is part of the VLFeat library and is made available under
% the terms of the BSD license (see the COPYING file).

if ~isempty(varargin) && ~ischar(varargin{1})  % passed in dzdy
    dzdy = varargin{1} ;
    varargin(1) = [] ;
else
    dzdy = [] ;
end

opts.marginMat = [] ;
opts.marginAlpha_ = 1 ;
opts.instanceWeights = [] ;
opts.classWeights = [] ;
opts.threshold = 0 ;
opts.loss = 'softmaxlog' ;
opts.topK = 5 ;
opts = vl_argparse(opts, varargin, 'nonrecursive') ;

inputSize = [size(x,1) size(x,2) size(x,3) size(x,4)] ;

% Form 1: C has one label per image. In this case, get C in form 2 or
% form 3.
c = gather(c) ;
if numel(c) == inputSize(4)
    c = reshape(c, [1 1 1 inputSize(4)]) ;
    c = repmat(c, inputSize(1:2)) ;
end

switch lower(opts.loss)
    case {'cosinesimilaritylogloss', 'cosinesimilaritymmloss', 'cosinesimilarityregloss', 'cosinesimilarityregloss_full', ...
            'cosinesimlaritylogloss', 'cosinesimlaritymmloss', 'cosinesimlarityregloss', 'cosinesimlarityregloss_full', ...
           'cosinesimlairtylogloss', 'cosinesimlairtymmloss', 'cosinesimlairtyregloss', 'cosinesimlairtyregloss_full', ...
           'cosinesimilarityabsregloss', 'cosinesimilarityabsregloss_full'}
        hasIgnoreLabel = 0;
    otherwise
        hasIgnoreLabel = any(c(:) == 0);
end

% --------------------------------------------------------------------
% Spatial weighting
% --------------------------------------------------------------------

% work around a bug in MATLAB, where native cast() would slow
% progressively
if isa(x, 'gpuArray')
    switch classUnderlying(x) ;
        case 'single', cast = @(z) single(z) ;
        case 'double', cast = @(z) double(z) ;
    end
    c = gpuArray(c);
else
    switch class(x)
        case 'single', cast = @(z) single(z) ;
        case 'double', cast = @(z) double(z) ;
    end
end

labelSize = [size(c,1) size(c,2) size(c,3) size(c,4)] ;
assert(isequal(labelSize(1:2), inputSize(1:2))) ;
assert(labelSize(4) == inputSize(4)) ;
instanceWeights = [] ;
switch lower(opts.loss)
    case {'classerror', 'topkerror', 'log', 'softmaxlog', 'mhinge', 'mshinge'}
        % there must be one categorical label per prediction vector
        assert(labelSize(3) == 1) ;

        if hasIgnoreLabel
            % null labels denote instances that should be skipped
            instanceWeights = cast(c(:,:,1,:) ~= 0) ;
        end

    case {'binaryerror', 'binarylog', 'logistic', 'hinge'}

        % there must be one categorical label per prediction scalar
        assert(labelSize(3) == inputSize(3)) ;

        if hasIgnoreLabel
            % null labels denote instances that should be skipped
            instanceWeights = cast(c ~= 0) ;
        end
    case {'cosinesimilaritylogloss', 'cosinesimilaritymmloss', 'cosinesimilarityregloss', 'cosinesimilarityregloss_full', ...
           'cosinesimlaritylogloss', 'cosinesimlaritymmloss', 'cosinesimlarityregloss', 'cosinesimlarityregloss_full', ...
           'cosinesimlairtylogloss', 'cosinesimlairtymmloss', 'cosinesimlairtyregloss', 'cosinesimlairtyregloss_full',...
           'cosinesimilarityabsregloss', 'cosinesimilarityabsregloss_full'}
        assert(labelSize(1) == inputSize(1)) ;
        assert(labelSize(2) == inputSize(2)) ;
        assert(labelSize(3) == inputSize(3)) ;
        assert(labelSize(4) == inputSize(4)) ;
        instanceWeights = opts.instanceWeights;
    otherwise
        error('Unknown loss ''%s''.', opts.loss) ;
end

if ~isempty(opts.instanceWeights)
    % important: this code needs to broadcast opts.instanceWeights to
    % an array of the same size as c
    if isempty(instanceWeights) && isempty(strfind(lower(opts.loss), 'cosinesimilarity'))  && isempty(strfind(lower(opts.loss), 'cosinesimlarity'))
        instanceWeights = bsxfun(@times, onesLike(c), opts.instanceWeights) ;
%         instanceWeights = c;
%         instanceWeights(instanceWeights~=11) = 3;
%         instanceWeights(instanceWeights==11) = 1;
%         instanceWeights = bsxfun(@times, instanceWeights, opts.instanceWeights) ;
    elseif isempty(instanceWeights) && ~isempty(strfind(lower(opts.loss), 'cosinesimilarity'))  && ~isempty(strfind(lower(opts.loss), 'cosinesimlarity'))
        instanceWeights = bsxfun(@times, instanceWeights, opts.instanceWeights);
    else
    end
end

% --------------------------------------------------------------------
% Do the work
% --------------------------------------------------------------------

switch lower(opts.loss)
    case {'log', 'softmaxlog', 'mhinge', 'mshinge'}
        % from category labels to indexes
        numPixelsPerImage = prod(inputSize(1:2)) ;
        numPixels = numPixelsPerImage * inputSize(4) ;
        imageVolume = numPixelsPerImage * inputSize(3) ;

        n = reshape(0:numPixels-1,labelSize) ;
        offset = 1 + mod(n, numPixelsPerImage) + ...
            imageVolume * fix(n / numPixelsPerImage) ;
        ci = offset + numPixelsPerImage * max(c - 1,0) ;
end

if nargin <= 2 || isempty(dzdy)
    switch lower(opts.loss)
        case 'classerror'
            [~,chat] = max(x,[],3) ;
            t = cast(c ~= chat) ;
        case 'topkerror'
            [~,predictions] = sort(x,3,'descend') ;
            t = 1 - sum(bsxfun(@eq, c, predictions(:,:,1:opts.topK,:)), 3) ;
        case 'log'
            t = - log(x(ci)) ;
        case 'softmaxlog'
            Xmax = max(x,[],3) ;
            ex = exp(bsxfun(@minus, x, Xmax)) ;
            t = Xmax + log(sum(ex,3)) - x(ci) ;
        case 'mhinge'
            t = max(0, 1 - x(ci)) ;
        case 'mshinge'
            Q = x ;
            Q(ci) = -inf ;
            t = max(0, 1 - x(ci) + max(Q,[],3)) ;
        case 'binaryerror'
            t = cast(sign(x - opts.threshold) ~= c) ;
        case 'binarylog'
            t = -log(c.*(x-0.5) + 0.5) ;
        case 'logistic'
            %t = log(1 + exp(-c.*X)) ;
            a = -c.*x ;
            b = max(0, a) ;
            t = b + log(exp(-b) + exp(a-b)) ;
        case 'hinge'
            t = max(0, 1 - c.*x) ;
        case 'cosinesimilaritylogloss'
            t = - (log(x).*c + log(1-x).*(1-c));
        case {'cosinesimilaritymmloss', 'cosinesimlaritymmloss', 'cosinesimlairtymmloss'}
            if isempty(opts.marginMat)
                t = max(0, x-opts.marginAlpha_) .* (1-c) ;
            else
                t = max(0, x-opts.marginMat) .* (1-c) ;
            end
        case {'cosinesimilarityregloss', 'cosinesimlarityregloss', 'cosinesimlairtyregloss'}
            t = c.* ((x-c).^2); % including positive pairs only
        case {'cosinesimilarityregloss_full', 'cosinesimlarityregloss_full', 'cosinesimlairtyregloss_full'}
            t = ((x-c).^2); % including positive and negative pairs
        case {'cosinesimilarityabsregloss', 'cosinesimlarityabsregloss', 'cosinesimlairtyabsregloss'}
            t = c.* (c-x); % including positive pairs only
        case {'cosinesimilarityabsregloss_full', 'cosinesimlarityabsregloss_full', 'cosinesimlairtyabsregloss_full'}
            t =  (c-x); % including positive pairs only
    end
    if ~isempty(instanceWeights)
        y = instanceWeights(:)' * t(:) ;
    else
        y = sum(t(:));
    end
else
    if ~isempty(instanceWeights)
        dzdy = dzdy * instanceWeights ;
    end
    switch lower(opts.loss)
        case {'classerror', 'topkerror'}
            y = zerosLike(x) ;
        case 'log'
            y = zerosLike(x) ;
            y(ci) = - dzdy ./ max(x(ci), 1e-8) ;
        case 'softmaxlog'
            Xmax = max(x,[],3) ;
            ex = exp(bsxfun(@minus, x, Xmax)) ;
            y = bsxfun(@rdivide, ex, sum(ex,3)) ;
            y(ci) = y(ci) - 1 ;
            y = bsxfun(@times, dzdy, y) ;
        case 'mhinge'
            % t = max(0, 1 - x(ci)) ;
            y = zerosLike(x) ;
            y(ci) = - dzdy .* (x(ci) < 1) ;
        case 'mshinge'
            Q = x ;
            Q(ci) = -inf ;
            [~, q] = max(Q,[],3) ;
            qi = offset + numPixelsPerImage * (q - 1) ;
            W = dzdy .* (x(ci) - x(qi) < 1) ;
            y = zerosLike(x) ;
            y(ci) = - W ;
            y(qi) = + W ;
        case 'binaryerror'
            y = zerosLike(x) ;
        case 'binarylog'
            y = - dzdy ./ (x + (c-1)*0.5) ;
        case 'logistic'
            % t = exp(-Y.*X) / (1 + exp(-Y.*X)) .* (-Y)
            % t = 1 / (1 + exp(Y.*X)) .* (-Y)
            y = - dzdy .* c ./ (1 + exp(c.*x)) ;
        case 'hinge'
            y = - dzdy .* c .* (c.*x < 1) ;
        case {'cosinesimilaritylogloss', 'cosinesimlaritylogloss', 'cosinesimlairtylogloss'}
            y = -dzdy .* ( c./x - (1-c) ./ (1-x) );
        case {'cosinesimilaritymmloss', 'cosinesimlaritymmloss', 'cosinesimlairtymmloss'}

            if isempty(opts.marginMat)
                % t = max(0, x-opts.marginAlpha_) .* (1-c) ;
                y = x - opts.marginAlpha_;
                y(y<0) = 0;
                y = y.* (1-c);
                y = dzdy .* y;
            else
                y = x - opts.marginMat;
                y(y<0) = 0;
                y = y.* (1-c);
                y = dzdy .* y;
            end
        case {'cosinesimilarityregloss', 'cosinesimlarityregloss', 'cosinesimlairtyregloss'}
            y = c.* dzdy .* (x-c);     % including positive pairs only
        case {'cosinesimilarityregloss_full', 'cosinesimlarityregloss_full', 'cosinesimlairtyregloss_full'}
            y = dzdy .* (x-c); % including positive and negative pairs
        case {'cosinesimilarityabsregloss', 'cosinesimlarityabsregloss', 'cosinesimlairtyabsregloss'}
            %t = c.* (c-x); % including positive pairs only
            y = c.* dzdy .* (-x);     % including positive pairs only
        case {'cosinesimilarityabsregloss_full', 'cosinesimlarityabsregloss_full', 'cosinesimlairtyabsregloss_full'}
            %t =  (c-x); % including positive pairs only
            y = dzdy .* (-x);     % including positive pairs only
    end
end

% --------------------------------------------------------------------
function y = zerosLike(x)
% --------------------------------------------------------------------
if isa(x,'gpuArray')
    y = gpuArray.zeros(size(x),classUnderlying(x)) ;
else
    y = zeros(size(x),'like',x) ;
end

% --------------------------------------------------------------------
function y = onesLike(x)
% --------------------------------------------------------------------
if isa(x,'gpuArray')
    y = gpuArray.ones(size(x),classUnderlying(x)) ;
else
    y = ones(size(x),'like',x) ;
end

author: aimerykong


text:
function [v0, vl, vu, vt] = params_var(om, name, idx)
%PARAMS_VAR  Returns initial value, lower bound and upper bound for opt variables.
%   [V0, VL, VU] = OM.PARAMS_VAR()
%   [V0, VL, VU] = OM.PARAMS_VAR(NAME)
%   [V0, VL, VU] = OM.PARAMS_VAR(NAME, IDX_LIST)
%   [V0, VL, VU, VT] = PARAMS_VAR(...)
%   Returns the initial value V0, lower bound VL and upper bound VU for
%   the full optimization variable vector, or for a specific named or named
%   and indexed variable set. Optionally also returns a corresponding char
%   vector VT of variable types, where 'C', 'I' and 'B' represent continuous
%   integer and binary variables, respectively.
%
%   Examples:
%       [x0, xmin, xmax] = om.params_var();
%       [Pg0, Pmin, Pmax] = om.params_var('Pg');
%       [zij0, zijmin, zijmax, ztype] = om.params_var('z', {i, j});
%
%   See also OPT_MODEL.

%   MP-Opt-Model
%   Copyright (c) 2008-2020, Power Systems Engineering Research Center (PSERC)
%   by Ray Zimmerman, PSERC Cornell
%
%   This file is part of MP-Opt-Model.
%   Covered by the 3-clause BSD License (see LICENSE file for details).
%   See https://github.com/MATPOWER/mp-opt-model for more info.

if nargout > 3
    have_vt = 1;
else
    have_vt = 0;
end
if nargin < 2
    v0 = []; vl = []; vu = []; vt = char([]);
    %% calls to substruct() are relatively expensive, so we pre-build the
    %% structs for addressing cell and numeric array fields, updating only
    %% the subscripts before use
    sc = struct('type', {'.', '{}'}, 'subs', {'', 1});  %% cell array field
    for k = 1:om.var.NS
        name = om.var.order(k).name;
        idx = om.var.order(k).idx;
        if isempty(idx)
            v0 = [ v0; om.var.data.v0.(name) ];
            vl = [ vl; om.var.data.vl.(name) ];
            vu = [ vu; om.var.data.vu.(name) ];
            if have_vt
                N = om.var.idx.N.(name);
                vt0 = om.var.data.vt.(name);
                if isscalar(vt0) && N > 1
                    vt = [ vt char(vt0 * ones(1, N)) ];
                else
                    vt = [ vt vt0 ];
                end
            end
        else
            % (calls to substruct() are relatively expensive ...
            % sc = substruct('.', name, '{}', idx);
            % ... so replace it with these more efficient lines)
            sc(1).subs = name;
            sc(2).subs = idx;
            v0 = [ v0; subsref(om.var.data.v0, sc) ];
            vl = [ vl; subsref(om.var.data.vl, sc) ];
            vu = [ vu; subsref(om.var.data.vu, sc) ];
            if have_vt
                % (calls to substruct() are relatively expensive ...
                % sn = substruct('.', name, '()', idx);
                % ... so replace it with these more efficient lines)
                sn = sc; sn(2).type = '()';
                N = subsref(om.var.idx.N, sn);
                vt0 = subsref(om.var.data.vt, sc);
                if isscalar(vt0) && N > 1
                    vt = [ vt char(vt0 * ones(1, N)) ];
                else
                    if ~isempty(vt0)
                        vt = [ vt vt0 ];
                    end
                end
            end
        end
    end
else
    if isfield(om.var.idx.N, name)
        if nargin < 3 || isempty(idx)
            v0 = om.var.data.v0.(name);
            vl = om.var.data.vl.(name);
            vu = om.var.data.vu.(name);
            if have_vt
                N = om.var.idx.N.(name);
                vt0 = om.var.data.vt.(name);
                if isscalar(vt0) && N > 1
                    vt = char(vt0 * ones(1, N));
                else
                    vt = vt0;
                end
            end
        else
            % (calls to substruct() are relatively expensive ...
            % sc = substruct('.', name, '{}', idx);
            % ... so replace it with these more efficient lines)
            sc = struct('type', {'.', '{}'}, 'subs', {name, idx});
            v0 = subsref(om.var.data.v0, sc);
            vl = subsref(om.var.data.vl, sc);
            vu = subsref(om.var.data.vu, sc);
            if have_vt
                % (calls to substruct() are relatively expensive ...
                % sn = substruct('.', name, '()', idx);
                % ... so replace it with these more efficient lines)
                sn = sc; sn(2).type = '()';
                N = subsref(om.var.idx.N, sn);
                vt0 = subsref(om.var.data.vt, sc);
                if isscalar(vt0) && N > 1
                    vt = char(vt0 * ones(1, N));
                else
                    vt = vt0;
                end
            end
        end
    else
        v0 = [];
        vl = [];
        vu = [];
        if have_vt
            vt = [];
        end
    end
end

author: MATPOWER

repo: matpower

save_path: github-repos/MATLAB/MATPOWER-matpower


text:
%USERKERNEL Construct user defined kernel mapping
%
%  K = USERKERNEL(B,R,FUNC,P1,P2, ...)
%  K = B*USERKERNEL([],R,FUNC,P1,P2, ...)
%  W = USERKERNEL([],R,FUNC,P1,P2, ...)
%  W = R*USERKERNEL([],[],FUNC,P1,P2, ...)
%  K = B*W
%
% INPUT
%   R     Dataset, representation set, default B
%   B     Dataset
%   FUNC  String with function name to compute kernels (proximities) by
%         K = FEVAL(FUNC,B,R,P1,P2, ...)
%
% OUTPUT
%   W     Trained kernel mapping
%   K     Kernel (proximity) matrix
%
% DESCRIPTION
% The kernel matrix K is computed according to the definition given by the
% user supplied function FUNC. The size of K is [SIZE(B,1) SIZE(R,1)]

% Copyright: R.P.W. Duin, r.p.w.duin@37steps.com
% Faculty EWI, Delft University of Technology
% P.O. Box 5031, 2600 GA Delft, The Netherlands

function w = userkernel(a,r,kernel,varargin)

	if nargin < 3, kernel = []; end
	if nargin < 2, r = []; end
	if nargin < 1, a = []; end

	if isempty(a) & isempty(r)
		w = prmapping(mfilename,'untrained',{[],kernel,varargin{:}});
		w = setname(w,'userkernel mapping');
	elseif isempty(r)
		[m,k] = size(a);
		w = prmapping(mfilename,'trained',{a,kernel,varargin},getlab(a),k,m);
	elseif isempty(a)
		[m,k] = size(r);
		w = prmapping(mfilename,'trained',{r,kernel,varargin},getlab(r),k,m);
	elseif ismapping(r) % execution of a*userkernel (trained)
		u = getdata(r);
		[r,kernel,pars] = deal(u{:});
		w = compute_kernel(kernel,a,r,pars);
	elseif isdataset(a) % execution
		w = compute_kernel(kernel,a,r,varargin);
	end

	function k = compute_kernel(kernel,a,r,pars)

		if isempty(kernel)
			error('No kernel function or kernel mapping supplied')
		elseif ~exist(kernel)
			error('kernel function not found')
		end

		if isempty(pars)
			k = feval(kernel,+a,+r);
		else
			k = feval(kernel,+a,+r,pars{:});
		end
		k = setdat(a,k);
		if isdataset(r)
			k = setfeatlab(k,getlabels(r));
		end
	return


text:
function vectors=VectorIndexing3D(tensor_in,idx)
% indexing vectors from a 3D tensor,
%2016-10-15 jlfeng
if  ndims(tensor_in)~=3
    return
end
dims=size(tensor_in);
data=reshape(tensor_in,[dims(1)*dims(2), dims(3)]);
vectors=data(idx(:),:);

author: zhangqianqianQQ

repo: MachineVisionAlgorithm


text:
function [ key ] = generate_checksum( X, k )
    sum_ = sum(X(:));
    key = [num2str(sum_) '_' num2str(k)];
end

author: jwyang


text:
function chebyshev_discrete_test ( )

%*****************************************************************************80
%
%% CHEBYSHEV_DISCRETE_TEST tests CHEBYSHEV_DISCRETE.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    18 March 2009
%
%  Author:
%
%    John Burkardt
%
  OFFSET = 1;

  fprintf ( 1, '\n' );
  fprintf ( 1, 'CHEBYSHEV_DISCRETE_TEST:\n' );
  fprintf ( 1, '  CHEBYSHEV_DISCRETE evaluates discrete Chebyshev polynomials.\n' );
  fprintf ( 1, '\n' );
  fprintf ( 1, '       N      M         X        T(N,M,X)\n' );

  n = 5;
  m = 5;

  for j = 0 : 5

    x = j / 2.0;

    value = chebyshev_discrete ( n, m, x );

    fprintf ( 1, '\n' );

    for i = 0 : n

      fprintf ( 1, '  %8d  %8d  %8f  %14f\n', i, m, x, value(i+OFFSET) );

    end

  end

  return
end


text:
% LFWriteMetadata - saves variables to a file in JSON format
%
% Usage:
%     LFWriteMetadata( JsonFileFname, DataToSave )
%
% This function saves data to a JSON file.
%
% User guide: LFToolbox.pdf
% See also: LFReadMetadata

% Copyright (c) 2013-2020 Donald G. Dansereau
%
% Based on code from JSONlab by Qianqian Fang,
% http://www.mathworks.com.au/matlabcentral/fileexchange/33381-jsonlab-a-toolbox-to-encodedecode-json-files-in-matlaboctave
%
% Minor modifications to simplify the interface as appropriate for use in the Light Field Toolbox,
% and to tolerate the sha keys appearing at the end of some Lytro JSON files.
% 2013, Donald G. Dansereau
%
% Enabled writing of class objects
% 2020, Nuno Monteiro

function LFWriteMetadata( JsonFileFname, DataToSave )

rootname = [];

varname=inputname(1);
opt.FileName = JsonFileFname;
opt.IsOctave=exist('OCTAVE_VERSION');
rootisarray=0;
rootlevel=1;
forceroot=jsonopt('ForceRootName',0,opt);
if((isnumeric(DataToSave) || islogical(DataToSave) || ischar(DataToSave) || isstruct(DataToSave) || iscell(DataToSave)) && isempty(rootname) && forceroot==0)
    rootisarray=1;
    rootlevel=0;
else
    if(isempty(rootname))
        rootname=varname;
    end
end
if((isstruct(DataToSave) || iscell(DataToSave))&& isempty(rootname) && forceroot)
    rootname='root';
end
json=obj2json(rootname,DataToSave,rootlevel,opt);
if(rootisarray)
    json=sprintf('%s\n',json);
else
    json=sprintf('{\n%s\n}\n',json);
end

jsonp=jsonopt('JSONP','',opt);
if(~isempty(jsonp))
    json=sprintf('%s(%s);\n',jsonp,json);
end

% save to a file if FileName is set, suggested by Patrick Rapin
if(~isempty(jsonopt('FileName','',opt)))
    fid = fopen(opt.FileName, 'wt');
    fwrite(fid,json,'char');
    fclose(fid);
end

%%-------------------------------------------------------------------------
function txt=obj2json(name,item,level,varargin)

if(iscell(item))
    txt=cell2json(name,item,level,varargin{:});
elseif(isstruct(item))
    txt=struct2json(name,item,level,varargin{:});
elseif(ischar(item))
    txt=str2json(name,item,level,varargin{:});
elseif(isobject(item))
    txt=class2json(name,item,level,varargin{:});
else
    txt=mat2json(name,item,level,varargin{:});
end

%%-------------------------------------------------------------------------
function txt=cell2json(name,item,level,varargin)
txt='';
if(~iscell(item))
        error('input is not a cell');
end

dim=size(item);
len=numel(item); % let's handle 1D cell first
padding1=repmat(sprintf('\t'),1,level-1);
padding0=repmat(sprintf('\t'),1,level);
if(len>1)
    if(~isempty(name))
        txt=sprintf('%s"%s": [\n',padding0, checkname(name,varargin{:})); name='';
    else
        txt=sprintf('%s[\n',padding0);
    end
elseif(len==0)
    if(~isempty(name))
        txt=sprintf('%s"%s": null',padding0, checkname(name,varargin{:})); name='';
    else
        txt=sprintf('%snull',padding0);
    end
end
for i=1:len
    txt=sprintf('%s%s%s',txt,padding1,obj2json(name,item{i},level+(len>1),varargin{:}));
    if(i<len) txt=sprintf('%s%s',txt,sprintf(',\n')); end
end
if(len>1) txt=sprintf('%s\n%s]',txt,padding0); end

%%-------------------------------------------------------------------------
function txt=struct2json(name,item,level,varargin)
txt = '';
if(~isstruct(item))
	error('input is not a struct');
end
txt = object2json(name,item,level,varargin{:});

%%-------------------------------------------------------------------------
function txt=class2json(name,item,level,varargin)
txt = '';
if(~isobject(item))
	error('input is not a class');
end
txt = object2json(name,item,level,varargin{:});

%%-------------------------------------------------------------------------
function txt=object2json(name,item,level,varargin)
txt='';
len=numel(item);
padding1=repmat(sprintf('\t'),1,level-1);
padding0=repmat(sprintf('\t'),1,level);
sep=',';

if(~isempty(name))
    if(len>1) txt=sprintf('%s"%s": [\n',padding0,checkname(name,varargin{:})); end
else
    if(len>1) txt=sprintf('%s[\n',padding0); end
end
for e=1:len
  names = fieldnames(item(e));
  if(~isempty(name) && len==1)
        txt=sprintf('%s%s"%s": {\n',txt,repmat(sprintf('\t'),1,level+(len>1)), checkname(name,varargin{:}));
  else
        txt=sprintf('%s%s{\n',txt,repmat(sprintf('\t'),1,level+(len>1)));
  end
  if(~isempty(names))
    for i=1:length(names)
	txt=sprintf('%s%s',txt,obj2json(names{i},getfield(item(e),...
             names{i}),level+1+(len>1),varargin{:}));
        if(i<length(names)) txt=sprintf('%s%s',txt,','); end
        txt=sprintf('%s%s',txt,sprintf('\n'));
    end
  end
  txt=sprintf('%s%s}',txt,repmat(sprintf('\t'),1,level+(len>1)));
  if(e==len) sep=''; end
  if(e<len) txt=sprintf('%s%s',txt,sprintf(',\n')); end
end
if(len>1) txt=sprintf('%s\n%s]',txt,padding0); end

%%-------------------------------------------------------------------------
function txt=str2json(name,item,level,varargin)
txt='';
if(~ischar(item))
        error('input is not a string');
end
item=reshape(item, max(size(item),[1 0]));
len=size(item,1);
sep=sprintf(',\n');

padding1=repmat(sprintf('\t'),1,level);
padding0=repmat(sprintf('\t'),1,level+1);

if(~isempty(name))
    if(len>1) txt=sprintf('%s"%s": [\n',padding1,checkname(name,varargin{:})); end
else
    if(len>1) txt=sprintf('%s[\n',padding1); end
end
isoct=jsonopt('IsOctave',0,varargin{:});
for e=1:len
    if(isoct)
        val=regexprep(item(e,:),'\\','\\');
        val=regexprep(val,'"','\"');
        val=regexprep(val,'^"','\"');
    else
        val=regexprep(item(e,:),'\\','\\\\');
        val=regexprep(val,'"','\\"');
        val=regexprep(val,'^"','\\"');
    end
    if(len==1)
        DataToSave=['"' checkname(name,varargin{:}) '": ' '"',val,'"'];
	if(isempty(name)) DataToSave=['"',val,'"']; end
        txt=sprintf('%s%s%s%s',txt,repmat(sprintf('\t'),1,level),DataToSave);
    else
        txt=sprintf('%s%s%s%s',txt,repmat(sprintf('\t'),1,level+1),['"',val,'"']);
    end
    if(e==len) sep=''; end
    txt=sprintf('%s%s',txt,sep);
end
if(len>1) txt=sprintf('%s\n%s%s',txt,padding1,']'); end

%%-------------------------------------------------------------------------
function txt=mat2json(name,item,level,varargin)
if(~isnumeric(item) && ~islogical(item))
        error('input is not an array');
end

padding1=repmat(sprintf('\t'),1,level);
padding0=repmat(sprintf('\t'),1,level+1);

if(length(size(item))>2 || issparse(item) || ~isreal(item) || jsonopt('ArrayToStruct',0,varargin{:}))
    if(isempty(name))
    	txt=sprintf('%s{\n%s"_ArrayType_": "%s",\n%s"_ArraySize_": %s,\n',...
              padding1,padding0,class(item),padding0,regexprep(mat2str(size(item)),'\s+',',') );
    else
    	txt=sprintf('%s"%s": {\n%s"_ArrayType_": "%s",\n%s"_ArraySize_": %s,\n',...
              padding1,checkname(name,varargin{:}),padding0,class(item),padding0,regexprep(mat2str(size(item)),'\s+',',') );
    end
else
    if(isempty(name))
    	txt=sprintf('%s%s',padding1,matdata2json(item,level+1,varargin{:}));
    else
        if(numel(item)==1 && jsonopt('NoRowBracket',1,varargin{:})==1)
            numtxt=regexprep(regexprep(matdata2json(item,level+1,varargin{:}),'^\[',''),']','');
           	txt=sprintf('%s"%s": %s',padding1,checkname(name,varargin{:}),numtxt);
        else
    	    txt=sprintf('%s"%s": %s',padding1,checkname(name,varargin{:}),matdata2json(item,level+1,varargin{:}));
        end
    end
    return;
end
dataformat='%s%s%s%s%s';

if(issparse(item))
    [ix,iy]=find(item);
    data=full(item(find(item)));
    if(~isreal(item))
       data=[real(data(:)),imag(data(:))];
       txt=sprintf(dataformat,txt,padding0,'"_ArrayIsComplex_": ','1', sprintf(',\n'));
    end
    txt=sprintf(dataformat,txt,padding0,'"_ArrayIsSparse_": ','1', sprintf(',\n'));
    if(find(size(item)==1))
        txt=sprintf(dataformat,txt,padding0,'"_ArrayData_": ',...
           matdata2json([ix,data],level+2,varargin{:}), sprintf('\n'));
    else
        txt=sprintf(dataformat,txt,padding0,'"_ArrayData_": ',...
           matdata2json([ix,iy,data],level+2,varargin{:}), sprintf('\n'));
    end
else
    if(isreal(item))
        txt=sprintf(dataformat,txt,padding0,'"_ArrayData_": ',...
            matdata2json(item(:)',level+2,varargin{:}), sprintf('\n'));
    else
        txt=sprintf(dataformat,txt,padding0,'"_ArrayIsComplex_": ','1', sprintf(',\n'));
        txt=sprintf(dataformat,txt,padding0,'"_ArrayData_": ',...
            matdata2json([real(item(:)) imag(item(:))],level+2,varargin{:}), sprintf('\n'));
    end
end
txt=sprintf('%s%s%s',txt,padding1,'}');

%%-------------------------------------------------------------------------
function txt=matdata2json(mat,level,varargin)
if(size(mat,1)==1)
    pre='';
    post='';
    level=level-1;
else
    pre=sprintf('[\n');
    post=sprintf('\n%s]',repmat(sprintf('\t'),1,level-1));
end
if(isempty(mat))
    txt='null';
    return;
end
floatformat=jsonopt('FloatFormat','%.10g',varargin{:});
formatstr=['[' repmat([floatformat ','],1,size(mat,2)-1) [floatformat sprintf('],\n')]];

if(nargin>=2 && size(mat,1)>1 && jsonopt('ArrayIndent',1,varargin{:})==1)
    formatstr=[repmat(sprintf('\t'),1,level) formatstr];
end
txt=sprintf(formatstr,mat');
txt(end-1:end)=[];
if(islogical(mat) && jsonopt('ParseLogical',0,varargin{:})==1)
   txt=regexprep(txt,'1','true');
   txt=regexprep(txt,'0','false');
end
%txt=regexprep(mat2str(mat),'\s+',',');
%txt=regexprep(txt,';',sprintf('],\n['));
% if(nargin>=2 && size(mat,1)>1)
%     txt=regexprep(txt,'\[',[repmat(sprintf('\t'),1,level) '[']);
% end
txt=[pre txt post];
if(any(isinf(mat(:))))
    txt=regexprep(txt,'([-+]*)Inf',jsonopt('Inf','"$1_Inf_"',varargin{:}));
end
if(any(isnan(mat(:))))
    txt=regexprep(txt,'NaN',jsonopt('NaN','"_NaN_"',varargin{:}));
end

%%-------------------------------------------------------------------------
function val=jsonopt(key,default,varargin)
val=default;
if(nargin<=2) return; end
opt=varargin{1};
if(isstruct(opt) && isfield(opt,key))
    val=getfield(opt,key);
end
%%-------------------------------------------------------------------------
function newname=checkname(name,varargin)
isunpack=jsonopt('UnpackHex',1,varargin{:});
newname=name;
if(isempty(regexp(name,'0x([0-9a-fA-F]+)_','once')))
    return
end
if(isunpack)
    isoct=jsonopt('IsOctave',0,varargin{:});
    if(~isoct)
        newname=regexprep(name,'(^x|_){1}0x([0-9a-fA-F]+)_','${native2unicode(hex2dec($2))}');
    else
        pos=regexp(name,'(^x|_){1}0x([0-9a-fA-F]+)_','start');
        pend=regexp(name,'(^x|_){1}0x([0-9a-fA-F]+)_','end');
        if(isempty(pos)) return; end
        str0=name;
        pos0=[0 pend(:)' length(name)];
        newname='';
        for i=1:length(pos)
            newname=[newname str0(pos0(i)+1:pos(i)-1) char(hex2dec(str0(pos(i)+3:pend(i)-1)))];
        end
        if(pos(end)~=length(name))
            newname=[newname str0(pos0(end-1)+1:pos0(end))];
        end
    end
end

repo: LFToolbox


text:
function handle = bvhVisualise(channels, skel, padding)

% BVHVISUALISE For updating a bvh representation of 3-D data.
% FORMAT
% DESC draws a bvh representation in a 3-D plot.
% ARG channels : the channels to update the skeleton with.
% ARG skel : the skeleton structure.
% RETURN handle : a vector of handles to the plotted structure.
%
% SEEALSO : bvhModify, skelVisualise
%
% COPYRIGHT : Neil D. Lawrence, 2005, 2006

% MOCAP

if nargin < 3
  padding = 0;
end
handle = skelVisualise(channels, skel, padding);


text:
%% Clear and Close Figures
clear ; close all; clc

correctPredictedSpam = 85
incorrectPredictedSpam = 890
correctPredictedNoSpam = 10
incorrectPredictedNoSpam = 15

total = correctPredictedSpam + incorrectPredictedSpam + correctPredictedNoSpam + incorrectPredictedNoSpam

accuracy = ( correctPredictedSpam + correctPredictedNoSpam ) / total
precision = correctPredictedSpam / (correctPredictedSpam + incorrectPredictedSpam)
recall = correctPredictedSpam / (correctPredictedSpam + incorrectPredictedNoSpam)
F1 = (2 * precision * recall) / (precision + recall)
F1

author: vugsus

repo: coursera-machine-learning


text:
classdef TP8 < PROBLEM
% <multi> <real> <large/none> <robust>
% Test problem for robust multi-objective optimization
% delta --- 0.05 --- Maximum disturbance degree
% H     ---   50 --- Number of disturbances

%------------------------------- Reference --------------------------------
% A. Gaspar-Cunha, J. Ferreira, and G. Recio, Evolutionary robustness
% analysis for multi-objective optimization: benchmark problems, Structural
% and Multidisciplinary Optimization, 2014, 49: 771-793.
%------------------------------- 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".
%--------------------------------------------------------------------------

    properties
        delta;      % Maximum disturbance degree
        H;          % Number of disturbances
    end
    methods
        %% Default settings of the problem
        function Setting(obj)
            [obj.delta,obj.H] = obj.ParameterSet(0.05,50);
            obj.M = 2;
            if isempty(obj.D); obj.D = 5; end
            obj.lower    = [0,0,-ones(1,obj.D-2)];
            obj.upper    = [1,1, ones(1,obj.D-2)];
            obj.encoding = ones(1,obj.D);
        end
        %% Calculate objective values
        function PopObj = CalObj(~,PopDec)
            PopObj(:,1) = PopDec(:,1);
            h = 2 - 0.8*exp(-((PopDec(:,2)-0.35)/0.25).^2) - exp(-((PopDec(:,2)-0.85)/0.03).^2);
            g = 50*sum(PopDec(:,3:end).^2,2);
            S = 1 - sqrt(PopObj(:,1));
            PopObj(:,2) = h.*(g+S);
        end
        %% Generate points on the Pareto front
        function R = GetOptimum(obj,N)
            R(:,1) = linspace(0,1,N)';
            R(:,2) = 1 - sqrt(R(:,1));
        end
        %% Generate the image of Pareto front
        function R = GetPF(obj)
            R = obj.GetOptimum(100);
        end
        %% Calculate the metric value
        function score = CalMetric(obj,metName,Population)
            switch metName
                case {'Mean_IGD','Mean_HV','Worst_IGD','Worst_HV'}
                    score = feval(metName,Population,obj);
                otherwise
                    score = feval(metName,Population,obj.optimum);
            end
        end
        %% Perturb solutions multiple times
        function PopX = Perturb(obj,PopDec,N)
            if nargin < 3; N = obj.H; end
            Delta = repmat(obj.delta.*(obj.upper-obj.lower),N*size(PopDec,1),1);
            w     = UniformPoint(N,obj.D,'Latin');
            Dec   = 2*Delta.*w(reshape(repmat(1:end,size(PopDec,1),1),1,[]),:) + repmat(PopDec,N,1) - Delta;
            Dec   = obj.CalDec(Dec);
            PopX  = SOLUTION(Dec,obj.CalObj(Dec),obj.CalCon(Dec));
        end
    end
end

author: BIMK

repo: PlatEMO

save_path: github-repos/MATLAB/BIMK-PlatEMO


text:
function [stat] = doMAPForward(B, S, param, stat)

if isempty(B)
    fprintf('Empty proposal set.\n');
    stat = [];
    return;
end
nB = size(B,2);
if nargin == 3 || isempty(stat)
    % initialization
    stat.W = [];
    stat.Xp = []; % optimal w_{ij} given the output set
    stat.X = zeros(nB,1); % assignment
    % construct W
    [stat.W, stat.Xp] = getW(B, S, param);
    stat.BGp = stat.Xp;
    stat.nms = zeros(size(B,2),1);
    stat.f = sum(stat.Xp);
    stat.O = [];
end

%% loop
while numel(stat.O) < min(param.maxnum, nB)
    V = max(stat.W - repmat(stat.Xp, [1 nB]),0);
    [score, vote] = max(sum(V) + stat.nms(:)');
    if score == 0 % no supporters
        break
    end
    tmpf = stat.f + score + param.phi;

    if (tmpf > stat.f)
        mask = V(:,vote) > 0;
        stat.X(mask) = vote;
        stat.O(end + 1) = vote;
        stat.Xp(mask) = stat.W(mask,vote);
        stat.f = tmpf;
        stat.nms = stat.nms + ...
            param.gamma*getNMSPenalty(B, B(:,vote));
    else
        break
    end
end


text:
%mg_cd  GMG preconditioner for convection-diffusion problem
%   IFISS scriptfile: DJS, HCE; 15 April 2005.
% Copyright (c) 2005 D.J. Silvester, H.C. Elman, A. Ramage

nc=log2(length(x)-1);
%
% compute new MG data or reload existing data?
compute_mg = default('compute / load MG data? 1/2 (default 1)',1);
if compute_mg==2
   load mgdata_cd
else
   h=2^(1-nc);
   fprintf('Setting up MG data ...')
   % top level
   mgdata(nc).matrix=Asupg;
   mgdata(nc).prolong=mg_prolong(2^nc,2^nc,x,y);
   xn=x(1:2:end);yn=y(1:2:end);
   % loop over remaining levels
   for level = nc-1:-1:2;
      mgdata(level).matrix=mg_cd_setup(xn,yn,viscosity,outbc);
      mgdata(level).prolong=mg_prolong(2^level,2^level,xn,yn);
      xn=xn(1:2:end);yn=yn(1:2:end);
   end
   fprintf('done\n')
   gohome, cd datafiles, save mgdata_cd mgdata
end
%
% MG parameters
smooth = default('Jacobi / Gauss-Seidel / ILU smoother? 1/2/3 (default is Gauss-Seidel)',2);
if smooth==3 % point ILU
   sweeps=1;stype=1;
elseif smooth==2 % Gauss-Seidel
   stype = default('point / line Gauss-Seidel? 1/2 (default is line)',2);
   if stype==2
      sweeps = default('number of Gauss-Seidel directions? 1/2/3/4 (default is 2)',2);
   else
      sweeps=1;
   end
else % point Jacobi
   sweeps=1;stype=1;
end
npre = default('number of pre-smoothing steps? (default is 1)',1);
npost = default('number of post-smoothing steps? (default is 1)',1);
%
% construct smoother
smooth_data = mg_smooth(mgdata,nc,sweeps,smooth,stype);


text:
function [Iq, A, W, S, index]=icassoShow(sR,varargin)
%function icassoShow(sR,['identifier1',value1,'indentifier2',value2,...])
%
%PURPOSE
%
%To generate explorative visualizations for Icasso
%
%EXAMPLES OF BASIC USAGE
%
%  [Iq, A, W, S]=icassoShow(sR);
%
%shows results for as many estimate-clusters as there are (reduced)
%data dimensions. Also, return the estimates of independent
%components estimates (A,W,S) that correspond to the centroid of
%each estimate-cluster. The first output Iq contain the quality of
%the  estimates. You can rank the estimates according to this index.
%
%  icassoShow(sR,'colorlimit',[0.7 0.9],'estimate','demixing','L',9);
%
%changes the color scale
%0...0.7 (not shown), 0.7...0.9 (light red), 0.9...1 (bright red)
%and suppresses the graph lines for similarities under value 0.7 in
%general, and inside clusters that are dense (0.9...1). Shows rows
%of demixing matrix instead of sources in the estimate
%window. Aggregate results in 9 estimate-clusters.
%
%INPUTS
%
%sR (struct) Icasso result data structure
%
%Optional input arguments are given as argument identifier - value
%pairs: 'identifier1', value1, 'identifier2', value2,...
%(case insensitive)
%
% 'L' (string) 'rdim' (default) | (integer)
%   sets the number of estimate-clusters 'rdim' sets it equal to the
%   (reduced) data dimension
% 'estimate' (string) 'source' (default) | 'demixing' | 'mixing' | 'off'
%   whether to show the estimates of the
%   - independent components (sources),
%   - rows of the demixing matrix (W), or
%   - columns of the mixing matrix (A)
%    that are associated to the centrotype of each estimate-cluster.
%   Argument 'off' suppresses the window.
% 'colorlimit' (vector) default [0.5 0.75 0.9]
%   sets the thresholds for color of graph lines and clusters
%   in the 2D plot; if the cluster density (the average
%   intra-cluster similarity) exceeds the highest value, the
%   cluster/lines will be bright red, if it is below the minimum,
%   the cluster is white/lines are suppressed. The rest is colored
%   with shades of red.
% 'line' (string) 'on' (default) | 'off'
%   whether to show the similarity graph lines in the 2D plot or not
% 'hull' (string) 'on' (default) | 'off'
%   whether to show the "cluster hulls" in the 2D plot or not
% 'graphlimit' (scalar) in 0...1 | (string) 'auto' (default)
%   Controls the 2D plot: See function icassoGraph
% 'dense' (scalar) in 0...1 | (string) 'auto' (default)
%   Controls the 2D plot: See function icassoGraph
% 'quality' (string) 'simple' (default) | 'detailed'
%   'simple' shows the cluster quality index, 'detailed' shows also
%   more detailed info in a separate window (figure 6)
%
%OUTPUT
%
% Iq    (vector) stability index of each estimate (see function
%         icassoResult)
% A     (matrix) estimated columns of the mixing matrix (A) =
%         pinv(W) (see function icassoResult)
% W     (matrix) estimated rows of the demixing matrix (W) (see
%         function icassoResult)
% S     (matrix) estimated independent components (see function
%         icassoResult)
% index (vector) indices to the centrotypes (centroids) of each
%         estimate-cluster (see function icassoResult)
%
%DETAILS
%
%Detailed explanation of the resulting figures can be found in help
%texts of the functions mentioned below:
%
%[Figures in brackets can be suppressed or are optional by default]
%
%Figure 1: Relative clustering quality index for different number
%of clusters and additional information: number of clusters used in
%the rest of figures (L), number of ICs, and (reduced) data
%dimension. See function icassoRindex.
%
%Figure 2: Stability (reliability) indices  of the selected L
%estimate-clusters. See function icassoStability.
%
%Figure 3: Correlation structure as a matrix and a dendrogram
%representation for L clusters. See function icassoDendrogram.
%
%Figure 4: Graph of the correlations between all the estimates and
%estimate-clusters as convex hulls for L clusters
%(see function icassoGraph)
%
%[Figure 5: IC estimates, rows of W, or columns of A (depends on user
%selection) that correspond to the centroid (actually centrotype)
%of the selected estimate-clusters. See function signalplot.]
%
%[Figure 6: more detailed statistics on the estimate
%clusters: components of stability index Iq. See function
%icassoStability.]
%
%SEE ALSO
% icassoGet
% icassoResult
% signalplot
% icassoStability
% icassoViz
% icasso

%COPYRIGHT NOTICE
%This function is a part of Icasso software library
%Copyright (C) 2003-2005 Johan Himberg
%
%This program is free software; you can redistribute it and/or
%modify it under the terms of the GNU General Public License
%as published by the Free Software Foundation; either version 2
%of the License, or any later version.
%
%This program is distributed in the hope that it will be useful,
%but WITHOUT ANY WARRANTY; without even the implied warranty of
%MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
%GNU General Public License for more details.
%
%You should have received a copy of the GNU General Public License
%along with this program; if not, write to the Free Software
%Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.

% ver 1.21 040305 johan

if nargin<1|isempty(sR),
  error('At least one input argument expected');
end

if isempty(sR.projection.coordinates) | isempty(sR.cluster.partition) | ...
      isempty(sR.cluster.index) | isempty(sR.cluster.similarity),
  error('Missing similarity/projection/cluster information.');
end

% initiate output args.
index2centrotypes=[]; clusterquality=[]; partition=[];

%% Set defaults and process optional input
default={'line','on','estimate','source', 'quality','simple','L',icassoGet(sR,'rdim'),...
	 'graphlimit','auto','colorlimit',[0.5 0.75 0.9],...
	 'dense','auto','hull','on'};

% initiate arguments to icassoGraph
graphArgs=[];

varargin=processvarargin(varargin,default);
num_of_args=length(varargin);

for i=1:2:num_of_args,
  id=varargin{i}; value=varargin{i+1};
  switch lower(id)
   % Check first icassoShow
   case 'quality'
    switch lower(value)
     case 'simple'
      detailedrankplot=0;
     case 'detailed'
      detailedrankplot=1;
     otherwise
      error('Option ''quality'' must be ''simple'' or ''detailed''.');
    end
   case 'estimate'
    switch lower(value)
     case {'demixing','mixing','off','source'}
      est=lower(value);
     otherwise
      error(['Option ''estimate'' must be ''source'',''demixing'',' ...
	     ' ''mixing'', or ''off''.']);
    end
   case 'l'
    if isnumeric(value);
      level=value;
    else
      switch lower(value)
       case 'rdim'
	level=icassoGet(sR,'rdim');
       otherwise
	error(['Option ''L'' must be an integer or string ''rdim''']);
      end
    end
    % submit also to icassoGraph
    graphArgs{1,end+1}=id;
    graphArgs{1,end+1}=value;
    % submit the following to icassoGraph
   case {'line','graphlimit','colorlimit','dense','hull'}
    graphArgs{1,end+1}=id;
    graphArgs{1,end+1}=value;
   otherwise
    error(['Option ''' lower(id) ''' not available.' sprintf('\n') ...
	   'Available: ''L'', ''line'',''graphlimit'',''colorlimit'',''hull'',' ...
	   '''estimate'', and ''quality''.']);
  end
end

% Check if cluster level is valid
maxCluster=size(sR.cluster.partition,1);
if level<=0 | level>maxCluster,
  error('Cluster level out of range or not specified.');
end

%%%% Get main results %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

[Iq, A, W, S, index2centrotypes]=icassoResult(sR,level);

%%%%% Compute some cluster statistics %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% Get the partition

partition=sR.cluster.partition(level,:);
Ncluster=max(partition);

% cluster statistics
c=sR.cluster.similarity;
s=clusterstat(c,partition);

%%%%%%%%%%%%%%%%% Clustering validity index %%%%%%%%%%%%%%%%%%%%%%%%%%%

figure(1); clf reset;

icassoRindex(sR,level);

%%%%%%%%%%%%%%%% Ranking clusters & number of estimates %%%%%%%%%%%%%%%

figure(2); clf;

% compute & plot quality index

subplot(1,2,1);
Iq=icassoStability(sR,level,'plotindex');

% compute the rank order of estimates
clusterlabels=1:Ncluster;
[tmp,estimateOrder]=sort(-Iq);

% plot number of estimates in each cluster
subplot(1,2,2);
barh(s.N(estimateOrder));
set(gca,'ytick',1:Ncluster,'yticklabel', ...
	clusterlabels(estimateOrder),'ydir','reverse');
text(s.N(estimateOrder),1:Ncluster,cellstr(num2str(s.N(estimateOrder)')));
axis([0 max(s.N) 0.5 Ncluster+.5]);
title('Number of ICA estimates in the estimate-clusters');
ylabel('Label');
xlabel('Number of estimates');

set(2,'name','Icasso: Estimate Quality');

%%%%%%%%%%%%%%% Dendrogram %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Show dendrogram

figure(3);
clf reset;
icassoDendrogram(sR,level);

%%%%%%%%%%%%%%% Correlation graph %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

figure(4); clf reset;
if ~isempty(sR.projection.coordinates),
  icassoGraph(sR,graphArgs{:});
else
  warning(['Projection coordinates not computed, can''t start ' ...
           'icassoGraph.']);
end

%%%%%%%%%%%%% Source plots %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

switch est
 case 'source'
  figure(5); clf reset;
  set(5,'name','Icasso: Source Estimates (centrotypes)');
  signalplot(S(estimateOrder,:));
  set(gca,'yticklabel',estimateOrder);
  ylabel('Label');
  xlabel('Sample #');
  title(['Independent components (ranked according to I_q)']);
 case 'demixing'
  figure(5); clf reset;
  set(5,'name','Icasso: Demixing Matrix Estimate');
  signalplot(W(estimateOrder,:));
  set(gca,'yticklabel',estimateOrder);
  ylabel('Label');
  xlabel('Column #');
  title(['Demixing matrix rows (ranked according to' ...
	 ' I_q)']);
 case 'mixing'
  figure(5); clf reset;
  set(5,'name','Icasso: Mixing Matrix Estimate');
  signalplot(A(:,estimateOrder)');
  set(gca,'yticklabel',estimateOrder);
  ylabel('Label');
  xlabel('Row #');
  view([-90 90])
  title(['Mixing matrix columns (ranked according to' ...
	 ' I_q)']);
end

%%%%%%%%%%%%% Details of Iq %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

if detailedrankplot,
  figure(6); clf;
  set(6,'name','Icasso: Detaild Estimate Stability');
  icassoStability(sR,level,'plotstat');
end

%%% subfunctions

author: fieldtrip

repo: fieldtrip

save_path: github-repos/MATLAB/fieldtrip-fieldtrip


text:
function varargout = scan(varargin)
% SCAN(f, varargin) is not supported for ballfuns.

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

error('CHEBFUN:BALLFUN:SLICE:fail',...
       'Scan is not supported for ballfuns.')

end

author: chebfun

repo: chebfun

save_path: github-repos/MATLAB/chebfun-chebfun


text:
function res = reduce(im,bfilt)
%
% function res = reduce(im,bfilt)
%
% Separable convolution and subsampling by a factor of two
% im: input image
% bfilt: convolution kernel (vector).
%       Default is: [.0625 .25 .375 .25 .0625]'
% res: result image
%
% res is 1/2 the size of im.  Fills NaNs for invalid pixels near
% edges.

if ~exist('bfilt')
  bfilt=[.0625 .25 .375 .25 .0625]';
end

% Use standard Matlab convolution routines instead for ease of
% distribution, and to set edge values to NaNs
bsize = floor(length(bfilt)/2);
tmp1 = conv2sep(im,bfilt,bfilt,'valid');
tmp2 = NaN*ones(size(im));
tmp2(1+bsize:size(im,1)-bsize,1+bsize:size(im,2)-bsize)=tmp1;
res = tmp2(1:2:size(im,1),1:2:size(im,2));

return;

%%%%%%%%%
% Debug %
%%%%%%%%%

foo=pgmRead('einstein.pgm');
bar=reduce(foo);
bar=replaceValue(bar,NaN,0);
displayImage(bar);

in=ones(7,7)
res=reduce(in)

author: vistalab

repo: vistasoft

save_path: github-repos/MATLAB/vistalab-vistasoft


text:
function [s,d,S_k,S_si] = pinHoleSegment(k,si)

% PINHOLESEGMENT  Pin hole projection of a segment.
%   [S,D] = PINHOLESEGMENT(K,SI) projects into a pinhole camera with
%   intrinsic aprameters K the segment SI. It returns the projected segment
%   S and the non-observable depths D of the two endpoints.
%
%   SI is a 6-vector containing the two endpoints of the 3D segment.
%   S is a 4-vector conteining the two endpoints of the 2D segment.
%
%   SI can also be a 6-by-N matrix with N segments. In such case S is a
%   4-by-N matrix with N projected segments.
%
%   [S,D,S_k,S_si] = POINHOLESEGMENT(...) returns the Jacobians of S wrt K
%   and SI. It only works for single segments.
%
%   See also PINHOLE.

%   Copyright 2008-2009 Joan Sola @ LAAS-CNRS.

p1 = si(1:3,:);
p2 = si(4:6,:);

if nargout <= 2

    [e1,d1] = pinHole(p1,k);
    [e2,d2] = pinHole(p2,k);
    s       = [e1;e2];
    d       = [d1;d2];

else % Jacobians

    if size(si,2) == 1

        [e1,d1,E1_p1,E1_k] = pinHole(p1,k);
        [e2,d2,E2_p2,E2_k] = pinHole(p2,k);
        s       = [e1;e2];
        d       = [d1;d2];

        S_k = [...
            E1_k
            E2_k];

        Z23 = zeros(2,3);

        S_si = [...
            E1_p1 Z23
            Z23   E2_p2];

    else

        error('Jacobians not available for multiple segments.')

    end

end

return

%% test
syms p1 p2 p3 q1 q2 q3 u0 v0 au av real
k  = [u0;v0;au;av];
si = [p1;p2;p3;q1;q2;q3];

[s,d,S_k,S_si] = pinHoleSegment(k,si);

simplify(S_k - jacobian(s,k))
simplify(S_si - jacobian(s,si))

% ========== End of function - Start GPL license ==========

%   # START GPL LICENSE

%---------------------------------------------------------------------
%
%   This file is part of SLAMTB, a SLAM toolbox for Matlab.
%
%   SLAMTB 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.
%
%   SLAMTB 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 SLAMTB.  If not, see <http://www.gnu.org/licenses/>.
%
%---------------------------------------------------------------------

%   SLAMTB is Copyright:
%   Copyright (c) 2008-2010, Joan Sola @ LAAS-CNRS,
%   Copyright (c) 2010-2013, Joan Sola,
%   Copyright (c) 2014-2015, Joan Sola @ IRI-UPC-CSIC,
%   SLAMTB is Copyright 2009
%   by Joan Sola, Teresa Vidal-Calleja, David Marquez and Jean Marie Codol
%   @ LAAS-CNRS.
%   See on top of this file for its particular copyright.

%   # END GPL LICENSE

author: joansola

repo: slamtb

save_path: github-repos/MATLAB/joansola-slamtb


text:
classdef TestEM
    %TestEM

    properties (Constant)
        X = [randn(10,3)+1; randn(10,3)-1; randn(10,3)-2];
    end

    methods (Static)
        function test_clustering
            [N,d] = size(TestEM.X);
            C = 3;
            model = cv.EM();
            model.ClustersNumber = C;
            model.CovarianceMatrixType = 'Diagonal';
            model.TermCriteria.maxCount = 100;
            model.TermCriteria.epsilon = 1e-6;

            model.train(TestEM.X);
            assert(model.isTrained());

            [LL, IDX, probs] = model.trainEM(TestEM.X);
            validateattributes(LL, {'numeric'}, {'vector', 'numel',N});
            validateattributes(IDX, {'numeric'}, ...
                {'vector', 'integer', '<',C, 'numel',N});
            validateattributes(probs, {'numeric'}, ...
                {'size',[N C], 'real', '>=',0, '<=',1});

            weights = model.getWeights();
            validateattributes(weights, {'numeric'}, {'vector', 'numel',C});

            means = model.getMeans();
            validateattributes(means, {'numeric'}, {'size',[C d]});

            covs = model.getCovs();
            validateattributes(covs, {'cell'}, {'vector'});
            if ~mexopencv.isOctave()
                %HACK: https://savannah.gnu.org/bugs/index.php?46328
                cellfun(@(covar) validateattributes(covar, {'numeric'}, ...
                    {'square', 'diag', 'size',[d d]}), covs);
            end

            [LL, IDX, probs] = model.trainE(TestEM.X, means);
            [LL, IDX, probs] = model.trainE(TestEM.X, means, 'Covs0',covs);
            [LL, IDX, probs] = model.trainE(TestEM.X, means, ...
                'Covs0',covs, 'Weights0',weights);

            [LL, IDX, probs] = model.trainM(TestEM.X, probs);

            probs = model.predict(TestEM.X);
            validateattributes(probs, {'numeric'}, ...
                {'size',[N C], 'real', '>=',0, '<=',1});

            [LL, IDX, probs] = model.predict2(TestEM.X);
            validateattributes(LL, {'numeric'}, {'vector', 'numel',N});
            validateattributes(IDX, {'numeric'}, ...
                {'vector', 'integer', '<',C, 'numel',N});
            validateattributes(probs, {'numeric'}, ...
                {'size',[N C], 'real', '>=',0, '<=',1});
        end

        function test_data_options
            model = cv.EM();
            N = size(TestEM.X, 1);
            model.train(TestEM.X, 'Data',{'Layout','Row', 'VarType','NNN', ...
                'VarIdx',[], 'SampleIdx',[], 'SampleWeights',ones(N,1), ...
                'TrainTestSplitRatio',1/3, 'TrainTestSplitShuffle',true});
        end

        function test_storage
            fname = tempname();
            model = cv.EM();
            model.train(TestEM.X);

            model.save([fname '.xml']);
            cleanObj = onCleanup(@() delete([fname '.xml']));
            model1 = cv.EM();
            model1.load([fname '.xml']);
            %isequal(model, model1)

            model.save([fname '.yaml']);
            cleanObj = onCleanup(@() delete([fname '.yaml']));
            model2 = cv.EM([fname '.yaml']);
            %isequal(model, model2)

            model1.clear();
            model2.clear();
        end

        function test_serialization
            model = cv.EM();
            model.train(TestEM.X);
            strXML = model.save('.xml');
            strYML = model.save('.yml');

            model2 = cv.EM();
            model2.load(strXML, 'FromString',true);
            Yhat = model2.predict(TestEM.X);

            model3 = cv.EM(strYML, 'FromString',true);
            Yhat = model3.predict(TestEM.X);
        end
    end

end

author: kyamagu

repo: mexopencv

save_path: github-repos/MATLAB/kyamagu-mexopencv


text:
% Test file for @chebfun/sum.m.

function pass = test_sum(pref)

% Obtain preferences.
if ( nargin == 0 )
    pref = chebfunpref();
end

% Generate a few random points in [-1 1] to use as test values.
seedRNG(7681);
xr = 2 * rand(1000, 1) - 1;

%% SCALAR-VALUED

% Check the empty case.
pass(1) = sum(chebfun()) == 0;

% Check operation in the general case.
f = chebfun({@(x) exp(4*pi*1i*x), @exp, @exp}, [-1 0 0.5 1], pref);
pass(2) = abs(sum(f) - (exp(1) - 1)) < 10*vscale(f)*eps;

% Check operation for row chebfuns.
ft = f.';
pass(3) = abs(sum(ft) - (exp(1) - 1)) < 10*vscale(ft)*eps;

% Check sum over a subdomain.
pass(4) = abs(sum(f, [-1 1]) - (exp(1) - 1)) < 10*vscale(f)*eps;
pass(5) = abs(sum(f, [-1 0])) < 10*vscale(f)*eps;
pass(6) = abs(sum(f, [0 1]) - (exp(1) - 1)) < 10*vscale(f)*eps;

% Check sum between chebfun limits.
f = chebfun(@exp, [-1 -0.5 0 0.5 1], pref);
a = chebfun(@(x) x.^2 - 1, [-1 1]);
b = chebfun(@(x) -x.^2 + 1, [-1 1]);

F1 = sum(f, a, 1);
F1_exact = @(x) exp(1) - exp(x.^2 - 1);
pass(7) = norm(feval(F1, xr) - F1_exact(xr), inf) < 10*vscale(F1)*eps;

F2 = sum(f, -1, b);
F2_exact = @(x) exp(-x.^2 + 1) - exp(-1);
pass(8) = norm(feval(F2, xr) - F2_exact(xr), inf) < 10*vscale(F2)*eps;

F3 = sum(f, a, b);
F3_exact = @(x) exp(-x.^2 + 1) - exp(x.^2 - 1);
pass(9) = norm(feval(F3, xr) - F3_exact(xr), inf) < 10*vscale(F3)*eps;

%% ARRAY-VALUED
f = chebfun(@(x) [sin(x) cos(x) exp(x)], [-1 -0.5 0 0.5 1]);
pass(10) = norm(sum(f) - [0 2*sin(1) (exp(1) - exp(-1))], inf) < ...
    10*vscale(f)*eps;

ft = f.';
pass(11) = norm(sum(ft) - [0 2*sin(1) (exp(1) - exp(-1))].', inf) < ...
    10*vscale(ft)*eps;

pass(12) = norm(sum(f, [-1 1]) - [0 2*sin(1) (exp(1) - exp(-1))], inf) < ...
    10*vscale(f)*eps;

pass(13) = norm(sum(f, [-1 0]) - [(cos(-1) - 1) sin(1) (1 - exp(-1))], inf) ...
    < 10*vscale(f)*eps;

F1 = sum(f, a, 1);
F1_col1_exact = @(x) cos(x.^2 - 1) - cos(1);
F1_col2_exact = @(x) sin(1) - sin(x.^2 - 1);
F1_col3_exact = @(x) exp(1) - exp(x.^2 - 1);
F1_exact = @(x) [F1_col1_exact(x) F1_col2_exact(x) F1_col3_exact(x)];
err = feval(F1, xr) - F1_exact(xr);
pass(14) = norm(err(:), inf) < 10*vscale(F1)*eps;

F2 = sum(f, -1, b);
F2_col1_exact = @(x) cos(-1) - cos(-x.^2 + 1);
F2_col2_exact = @(x) sin(-x.^2 + 1) - sin(-1);
F2_col3_exact = @(x) exp(-x.^2 + 1) - exp(-1);
F2_exact = @(x) [F2_col1_exact(x) F2_col2_exact(x) F2_col3_exact(x)];
err = feval(F2, xr) - F2_exact(xr);
pass(15) = norm(err(:), inf) < 10*vscale(F2)*eps;

F3 = sum(f, a, b);
F3_col1_exact = @(x) -cos(-x.^2 + 1) + cos(x.^2 - 1);
F3_col2_exact = @(x) sin(-x.^2 + 1) - sin(x.^2 - 1);
F3_col3_exact = @(x) exp(-x.^2 + 1) - exp(x.^2 - 1);
F3_exact = @(x) [F3_col1_exact(x) F3_col2_exact(x) F3_col3_exact(x)];
err = feval(F3, xr) - F3_exact(xr);
pass(16) = norm(err(:), inf) < 1e2*vscale(F3)*eps;

%% Check dim argument.
g = sum(f, 2);
g_exact = @(x) sin(x) + cos(x) + exp(x);
pass(17) = norm(feval(g, xr) - g_exact(xr), inf) < 10*vscale(g)*eps;

g = sum(ft, 1);
g_exact = @(x) (sin(x) + cos(x) + exp(x)).';
pass(18) = norm(feval(g, xr) - g_exact(xr), inf) < 10*vscale(g)*eps;

%% Check error conditions.
try
    s = sum(f, -2, 2);
    pass(19) = false;
catch ME
    pass(19) = strcmp(ME.identifier, 'CHEBFUN:CHEBFUN:sum:sumSubDom:ab');
end

try
    s = sum(f, -2, b);
    pass(20) = false;
catch ME
    pass(20) = strcmp(ME.identifier, 'CHEBFUN:CHEBFUN:sum:sumSubDom:a');
end

try
    s = sum(f, a, 2);
    pass(21) = false;
catch ME
    pass(21) = strcmp(ME.identifier, 'CHEBFUN:CHEBFUN:sum:sumSubDom:b');
end

%% QUASIMATRICES:

f = quasimatrix(@(x) [sin(x) cos(x) exp(x)], [-1 -0.5 0 0.5 1]);
pass(22) = norm(sum(f) - [0 2*sin(1) (exp(1) - exp(-1))], inf) < ...
    10*vscale(f)*eps;

ft = f.';
pass(23) = norm(sum(ft) - [0 2*sin(1) (exp(1) - exp(-1))].', inf) < ...
    10*vscale(ft)*eps;

pass(24) = norm(sum(f, [-1 1]) - [0 2*sin(1) (exp(1) - exp(-1))], inf) < ...
    10*vscale(f)*eps;

pass(25) = norm(sum(f, [-1 0]) - [(cos(-1) - 1) sin(1) (1 - exp(-1))], inf) ...
    < 10*vscale(f)*eps;

%% Test on singular function: piecewise smooth chebfun - splitting on.

% Set a domain
dom = [-2 7];

pow = -0.5;
op = @(x) (x - dom(1)).^pow.*sin(100*x);
f = chebfun(op, dom, 'exps', [pow 0], 'splitting', 'on');
I = sum(f);
I_exact = 0.17330750941063138;
pass(26) = ( abs(I-I_exact) < 1e3*eps*abs(I_exact) );

%% Test for functions defined on unbounded domain:

% Functions on [-inf inf]:

% Set the domain:
dom = [-Inf 2 Inf];

op1 = @(x) x.^2.*exp(-x.^2);
op2 = @(x) (1-exp(-x.^2))./x.^2;
f = chebfun({op1 op2}, dom);
I = sum(f);
% The exact solution is obtained using Matlab symbolic toolbox:
IExact = 1.364971769155161;
err = abs(I - IExact);
pass(27) = err < 1e7*eps*get(f,'vscale');

% Functions on [0 inf]:

% Set the domain:
dom = [0 Inf];

op1 = @(x) x.*exp(-x);
f = chebfun(op1, dom);
I1 = sum(f);

IExact = 1;
err1 = abs(I1 - IExact);
pass(28) = err1 < 2e5*eps*get(f,'vscale');

x = chebfun('x', dom);
g = x.*exp(-x);
warning('off', 'CHEBFUN:UNBNDFUN:sum:slowDecay');
I2 = sum(g);
warning('on', 'CHEBFUN:UNBNDFUN:sum:slowdDecay');
err2 = abs(I2 - IExact);
tol = 2e10*eps*get(f,'vscale');
pass(29) = err2 < tol;

% Function on [-Inf Inf]:
f = chebfun('exp(-x.^2/16).*(1+.2*cos(10*x))',[-inf,inf]);

% Suppress expected warnings which may occur on certain machines:
warning('off','CHEBFUN:UNBNDFUN:sum:slowDecay');
I = sum(f);
% Re-enable warnings:
warning('on','CHEBFUN:UNBNDFUN:sum:slowDecay');
IExact = 7.0898154036220641;
err = abs(I - IExact);
pass(30) = err < 1e9*eps*get(f,'vscale');

% #496
k = 3;
op = @(t) exp(-t)./(1+k*t);
f = chebfun(op, [0 inf]);
I = sum(f);
% The following exact result is obtained using Matlab symbolic toolbox:
I_exact = 0.385602012136694;
err = abs(I-I_exact);
pass(31) = err < 2e1*eps*get(f,'vscale');

% #920: sum of array-valued chebfun defined on unbounded domain:
% (same function which decays fast enough to be integrable):
f = chebfun(@(x) exp(-[x x].^2), [0 Inf]);
I = sum(f);
% The following exact value is obtained by Mathematica.
I_exact = 0.886226925452758*ones(1, 2);
pass(32) = norm(I-I_exact, inf) < 1e1*eps*get(f,'vscale');

% #920: sum of array-valued chebfun defined on unbounded domain:
% (different functions but both decay fast so that integrable)
f = chebfun(@(x)[exp(-x.^2) 1./(x.^3)], [1 Inf]);
I = sum(f);
% The following exact value is obtained by Mathematica.
I_exact = [0.139402792640331 0.5];
pass(33) = norm(I-I_exact, inf) < 1e6*eps*get(f,'vscale');

% #920: sum of array-valued chebfun defined on unbounded domain:
% (two different functions one of which is integrable and the other one is not):
f = chebfun(@(x)[exp(-x.^2) 0*x+1], [1 Inf]);
I = sum(f);
% The following exact value is obtained by Mathematica.
I_exact = 0.139402792640331;
pass(34) = abs(I(1)-I_exact) < 1e1*eps*get(f,'vscale') && ...
    isinf(I(2));

%%

% Complex limits are not supported
f = chebfun(@sin);
try
    sum(f, [0 1i]);
    pass(35) = false;
catch
    pass(35) = true;
end

% Invalid subdomains in reverse order
t = chebfun(@(t) t, [0 1]);
try
    sum(t, 0, -1)
    pass(36) = false;
catch
    pass(36) = true;
end

end


text:
function [d1, d2, sel] = partition(d, sel, p)

% PARTITION - partition a DATA_GLOBAL object randomly or deterministically
%
% D1 = PARTITION(D, SEL)
% [D1, D2] = PARTITION(D, SEL)
% [D1, D2, SEL] = PARTITION(D, 'rand', P)
%
% The same as the DATA/PARTITION method, except that it operates on
% DATA_GLOBAL objects, manipulating D.index, and D.myX and D.myY if
% present.

n = length(d.index);

sel = sel(:);
if ischar(sel)
    if ~strcmp(lower(sel'), 'rand'), error(['unknown option ''' sel' '''']), end
    if nargin < 3, p = 0.5; end
    if ~isnumeric(p) | prod(size(p)) ~= 1, error('P must be a numeric scalar'), end
    p = full(real(double(p)));
    if p < 0
        error('P may not be negative')
    elseif p > n
        error('P may not be larger than the number of data points')
    elseif p < 1
        sel = (rand(n, 1) < p);
    elseif p ~= round(p)
        error('P must either be a probability, or a positive integer')
    else
        sel = randperm(n)';
        sel = sel(1:p);
    end
end
if islogical(sel)
    if length(sel) ~= n
        error('logical array must have one element per data point')
    end
elseif ~isreal(sel) | ~isnumeric(sel)
    error('SEL should be a logical array or an array of valid indices')
else
    sel = full(double(sel));
    if any(sel < 1 | sel > n | sel ~= round(sel))
        error('invalid/out of range indices')
    end
    z = logical(zeros(n,1));
    z(sel) = 1;
    sel = z;
end

d1 = d;
d1.index(~sel) = [];
if ~isempty(d1.myX), d1.myX(~sel, :) = []; end
if ~isempty(d1.myY), d1.myY(~sel, :) = []; end

if nargout < 2, return, end

d2 = d;
d2.index(sel) = [];
if ~isempty(d2.myX), d2.myX(sel, :) = []; end
if ~isempty(d2.myY), d2.myY(sel, :) = []; end

author: canlab

repo: CanlabCore

save_path: github-repos/MATLAB/canlab-CanlabCore


text:
% LSADM (Goldfarb et al. 2010)
% process_video('RPCA', 'LSADM', 'dataset/demo.avi', 'output/demo_LSADM.avi');
opts.D = M;
opts.mu = norm(M)/1.25;
[n1,n2] = size(M);
opts.Xs = M;
opts.Ys = M;
opts.n1 = n1;
opts.n2 = n2;
opts.sigma = 1e-6;
opts.maxitr = 500;
opts.rho = 1/sqrt(n1);
opts.eta_mu = 2/3;
opts.eta_sigma = 2/3;
opts.muf = 1e-6;
opts.sigmaf = 1e-6;
opts.epsilon = 1e-7;
opts.sv = 100;
out_ALM = ALM_SADAL_smoothed(opts.D,opts);
L = out_ALM.X;
S = out_ALM.Y;

author: andrewssobral

repo: lrslibrary