text:
function a = lock ( n )

%*****************************************************************************80
%
%% LOCK returns the number of codes for a lock with N buttons.
%
%  Discussion:
%
%    A button lock has N numbered buttons.  To open the lock, groups
%    of buttons must be pressed in the correct order.  Each button
%    may be pushed no more than once.  Thus, a code for the lock is
%    an ordered list of the groups of buttons to be pushed.
%
%    For this discussion, we will assume that EVERY button is pushed
%    at some time, as part of the code.  To count the total number
%    of codes, including those which don't use all the buttons, then
%    the number is 2 * A(N), or 2 * A(N) - 1 if we don't consider the
%    empty code to be valid.
%
%  Examples:
%
%    If there are 3 buttons, then there are 13 possible "full button" codes:
%
%      (123)
%      (12) (3)
%      (13) (2)
%      (23) (1)
%      (1) (23)
%      (2) (13)
%      (3) (12)
%      (1) (2) (3)
%      (1) (3) (2)
%      (2) (1) (3)
%      (2) (3) (1)
%      (3) (1) (2)
%      (3) (2) (1)
%
%    and, if we don't need to push all the buttons, every "full button" code above
%    yields a distinct "partial button" code by dropping the last set of buttons:
%
%      ()
%      (12)
%      (13)
%      (23)
%      (1)
%      (2)
%      (3)
%      (1) (2)
%      (1) (3)
%      (2) (1)
%      (2) (3)
%      (3) (1)
%      (3) (2)
%
%  First values:
%
%     N         A(N)
%     0           1
%     1           1
%     2           3
%     3          13
%     4          75
%     5         541
%     6        4683
%     7       47293
%     8      545835
%     9     7087261
%    10   102247563
%
%  Recursion:
%
%    A(I) = sum ( 0 <= J < I ) Binomial ( I, N-J ) * A(J)
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    02 June 2007
%
%  Author:
%
%    John Burkardt
%
%  Reference:
%
%    Daniel Velleman, Gregory Call,
%    Permutations and Combination Locks,
%    Mathematics Magazine,
%    Volume 68, Number 4, October 1995, pages 243-253.
%
%  Parameters:
%
%    Input, integer N, the maximum number of lock buttons.
%
%    Output, integer A(1:N+1), the number of lock codes.
%
  if ( n < 0 )
    a = [];
    return
  end

  a(1) = 1;

  for i = 1 : n
    a(i+1) = 0;
    for j = 0 : i-1
      a(i+1) = a(i+1) + i4_choose ( i, i-j ) * a(j+1);
    end
  end

  return
end


text:
function [document, scores] = report_failures(context, experiment, trackers, sequences, varargin)
% report_failures Generate a statistic overview of occurences of failures
%
% Generate a statistic overview of occurences of failures in individual sequences.
%
% Input:
% - context (structure): Report context structure.
% - experiment (struct): An experiment structure.
% - trackers (cell): An array of tracker structures.
% - sequences (cell): An array of sequence structures.
%
% Output:
% - document (structure): Resulting document structure.
%

scores = [];

document = create_document(context, 'failures', 'title', 'Failure overview');

for i = 1:2:length(varargin)
    switch lower(varargin{i})
        otherwise
            error(['Unknown switch ', varargin{i},'!']) ;
    end
end

tracker_labels = cellfun(@(x) iff(isfield(x.metadata, 'verified') && x.metadata.verified, [x.label, '*'], x.label), trackers, 'UniformOutput', 0);

sequences_hash = md5hash(strjoin((cellfun(@(x) x.name, sequences, 'UniformOutput', false)), '-'), 'Char', 'hex');
trackers_hash = md5hash(strjoin((cellfun(@(x) x.identifier, trackers, 'UniformOutput', false)), '-'), 'Char', 'hex');

if ~strcmp(experiment.type, 'supervised')
   error('Failure analysis only suitable for supervised experiments!');
end

cache_identifier = sprintf('failures_%s_%s_%s.mat', experiment.name, trackers_hash, sequences_hash);

histograms = report_cache(context, cache_identifier, @analyze_failures, ...
    experiment, trackers, sequences);

document.subsection('Experiment %s', experiment.name);

experiment_sequences = convert_sequences(sequences, experiment.converter);

for s = 1:length(experiment_sequences)

    if numel(experiment_sequences{s}.tags.names) == 0
        continue
    end

    document.subsection('Sequence %s',  experiment_sequences{s}.name);

    document.raw('
\n');

    hf = figure('Visible', 'off');
    imagesc(histograms{s});
    set(gca,'ytick',(1:numel(trackers)),'yticklabel', tracker_labels);

    document.figure(hf, sprintf('failures_%s_%s_individual', experiment.name, experiment_sequences{s}.name), ...
        sprintf('Sequence %s, individual trackers', experiment_sequences{s}.name));

    hf = figure('Visible', 'off');
    tags = experiment_sequences{s}.tags.data;
    tagsplit = mat2cell(tags, size(tags, 1), ones(1, size(tags, 2)));
    starts = cellfun(@(x) find(diff([0; x; 0]) > 0), tagsplit, 'UniformOutput', 0);
    ends = cellfun(@(x) find(diff([0; x; 0]) < 0), tagsplit, 'UniformOutput', 0);

    hold on;
    plot_timeline(experiment_sequences{s}.tags.names, starts, ends);
    combined_histogram = sum(histograms{s}, 1);
    combined_histogram = (combined_histogram / max(combined_histogram)) * numel(experiment_sequences{s}.tags.names);
    plot(combined_histogram);
    set(gca, 'XLim', [1, experiment_sequences{s}.length]);
    hold off;

    document.figure(hf, sprintf('failures_%s_%s_combined', experiment.name, experiment_sequences{s}.name), ...
        sprintf('Sequence %s, combined failures with properties', experiment_sequences{s}.name));

    document.raw('
\n');
end;

document.write();

author: votchallenge

repo: toolkit-legacy


text:
function varargout = drawCurve3d(varargin)
%DRAWCURVE3D draw a 3D curve specified by a list of points.
%
%   drawCurve3d(COORD) packs coordinates in a single [N*3] array.
%
%   drawCurve3d(PX, PY, PZ) specify coordinates in separate arrays.
%
%   H = drawCurve3d(...) also return a handle to the list of line objects.
%
%   See Also :
%   drawPolygon

% ------
% Author: David Legland
% e-mail: david.legland@inrae.fr
% Created: 2005-02-18
% Copyright 2005 INRA - TPV URPOI - BIA IMASTE

% deprecation warning
warning('geom3d:deprecated', ...
    '''drawCurve3d'' is deprecated, use ''drawPolyline3d'' instead');

% default value for closed or open curve
closed = false;

% check case we want to draw several curves, stored in a cell array
var = varargin{1};
if iscell(var)
    hold on;
    h = [];
    for i = 1:length(var(:))
        h = [h; drawCurve3d(var{i}, varargin{2:end})]; %#ok<AGROW>
    end
    if nargout>0
        varargout{1}=h;
    end
    return;
end

% extract curve coordinate
if size(var, 2)==1
    % first argument contains x coord, second argument contains y coord
    % and third one the z coord
    px = var;
    if length(varargin)<3
        error('Wrong number of arguments in drawCurve3d');
    end
    py = varargin{2};
    pz = varargin{3};
    varargin = varargin(4:end);
else
    % first argument contains both coordinate
    px = var(:, 1);
    py = var(:, 2);
    pz = var(:, 3);
    varargin = varargin(2:end);
end

% check if curve is closed or open
if ~isempty(varargin)
    var = varargin{1};
    if strncmpi(var, 'close', 5)
        closed = true;
        varargin = varargin(2:end);
    elseif strncmpi(var, 'open', 4)
        closed = false;
        varargin = varargin(2:end);
    end
end

% for closed curve, add the first point at the end to close curve
if closed
    px = [px; px(1)];
    py = [py; py(1)];
    pz = [pz; pz(1)];
end

%% draw the curve ! !! ! !
h = plot3(px, py, pz, varargin{:});

if nargout>0
    varargout{1} = h;
end

author: mattools

repo: matGeom

save_path: github-repos/MATLAB/mattools-matGeom


text:
function s = json2struct(j)

j = j(2:end-1); % remove the { and }
t = tokenize(j, ',');

s = struct();
for i=1:numel(t)
  kv = tokenize(t{i}, ':');
  key = kv{1}(2:end-1); % remove the " and "
  val = kv{2};
  if any(val(1)=='0123456789') && ~isnan(str2double(val))
    val = str2double(val);
  elseif isequal(val, 'true')
    val = true;
  elseif isequal(val, 'false')
    val = false;
  else
    val = val(2:end-1); % remove the " and "
  end
  key = fixname(key);
  s.(key) = val;
end

author: fieldtrip

repo: fieldtrip

save_path: github-repos/MATLAB/fieldtrip-fieldtrip


text:
function nearest = find_closest2 ( m, nr, r, ns, s )

%*****************************************************************************80
%
%% FIND_CLOSEST2 finds the nearest R point to each S point.
%
%  Discussion:
%
%    We are given R, a set of NR points in M dimensions.
%
%    We are given S, a set of NS points in M dimensions.
%
%    For each S(I) in S, we seek the index J of the point R(J)
%    which is nearest to S(I) over all points in R.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    08 July 2010
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, integer M, the spatial dimension.
%
%    Input, integer NR, the number of cell generators.
%
%    Input, real R(M,NR), the cell generators.
%
%    Input, integer SN, the number of sample points.
%
%    Input, real S(M,NS), the points to be checked.
%
%    Output, integer NEAREST(NS), the index of the nearest cell generators.
%

%
%  Find nearest R to each S.
%
  nearest = zeros ( ns, 1 );

  for js = 1 : ns

    rs = sum ( ( r - repmat ( s(:,js), 1, nr ) ).^2 );

    [ dummy, nearest(js) ] = min ( rs );

  end

  return
end

author: johannesgerer

repo: jburkardt-m


text:
function x = p22_start ( n )

%*****************************************************************************80
%
%% P22_START returns a starting point for optimization for problem 22.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    30 December 2000
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, integer N, the number of variables X.
%
%    Output, real X(N), a starting point for the optimization.
%
  x = ( linspace ( - 5.12, + 5.12, n ) )';

  return
end


text:
function []=plot_reprojectVSreal_points_stereo(I1,I2,imagePoints,stereoParams,IndpairsUsed)
%% function for plot each image with reprojected points vs. real points and straight lines, and reprojection error statistics in STEP0
%
% INPUTS:
% * CBimagesInfo
% * cameraCBparameters:  a structure containing all the calibration parameters created in STEP0
%
%%

% plot each image with reprojected points vs. real points and straight lines
Nimages=length(IndpairsUsed);
reprojectedPoints1=stereoParams.CameraParameters1.ReprojectedPoints;
reprojectedPoints2=stereoParams.CameraParameters2.ReprojectedPoints;
reprojectionErrors1=stereoParams.CameraParameters1.ReprojectionErrors;
reprojectionErrors2=stereoParams.CameraParameters2.ReprojectionErrors;

if (stereoParams.CameraParameters1.ImageSize(1)>=stereoParams.CameraParameters1.ImageSize(2))
    n=1; m=8; p1=1:3; p2=4; p3=5:7; p4=8;
else
   n=2; m=6; p1=1:5; p2=6; p3=7:11; p4=12;
end

% PLOT
% define tabs
f = figure('name','Reprojected points on all images. Scroll between tabs to view the different images','units','normalized','outerposition',[.1 .1 .8 .8]);
tabgp = uitabgroup(f);

icount=0;
for iplot=IndpairsUsed'
    icount=icount+1;
    %plot
    tab(icount) = uitab(tabgp,'Title',['IM' num2str(iplot)]);

    axes('Parent',tab(icount)); % somewhere to plot

    subplot(n,m,p1)
    imshow(I1(:,:,:,iplot)); hold all;
    plot(imagePoints(:,1,icount,1), imagePoints(:,2,icount,1),'go','linewidth',1.5);
    plot(reprojectedPoints1(:,1,icount),reprojectedPoints1(:,2,icount),'r+','linewidth',1.5);
    title(['Camera left, Image ' num2str(iplot)]);
    drawnow
    legend('Detected Points','Reprojected Points');
    hold off;

    reprojectionErrorsNow1=reprojectionErrors1(:,:,icount);
    reprojectionErrorsMgnNow1=sqrt(sum(reprojectionErrorsNow1.^2,2));
    reprojectionErrorsNow1=[reprojectionErrorsNow1 reprojectionErrorsMgnNow1];

    subplot(n,m,p2)
    boxplot(reprojectionErrorsNow1,'Labels',{'X','Y','Mgn'});
    ylim([-max(abs(reprojectionErrorsNow1(:))) max(abs(reprojectionErrorsNow1(:)))]);
    title({'Reprojection error'; 'statistics [pix]'});

    subplot(n,m,p3)
    imshow(I2(:,:,:,iplot)); hold all;
    plot(imagePoints(:,1,icount,2), imagePoints(:,2,icount,2),'go','linewidth',1.5);
    plot(reprojectedPoints2(:,1,icount),reprojectedPoints2(:,2,icount),'r+','linewidth',1.5);
    title(['Camera right, Image ' num2str(iplot)]);
    drawnow
    legend('Detected Points','Reprojected Points');
    hold off;

    reprojectionErrorsNow2=reprojectionErrors2(:,:,icount);
    reprojectionErrorsMgnNow2=sqrt(sum(reprojectionErrorsNow2.^2,2));
    reprojectionErrorsNow2=[reprojectionErrorsNow2 reprojectionErrorsMgnNow2];

    subplot(n,m,p4)
    boxplot(reprojectionErrorsNow2,'Labels',{'X','Y','Mgn'});
    ylim([-max(abs(reprojectionErrorsNow2(:))) max(abs(reprojectionErrorsNow2(:)))]);
    title({'Reprojection error'; 'statistics [pix]'});

end

end

%%
% MultiDIC: a MATLAB Toolbox for Multi-View 3D Digital Image Correlation
%
% License: <https://github.com/MultiDIC/MultiDIC/blob/master/LICENSE.txt>
%
% Copyright (C) 2018  Dana Solav
%
% If you use the toolbox/function for your research, please cite our paper:
% <https://engrxiv.org/fv47e>

author: MultiDIC

repo: MultiDIC

save_path: github-repos/MATLAB/MultiDIC-MultiDIC


text:
% DIRECT JACOBIAN DEMO
function cinematica_inversa_3gdl
% definimos la posición y orientación del extremo en 2D
phi = pi/8
px = 2
py = 2
% definimos una matriz T
T = [cos(phi) -sin(phi) 0 px;
     sin(phi) cos(phi)  0 py;
     0 0 1 0;
     0 0 0 1];

 q = ikine3DOF(T)

 function q = ikine3DOF(T)
 L1 = 1;
 L2 = 1;
 L3 = 1; % 1m, parámetro del robot

 % hallamos pm
 p = T(1:3, 4);

%find angle Phi
x3 = T(1:3,1);

cphi = x3'*[1 0 0]';
sphi = x3'*[0 1 0]';
phi = atan2(sphi, cphi);

% pm
pm = p - L3*x3

%Distance of the point to the origin.
R= sqrt(pm(1)^2+pm(2)^2);

if R > (L1+L2)
   disp('\ninversekinematic_3dofplanar: unfeasible solution. The point cannot be reached');
end

%compute geometric solution
beta = atan2(pm(2),pm(1));
gamma = real(acos((L1^2+R^2-L2^2)/(2*R*L1)));
delta = real(acos((L1^2+L2^2-R^2)/(2*L1*L2)));

%arrange possible combinations for q(1) and q(2)
%elbow down     elbow up solutions
q =[beta+gamma beta-gamma;
    delta-pi   pi-delta];

%in this case, phi = q(1) + q(2) + q(3) and
%q(3) can be computed as q(3) = phi - q(1) - q(2)
%corresponding to each of the previous solutions, a unique q(3) can be
%computed for each case
for i=1:2 %iterate through columns
    q(3,i) = phi - q(1,i) - q(2,i);
end

repo: ARTE


text:
% IN ORDER TO SIMULATE THE PROGRAM:
%   A) FIRST, LOAD A ROBOT
%       robot = load_robot('abb','irb140');
%   B) NEXT, LOAD SOME EQUIPMENT.
%       robot.equipment = load_robot('equipment','tables/table_small');
%       OR
%       robot.equipment = load_robot('equipment','bumper_cutting');
%   C) NOW, LOAD AN END TOOL
%       robot.tool= load_robot('equipment','end_tools/parallel_gripper_0');
%   D) FINALLY, LOAD A PIECE TO GRAB BY THE ROBOT
%       robot.piece{1}=load_robot('equipment','cylinders/cylinder_tiny');
%
%   E) IF NECESSARY, CHANGE THE POSITION AND ORIENTATION OF THE ROBOT'S
%   BASE
%       robot.piece{1}.T0= [1 0 0 -0.35;
%                        0 1 0 -0.55;
%                        0 0 1 0.2;
%                        0 0 0 1];

function kuka_KR60_painting_jaguar_ftype

global TD_tool0 RT_tp1 RT_tp2 RT_tp3 RT_tp4 RT_tp5 RT_tp6 RT_tp7 RT_tp8 RT_tp9 RT_tp10 RT_tp11 RT_tp12 RT_tp13 RT_tp14 RT_tp15 RT_tp16

TD_tool0=[1,[[0,0,0],[1,0,0,0]],[0,[0,0,0],[1,0,0,0],0,0,0]];

RT_tp1=[[1.3400, -0.0000, 1.8100],[0.7006, -0.0954, 0.7006, -0.0954], [0, 0, 1, 0], [9E+09,9E+09,9E+09,9E+09,9E+09,9E+09]];
RT_tp2=[[1.1648, 0.6617, 1.8206],[0.7233, -0.2626, 0.6324, 0.0894], [-1, 0, 1, 1], [9E+09,9E+09,9E+09,9E+09,9E+09,9E+09]];
RT_tp3=[[1.4924, 0.9803, 1.1556],[0.5816, -0.3339, 0.7362, 0.0911], [-1, 0, 1, 1], [9E+09,9E+09,9E+09,9E+09,9E+09,9E+09]];
RT_tp4=[[1.5020, 0.8533, 1.3436],[0.6364, -0.2942, 0.7086, 0.0787], [-1, 0, 1, 1], [9E+09,9E+09,9E+09,9E+09,9E+09,9E+09]];
RT_tp5=[[1.5207, 0.6786, 1.4836],[0.6784, -0.2480, 0.6898, 0.0507], [-1, -1, 1, 1], [9E+09,9E+09,9E+09,9E+09,9E+09,9E+09]];
RT_tp6=[[1.6593, 0.3815, 1.4044],[0.6580, -0.1850, 0.7298, -0.0146], [-1, -1, 1, 1], [9E+09,9E+09,9E+09,9E+09,9E+09,9E+09]];
RT_tp7=[[1.6990, 0.1099, 1.4044],[0.6547, -0.1255, 0.7423, -0.0677], [-1, 0, 1, 1], [9E+09,9E+09,9E+09,9E+09,9E+09,9E+09]];
RT_tp8=[[1.6707, -0.3278, 1.4044],[0.6406, -0.0289, 0.7523, -0.1514], [0, -1, 1, 1], [9E+09,9E+09,9E+09,9E+09,9E+09,9E+09]];
RT_tp9=[[1.6064, -0.8202, 1.1348],[0.5463, 0.0815, 0.8060, -0.2127], [0, -1, 1, 1], [9E+09,9E+09,9E+09,9E+09,9E+09,9E+09]];
RT_tp10=[[1.4123, -1.1863, 0.9392],[0.5177, 0.1723, 0.7925, -0.2725], [0, -1, 1, 1], [9E+09,9E+09,9E+09,9E+09,9E+09,9E+09]];
RT_tp11=[[1.2529, -1.3639, 0.9736],[0.5680, 0.2052, 0.7172, -0.3476], [0, 0, 1, 1], [9E+09,9E+09,9E+09,9E+09,9E+09,9E+09]];
RT_tp12=[[1.2075, -0.7363, 1.1201],[0.5774, 0.1073, 0.7697, -0.2503], [0, -1, 1, 1], [9E+09,9E+09,9E+09,9E+09,9E+09,9E+09]];
RT_tp13=[[1.2048, -0.5831, 1.2640],[0.5698, 0.0714, 0.7900, -0.2149], [0, 0, 1, 1], [9E+09,9E+09,9E+09,9E+09,9E+09,9E+09]];
RT_tp14=[[1.2754, -0.1696, 1.3640],[0.6362, -0.0526, 0.7587, -0.1299], [0, -1, 1, 1], [9E+09,9E+09,9E+09,9E+09,9E+09,9E+09]];
RT_tp15=[[1.2481, 0.3127, 1.3640],[0.6493, -0.1940, 0.7354, -0.0081], [-1, 0, 1, 1], [9E+09,9E+09,9E+09,9E+09,9E+09,9E+09]];
RT_tp16=[[1.1268, 0.6211, 1.3640],[0.6449, -0.2871, 0.7042, 0.0756], [-1, -1, 1, 1], [9E+09,9E+09,9E+09,9E+09,9E+09,9E+09]];
main
end

function main()

global TD_tool0 RT_tp1 RT_tp2 RT_tp3 RT_tp4 RT_tp5 RT_tp6 RT_tp7 RT_tp8 RT_tp9 RT_tp10 RT_tp11 RT_tp12 RT_tp13 RT_tp14 RT_tp15 RT_tp16

 MoveL(RT_tp1, 'vmax' , 'fine' , TD_tool0, 'wobj0');
 MoveL(RT_tp2, 'vmax' , 'fine' , TD_tool0, 'wobj0');
 MoveL(RT_tp3, 'vmax' , 'fine' , TD_tool0, 'wobj0');
 MoveL(RT_tp4, 'vmax' , 'fine' , TD_tool0, 'wobj0');
 MoveL(RT_tp5, 'vmax' , 'fine' , TD_tool0, 'wobj0');
 MoveL(RT_tp6, 'vmax' , 'fine' , TD_tool0, 'wobj0');
 MoveL(RT_tp7, 'vmax' , 'fine' , TD_tool0, 'wobj0');
 MoveL(RT_tp8, 'vmax' , 'fine' , TD_tool0, 'wobj0');
 MoveL(RT_tp9, 'vmax' , 'fine' , TD_tool0, 'wobj0');
 MoveL(RT_tp10, 'vmax' , 'fine' , TD_tool0, 'wobj0');
 MoveL(RT_tp11, 'vmax' , 'fine' , TD_tool0, 'wobj0');
 MoveL(RT_tp12, 'vmax' , 'fine' , TD_tool0, 'wobj0');
 MoveL(RT_tp13, 'vmax' , 'fine' , TD_tool0, 'wobj0');
 MoveL(RT_tp14, 'vmax' , 'fine' , TD_tool0, 'wobj0');
 MoveL(RT_tp15, 'vmax' , 'fine' , TD_tool0, 'wobj0');
 MoveL(RT_tp16, 'vmax' , 'fine' , TD_tool0, 'wobj0');

end


text:
function triangle_fekete_rule_test01 ( )

%*****************************************************************************80
%
%% TRIANGLE_FEKETE_RULE_TEST01 tests FEKETE_RULE_NUM, FEKETE_DEGREE, FEKETE_ORDER_NUM.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    16 March 2014
%
%  Author:
%
%    John Burkardt
%
  fprintf ( 1, '\n' );
  fprintf ( 1, 'TRIANGLE_FEKETE_RULE_TEST01\n' );
  fprintf ( 1, '  FEKETE_RULE_NUM returns the number of rules;\n' );
  fprintf ( 1, '  FEKETE_DEGREE returns the degree of a rule;\n' );
  fprintf ( 1, '  FEKETE_ORDER_NUM returns the order of a rule.\n' );

  rule_num = fekete_rule_num ( );

  fprintf ( 1, '\n' );
  fprintf ( 1, '  Number of available rules = %d\n', rule_num );
  fprintf ( 1, '\n' );
  fprintf ( 1, '      Rule    Degree     Order\n' );
  fprintf ( 1, '\n' );

  for rule = 1 : rule_num
    order_num = fekete_order_num ( rule );
    degree = fekete_degree ( rule );
    fprintf ( 1, '  %8d  %8d  %8d\n', rule, degree, order_num );
  end

  return
end


text:
function conv_feat_data = extract_image_activation_maps(CNN, image, scales, mean_pix, ...
    Semantic_Aware_CNN, semantic_scales)
% extract_image_activation_maps(CNN, image, scales, mean_pix): extract the
% activation maps of one image (section 3 of technical report) for the
% specified scales using the convolutional neural network CNN.
% extract_image_activation_maps(CNN, image, scales, mean_pix, Semantic_Aware_CNN, semantic_scales):
% In this case the current function it also extracts the semantic
% segmentation aware activation maps (see section 4 of the technical
% report).
%
% INPUTS:
% 1) CNN: the caffe net struct with the convolutional neural network that
% implements the activation mas module (section 3)
% 2) image: a Height x Width x 3 uint8 array that represents the image
% pixels
% 3) scales: NumScales x 1 or 1 x NumScales vector with the images scales
% that will be used. The i-th value should be the size in pixels of the
% smallest dimension of the image in the i-th scale.
% 4) mean_pix: is a 3 x 1 or 1 x 3 vector with the mean pixel value per
% color channel that is subtracted from the scaled image before is being
% fed to the CNN
% 5) Semantic_Aware_CNN (OPTIONAL): the caffe net struct with the
% convolutional neural that implements the activation mas module for the
% semantic segmentation aware CNN features (section 4). The Semantic_Aware_CNN
% network gets as input the convolutional feature maps that the CNN network
% yields and outputs semantic segmentation aware activation maps.
% 6) semantic_scales: a NumScales2 x 1 or 1 x NumScales2 vector with the
% images scales that will be used for the semantic segmentation aware
% features. The elements of this vector should be a subset of the scales vector.
%
% OUTPUTS:
% 1) conv_feat_data: a struct that includes the activation maps of the
% image. Its field is:
%    conv_feat_data.feat:
%    1.a) In case the function is called with the arguments:
%    extract_image_activation_maps(CNN, image, scales, mean_pix)
%    then conv_feat_data.feat is a struct that includes 1) the convolutional
%    feature maps (field rsp) that the CNN network yields, 2) the image scales
%    from which they were extracted (field scale), and 3) the original size
%    of the image (fields im_height and im_width)
%    1.b) In case the function is called with the arguments:
%    extract_image_activation_maps(CNN, image, scales, mean_pix, Semantic_Aware_CNN, semantic_scales):
%    then it is a 1 x 2 cell array where 1st element is a struct with the
%    convolutional feature maps of the CNN network (like in the 1.a case)
%    and the 2nd element is a struct with the convolutional feature maps of
%    the Semantic_Aware_CNN network (like in the 1.a case).
%
% This file is part of the code that implements the following ICCV2015 accepted paper:
% title: "Object detection via a multi-region & semantic segmentation-aware CNN model"
% authors: Spyros Gidaris, Nikos Komodakis
% institution: Universite Paris Est, Ecole des Ponts ParisTech
% Technical report: http://arxiv.org/abs/1505.01749
% code: https://github.com/gidariss/mrcnn-object-detection
%
%
% AUTORIGHTS
% --------------------------------------------------------
% Copyright (c) 2015 Spyros Gidaris
%
% "Object detection via a multi-region & semantic segmentation-aware CNN model"
% Technical report: http://arxiv.org/abs/1505.01749
% Licensed under The MIT License [see LICENSE for details]
% ---------------------------------------------------------

conv_feat_data = init_feat_data();
conv_feat_data.feat.im_height = size(image,1);
conv_feat_data.feat.im_width  = size(image,2);
% extract the activation maps of an image for a given set of scales
[conv_feat_data.feat.rsp, conv_feat_data.feat.scale] = extract_conv_features(...
    CNN, image, scales, mean_pix);

if exist('Semantic_Aware_CNN','var')>0
    assert(exist('semantic_scales','var')>0)
    % extract the semantic segmentation aware activation maps of an image
    % of a given set of scales and the convolutional feature maps (activation
    % maps) that were previously extracted from the image using the CNN
    % network
    semantic_conv_feat_data      = conv_feat_data;
    semantic_conv_feat_data.feat = pick_scales_if_there(...
        semantic_conv_feat_data.feat, semantic_scales);
    conf.do_interleave        = true; % if set to true then the
    % resolution augmentation technique described on the OverFeat paper: http://arxiv.org/abs/1312.6229
    % (section 3.3 of OverFeat technical report) is being used
    conf.interleave_num_steps = 2; % a scalar value for the number of steps
    % that are being using on the above resolution augmentation technique

    % extract the semantic segmentation aware activation maps
    semantic_conv_feat_data.feat.rsp = extract_semantic_seg_features_from_conv5(...
        Semantic_Aware_CNN, semantic_conv_feat_data.feat.rsp, conf);
    conv_feat_data.feat = {conv_feat_data.feat, semantic_conv_feat_data.feat};
end

end

function d = init_feat_data()
d.feat     = [];
end

function feat = pick_scales_if_there(feat, scales)
num_scales = length(scales);

found_scales = zeros(size(scales));
found_rsp = {};
c = 0;
for s = 1:num_scales
    scale_index = find(feat.scale == scales(s));
    if ~isempty(scale_index)
        assert(numel(scale_index) == 1);
        c = c + 1;
        found_scales(c) = scales(s);
        found_rsp{c}    = feat.rsp{scale_index};
    end
end
feat.scale = found_scales(1:c);
feat.rsp   = found_rsp;

end

author: gidariss

repo: mrcnn-object-detection


text:
% Converts a MatLab array into a java.util.ArrayList.
%
% Input arguments:
% array:
%    a MatLab row or column vector (with elements of any type)
%
% Output arguments:
% list:
%    a java.util.ArrayList instance
%
% See also: javaArray

% Copyright 2010 Levente Hunyadi
function list = javaArrayList(array)
list = java.util.ArrayList;
if ~isempty(array)
    assert(isvector(array), 'javaArrayList:DimensionMismatch', ...
        'Row or column vector expected.');
    if iscell(array)  % convert cell array into ArrayList
        for k = 1 : numel(array)
            list.add(array{k});
        end
    else  % convert (numeric) array into ArrayList
        for k = 1 : numel(array)
            list.add(array(k));
        end
    end
end

repo: REST


text:
function printStack(ME)
global logger %#ok<NUSED>

if ~exist('ME','var')
    ME = [];
end

stdoutFuncName = getStandardOutputFuncName();

if isempty(ME)
    s = dbstack;
else
    s = ME.stack;
end

eval( sprintf('%s(''----------------------------------------------\\n'')', stdoutFuncName) );
if ~isempty(ME)
    eval( sprintf('%s(''ERROR:    %%s\\n'', ME.message)', stdoutFuncName) );
    eval( sprintf('%s(''Current Folder :  %%s\\n'', filesepStandard(pwd))', stdoutFuncName) );
end
eval( sprintf('%s(''Call stack:\\n'')', stdoutFuncName) );
for ii = 1:length(s)
    [~,f,e] = fileparts(s(ii).file); %#ok<*ASGLU>
    eval( sprintf('%s(''    Error in %%s > %%s (line %%d)\\n'', [f, e], s(ii).name, s(ii).line)', stdoutFuncName) );
end
eval( sprintf('%s(''----------------------------------------------\\n\\n'')', stdoutFuncName) );

author: BUNPC


text:
function a = subsref(t,s)
%SUBSREF Subscripted reference for a sparse tensor.
%
%   We can extract elements or subtensors from a sparse tensor in the
%   following ways.
%
%   Case 1a: y = X(i1,i2,...,iN), where each in is an index, returns a
%   scalar.
%
%   Case 1b: Y = X(R1,R2,...,RN), where one or more Rn is a range and
%   the rest are indices, returns a sparse tensor. The elements are
%   renumbered here as appropriate.
%
%   Case 2a: V = X(S) or V = X(S,'extract'), where S is a p x n array
%   of subscripts, returns a vector of p values.
%
%   Case 2b: V = X(I) or V = X(I,'extract'), where I is a set of p
%   linear indices, returns a vector of p values.
%
%   Any ambiguity results in executing the first valid case. This
%   is particularily an issue if ndims(X)==1.
%
%   S = X.subs returns the subscripts of the nonzero entries in X.
%
%   V = X.vals returns the values of the nonzero entries in X.
%
%   Examples
%   X = sptensor([4,4,4;2,2,1;2,3,2],[3;5;1],[4 4 4]);
%   X(1,2,1) %<-- returns zero
%   X(4,4,4) %<-- returns 3
%   X(3:4,:,:) %<-- returns 2 x 4 x 4 sptensor
%   X(2,:,:) %<-- returns a 2 x 2 tensor
%   X([1,1,1;2,2,1]) %<-- returns a vector of 2 elements
%   X = sptensor([6;16;26],[1;1;1],30);
%   X([1:6]') %<-- extracts a subtensor
%   X([1:6]','extract') %<-- extracts a vector of 6 elements
%
%   See also SPTENSOR, SPTENSOR/FIND, TENSOR/SUBSREF.
%
%MATLAB Tensor Toolbox.
%Copyright 2012, Sandia Corporation.

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

switch s(1).type

    case '{}'
        error('Subscript cell reference cannot be used with sptensor.')

    case '.'
        switch s(1).subs
            case {'subs','indices'}
                a = tt_subsubsref(t.subs, s);
            case {'vals','values'}
                a = tt_subsubsref(t.vals, s);
            case 'size'
                a = tt_subsubsref(t.size, s);
            otherwise
                error(['No such field: ', s(1).subs]);
        end
        return;

    case '()'

        % *** CASE 1: Rectangular Subtensor ***
        if (numel(s(1).subs) == ndims(t)) && ...
                ~isequal(s(1).subs{end},'extract')

            % Extract the subdimensions to be extracted from t
            region = s(1).subs;

            % Pare down the list of subscripts (and values) to only
            % those within the subdimensions specified by region.
            loc = subdims(region, t);
            subs = t.subs(loc,:);
            vals = t.vals(loc);

            % Find the size of the subtensor and renumber the
            % subscripts
            [subs, sz] = renumber(subs, size(t), region);

            % Determine the subscripts
            newsiz = [];                % (future) new size
            kpdims = [];                % dimensions to keep
            rmdims = [];                % dimensions to remove

            % Determine the new size and what dimensions to keep
            for i = 1:length(region)
                if ischar(region{i}) && (region{i} == ':')
                    newsiz = [newsiz size(t,i)];
                    kpdims = [kpdims i];
                elseif numel(region{i}) > 1
                    newsiz = [newsiz numel(region{i})];
                    kpdims = [kpdims i];
                else
                    rmdims = [rmdims i];
                end
            end

            % Return a single double value for a zero-order sub-tensor
            if isempty(newsiz)
                if isempty(vals)
                    a = 0;
                else
                    a = vals;
                end
                return;
            end

            % Assemble the resulting sparse tensor
            if isempty(subs)
                a = sptensor([],[],sz(kpdims));
            else
                a = sptensor(subs(:,kpdims), vals, sz(kpdims));
            end
            return;
        end

        % Case 2: EXTRACT

        % *** CASE 2a: Subscript indexing ***
        if size(s(1).subs{1},2) == ndims(t)

            % extract array of subscripts
            srchsubs = s(1).subs{1};

       % *** CASE 2b: Linear indexing ***
        else

            % Error checking
            if numel(s(1).subs) ~= 1
                error('Invalid indexing');
            end

            idx = s(1).subs{1};
            if ndims(idx) ~=2 || size(idx,2) ~= 1
                error('Expecting a column index');
            end

            % extract linear indices and convert to subscripts
            srchsubs = tt_ind2sub(size(t),idx);

        end

        a = extract(t,srchsubs);
        a = tt_subsubsref(a,s);

        return;

    otherwise
        error('Incorrect indexing into sptensor.')
end

author: andrewssobral

repo: mtt

save_path: github-repos/MATLAB/andrewssobral-mtt


text:
function m = sqdist(p, q, A)
% SQDIST      Squared Euclidean or Mahalanobis distance.
% SQDIST(p,q)   returns m(i,j) = (p(:,i) - q(:,j))'*(p(:,i) - q(:,j)).
% SQDIST(p,q,A) returns m(i,j) = (p(:,i) - q(:,j))'*A*(p(:,i) - q(:,j)).
% The Lightspeed Matlab toolbox
% Written by Tom Minka

[d, pn] = size(p);
[d, qn] = size(q);

if pn == 0 || qn == 0
  m = zeros(pn,qn);
  return
end

if nargin == 2

  pmag = col_sum(p .* p);
  qmag = col_sum(q .* q);
  m = repmat(qmag, pn, 1) + repmat(pmag', 1, qn) - 2*p'*q;
  %m = ones(pn,1)*qmag + pmag'*ones(1,qn) - 2*p'*q;

else

  Ap = A*p;
  Aq = A*q;
  pmag = col_sum(p .* Ap);
  qmag = col_sum(q .* Aq);
  m = repmat(qmag, pn, 1) + repmat(pmag', 1, qn) - 2*p'*Aq;

end

author: layumi

repo: Image-Text-Embedding


text:
function Q = planarize(V,F,varargin)
  % PLANARIZE Planarize the quads of quad mesh using the local-local technique
  % of "Interactive Planarization and Optimization of 3D Meshes" [Poranne et
  % al. 2012]
  %
  % Q = planarize(V,F,varargin)
  %
  % Inputs:
  %   V  #V by dim list of mesh vertex positions
  %   F  #F by 4 list of quad indices into V
  %     Optional:
  %       'MaxIter' followed by maximum number of iterations {100}
  %       'mu' followed by "mu" penalty value {0.9}
  %       'r' followed by "r" penalty shrink factor {0.9}
  % Outputs:
  %   Q  #V by dim list of new mesh vertex positions
  %

  % default values
  max_iter = 100;
  mu = 0.9;
  r = 0.9;
  % Map of parameter names to variable names
  params_to_variables = containers.Map( ...
    {'MaxIter','mu','r'}, ...
    {'max_iter','mu','r'});
  v = 1;
  while v <= numel(varargin)
    param_name = varargin{v};
    if isKey(params_to_variables,param_name)
      assert(v+1<=numel(varargin));
      v = v+1;
      % Trick: use feval on anonymous function to use assignin to this workspace
      feval(@()assignin('caller',params_to_variables(param_name),varargin{v}));
    else
      error('Unsupported parameter: %s',varargin{v});
    end
    v=v+1;
  end

  % number of quads
  m = size(F,1);
  % number of vertices
  n = size(V,1);
  dim = 3;
  T = [F(:,[1 2 3]);F(:,[1 3 4])];
  N = normalizerow(normals(V,T));
  % Initial guess for N
  N = 0.5*(N(1:m,:)+N(m+1:end,:));
  A = doublearea(V,T);
  % A sort of barycentric massmatrix for quad meshes
  M = repdiag(sparse(T(:),T(:),0.25*repmat(A,3,1),n,n),dim);
  E = [];

  % nsp = 3;
  % subplot(1,nsp,1);
  % BC = barycenter(V,F);
  % trisurf(F,V(:,1),V(:,2),V(:,3));
  % hold on;
  % quiver3(BC(:,1),BC(:,2),BC(:,3),N(:,1),N(:,2),N(:,3));
  % hold off;

  % "corner" incidence
  C = sparse(F(:),repmat(1:m,1,4)',1,n,m);

  D = zeros(m,1);

  Q = V;
  iter = 1;
  while true

    % LOCAL
    for ni = 1:m
      Qn = Q(F(ni,:),:);
      Mn = mean(Qn,1);
      Qnc = bsxfun(@minus,Qn,Mn);
      [EV,ED] = eig(Qnc'*Qnc);
      [~,mi] = min(diag(ED));
      N(ni,:) = EV(:,mi);
      D(ni) = -Mn*N(ni,:)';
    end

    % LOCAL
    for v = 1:n
      % faces incident on v
      J = find(C(v,:));
      NJ = N(J,:);
      DJ = D(J);
      AJ = mu*eye(3) + (1-mu)*(NJ'*NJ);
      bJ = mu*V(v,:)' - (1-mu)*NJ'*DJ;
      Q(v,:) = AJ\bJ;
    end

    % GLOBAL
    %% "corner" incidence for each dimension with Nj
    %This is wrong becase D should be scalar...
    %CN = sparse( ...
    %  reshape(bsxfun(@plus,F(:),(0:dim-1)*n),dim*m*4,1), ...
    %  reshape(bsxfun(@plus,repmat((1:m)',4,1),(0:dim-1)*m),dim*m*4,1), ...
    %  reshape(repmat(N,4,1),4*m*dim,1), ...
    %  n*dim,m);
    %% μ ‖P-Q‖² + (1-μ)∑_ij(n_j*q_i-d_j)²
    %% Identity
    %I = repdiag(speye(m,m),dim);
    %H = [          mu*M        (1-mu)*CN; ...
    %          (1-mu)*CN'       (1-mu)*I];
    %l = [-mu*2*M*V(:);zeros(m*dim,1)];
    %QD = min_quad_with_fixed(H,l,[],[]);
    %Q = reshape(QD(1:n*dim),n,dim);
    %D = reshape(QD(n*dim+1:end),m,dim);
    %% (1-μ)∑_ij(n_j*q_i-d_j)²
    %% ‖n_i‖² = 1
    %%
    %CQ = sparse( ...
    %  reshape(bsxfun(@plus,F(:),(0:dim-1)*n),dim*m*4,1), ...
    %  reshape(bsxfun(@plus,repmat((1:m)',4,1),(0:dim-1)*m),dim*m*4,1), ...
    %  reshape(Q(F(:),:),m*dim*4,1), ...
    %  n*dim,m*dim);
    %%H = (1-mu)*diag(sum(CQ'.^2,2));
    %%l = (1-mu)*-2*diag(sum(CQ',2))*D(:);
    %%% min N'HN - 2N'l
    %%% ‖N‖² = m

    %subplot(1,nsp,2);
    %BC = barycenter(Q,F);
    %trisurf(F,Q(:,1),Q(:,2),Q(:,3));
    %hold on;
    %quiver3(BC(:,1),BC(:,2),BC(:,3),N(:,1),N(:,2),N(:,3));
    %hold off;

    %NT = normalizerow(normals(Q,T));
    %e = sum(sum((NT(1:m,:)-NT(m+(1:m),:)).^2,2));
    %E = [E(:);e];
    %subplot(1,nsp,3);
    %loglog(E);
    iter = iter+1;
    if iter>max_iter
      break;
    end

    mu = r*mu;
  end
end

author: alecjacobson

repo: gptoolbox


text:
function slice_regressors = tapas_physio_split_regressor_slices(...
    slice_regressors_concat, nSampleSlices)
% Splits a concatenated phase regressor into slice-wise regressors, one per
% column
%
%   slice_regressors = tapas_physio_split_regressor_slices(...
%       slice_regressors_concat, nSampleSlices)
%
% IN
%   slice_regressors_concat [nSamples*nSlices, nRegressors] - matrix which
%                           per column has entries of the form:
%                           regressor value volume 1, slice 1
%                                           volume 1, slice 2
%                                           ...
%                                           volume 1, slice nSampleSlices
%                                           volume 2, slice 1
%                                           volume 2, slice 2
%
% OUT
%
% EXAMPLE
%   tapas_physio_split_regressor_slices
%
%   See also

% Author: Lars Kasper
% Created: 2014-08-14
% Copyright (C) 2014 TNU, Institute for Biomedical Engineering, University of Zurich and ETH Zurich.
%
% This file is part of the TAPAS PhysIO Toolbox, which is released under the terms of the GNU General Public
% Licence (GPL), version 3. You can redistribute it and/or modify it under the terms of the GPL
% (either version 3 or, at your option, any later version). For further details, see the file
% COPYING or <http://www.gnu.org/licenses/>.

% make multiple columns for multiple slices
nRegressors = size(slice_regressors_concat,2);
slice_regressors = [];
for iRegressor = 1:nRegressors
    slice_regressors = [slice_regressors, ...
        reshape(reshape(slice_regressors_concat(:,iRegressor),1,[]), ...
        [], nSampleSlices)];
end

author: translationalneuromodeling

repo: tapas


text:
% SIMPLE ALGORITHM TO FOLLOW A LINE IN SPACE.
% Copyright (C) 2019, by Arturo Gil Aparicio
%
% This file is part of ARTE (A Robotics Toolbox for Education).
%
% ARTE is free software: you can redistribute it and/or modify
% it under the terms of the GNU Lesser General Public License as published by
% the Free Software Foundation, either version 3 of the License, or
% (at your option) any later version.
%
% ARTE is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
% GNU Lesser General Public License for more details.
%
% You should have received a copy of the GNU Leser General Public License
% along with ARTE.  If not, see <http://www.gnu.org/licenses/>.
function speed_demo_line
close all
% velocidad lineal entre puntos consecutivos
abs_linear_speed = 0.5; % (m/s)
delta_time = 0.05;

robot = load_robot('ABB', 'IRB52')

fprintf('\nPRESS ANY KEY TO CONTINUE...')
pause

%NOA matrix initial point
T1=[1 0 0 0.8;
    0 1 0 -0.3;
    0 0 1 0.9;
    0 0 0  1]
%NOA matrix end point
T2=[1 0 0 0.5;
    0 1 0 0.4;
    0 0 1 0.7;
    0 0 0  1]

punto_inicial = T1(1:3,4);
punto_final = T2(1:3,4);
% vector velocidad en la direcci�n de la trayectoria
v = (punto_final-punto_inicial);
v = abs_linear_speed*v/norm(v); %vector normalizado en la direcci�n de la recta
w = [0 0 0]';
xd = [v; w];

qinv = inversekinematic(robot, T1);

%Select arbitrarily the first solution
q = qinv(:,1);
qt = [q];
for i=1:1000
   J = manipulator_jacobian(robot, q);
   qd = inv(J)*xd;
   q = q + delta_time*qd;
   qt = [qt q];
end
animate(robot, qt)


text:
classdef TVPsettings < bear.settings.BASEsettings
    %TVPSETTINGS Panel VAR settings class
    %   The bear.settings.TVPsettings class is a class that creates a settings
    %   object to run a time-varying VAR. It can be created directly by
    %   running:
    %
    %   bear.settings.TVPsettings(ExcelFile, varargin)
    %
    %   or in its more convenient form:
    %
    %   BEARsettings('TVP', ExcelFile = 'path/To/file.xlsx')
    %
    % TVPsettings Properties:
    %    tvbvar          - choice of time-varying BVAR model
    %    It              - Gibbs sampler iterations
    %    Bu              - Gibbs sampler burn-in iterations
    %    pick            - retain only one post burn iteration
    %    pickf           - frequency of iteration picking
    %    ar              - auto-regressive coefficients
    %    PriorExcel      - Select individual priors
    %    priorsexogenous - Gibbs sampler burn-in iterations
    %    lambda4         - hyperparameter
    %    gamma           - hyperparameter
    %    alpha0          - hyperparameter
    %    delta0          - hyperparameter
    %    alltirf         - calculate IRFs for every sample period
    %    favar           - FAVAR Options

    properties
        % choice of time-varying BVAR model
        % 1 = time-varying coefficients (TVP)
        % 2 = general time-varying (TVP_SV)
        tvbvar (1,1) bear.TVPtype = 2;
        % switch to Excel interface
        PriorExcel (1,1) logical = false; % set to 1 if you want individual priors, 0 for default
        % switch to Excel interface for exogenous variables
        priorsexogenous (1,1) logical = true; % set to 1 if you want individual priors, 0 for default
        % hyperparameter: lambda4
        lambda4=100;
    end

    properties %Hyperparameters
        % hyperparameter: gamma
        gamma (1,1) double = 0.85;
        % hyperparameter: alpha0
        alpha0 (1,1) double = 0.001;
        % hyperparameter: delta0
        delta0 (1,1) double = 0.001;
        % calculate IRFs for every sample period (1=yes, 0=no)
        alltirf (1,1) logical = true;
    end

    properties
        % total number of iterations for the Gibbs sampler
        It (1,1) double {mustBeGreaterThanOrEqual(It,1)} = 2000;
        % number of burn-in iterations for the Gibbs sampler
        Bu (1,1) double = 1000;
        % choice of retaining only one post burn iteration over 'pickf' iterations (1=yes, 0=no)
        pick (1,1) logical = false;
        % frequency of iteration picking (e.g. pickf=20 implies that only 1 out of 20 iterations will be retained)
        pickf=20;
        % just for the code to run (do not touch)
        ar=0;
    end

    properties % FAVAR
        % FAVAR options
        favar (1,1) bear.settings.favar.FAVARsettings = bear.settings.favar.VARtypeSpecificFAVARsettings; % augment VAR model with factors (1=yes, 0=no)
    end

    methods

        function obj = TVPsettings(excelPath, varargin)

            obj@bear.settings.BASEsettings(6, excelPath)

            obj = parseBEARSettings(obj, varargin{:});

        end

        function obj = set.Bu(obj,value)
            if (value <= obj.It-1) %#ok<MCSUP>
                obj.Bu = value;
            else
                error('bear:settings:TVPsettings',"The maximum value of Bu is It-1: " + (obj.It-1)) %#ok<MCSUP>
            end
        end

        function obj = set.It(obj,value)
            if (value > obj.Bu-1) %#ok<MCSUP>
                obj.It = value;
            else
                error('bear:settings:TVPsettings',"The minimum value of It is Bu+1: " + (obj.Bu+1)) %#ok<MCSUP>
            end
        end

    end

end

author: european-central-bank

repo: BEAR-toolbox


text:
function [model, newHeight, newWidth] = vargplvmReduceVidModel(model, height, width, varargin)

% VARGPLVMREDUCEVIDMODEL Take a model computed on a video dataset and return a model which is computed on
% the same video but with lower resolution
% DESC Take a model computed on a video dataset and return a model which is computed on
% the same video but with lower resolution.
% SEEALSO: vargplvmReduceVideo.m
% VARGPLVM

% varargin may include "factor1" and "factor2" optionally

[model.m, newHeight, newWidth] = vargplvmReduceVideo(model.m, height, width, varargin{:});
[model.scale, newHeight, newWidth] = vargplvmReduceVideo(model.scale, height, width, varargin{:});
[model.bias, newHeight, newWidth] = vargplvmReduceVideo(model.bias, height, width, varargin{:});
model.d = newHeight * newWidth;

if ~isempty(model.y)
    [model.y, newHeight, newWidth] = vargplvmReduceVideo(model.y, height, width, varargin{:});
end

model.TrYY = sum(sum(model.m .* model.m));
model.P = model.P1 * (model.Psi1' * model.m);
model.B = model.P1' * model.P;

params = vargplvmExtractParam(model);
model = vargplvmExpandParam(model, params);

author: SheffieldML

repo: vargplvm


text:
function x= BesseliRatio(nu,kappa,maxIter)
%%BESSELIRATIO Evaluate the ratio of modified Bessel functions of the first
%              kind of the form x=I_{nu}(kappa)/I_{nu-1}(kappa).
%
%INPUTS: nu The positive integer (upper) order of the modified Bessel
%           function of the first kind in the ratio; nu>=1.
%     kappa The real argument of the modified Bessel function of the first
%           kind; kappa>=0.
%   maxIter An optional parameter specifying the maximum number of
%           iterations to use for computing the ratio using an iterative
%           method. If this parameter is omitted or an empty matrix is
%           passed, the default value of 2000 is used. Convergence usually
%           occurs long before the maximum number of iterations is reached.
%
%OUTPUTS: x The value of the ratio I_{nu}(kappa)/I_{nu-1}(kappa).
%
%Numerical precision limitations can make the evaluation of Bessel function
%radios difficult if one tries to explicitly evaluate the functions. Here,
%the algorithm of Perron described in [1] is implemented.
%
%REFERENCES:
%[1] W. Gautschi and J. Slavik, "On the computation of modified Bessel
%    function ratios," Mathematics of Computation, vol. 32, no. 143, pp.
%    865-875, Jul. 1978.
%
%September 2013 David F. Crouse, Naval Research Laboratory, Washington D.C.
%(UNCLASSIFIED) DISTRIBUTION STATEMENT A. Approved for public release.

    if(nargin<3||isempty(maxIter))
        %It usually converges long before the maximum number of iterations.
        maxIterations=2000;
    end
    %tolerance for convergence. This should be suitable for double-
    %precision computations.
    tol=1e-15;

    k=1;
    cumProd=1;
    pCur=0.5*kappa*(nu+0.5)/((nu+kappa/2)*(nu+kappa+0.5)-0.5*kappa*(nu+0.5));
    cumProd=cumProd*pCur;
    cumSum=1+cumProd;
    while(k<maxIterations)
        pPrev=pCur;
        k=k+1;
        pCur=0.5*kappa*(nu+k-0.5)*(1+pPrev)/((nu+kappa+(k-1)/2)*(nu+kappa+k/2)-0.5*kappa*(nu+k-0.5)*(1+pPrev));
        cumProd=cumProd*pCur;

        cumSumNew=cumSum+cumProd;
        %If we have hit a precision limit; this usually happens before k
        %gets too big.
        if(abs(cumSumNew-cumSum)<tol)
            break;
        end
        cumSum=cumSumNew;
    end
    a0=(kappa+2*nu)/kappa;
    x=cumSum/a0;
end

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

author: USNavalResearchLaboratory

repo: TrackerComponentLibrary


text:
function Select = subsetSelection(Obj,Objhat,N)

%------------------------------- 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".
%--------------------------------------------------------------------------

% This function is written by Tomoaki Takagi

    Len = length(Obj);
    Obj = [Obj;Objhat];

    %% Select the representative objective vector set
    LpNormD = pdist2(Obj,Obj);
    Select = false(1,size(Obj,1));
    Select(1:Len) = true;
    % Greedy inclusion distance-based subset slection
    while sum(Select) < N
        Remain   = find(~Select);
        [~, rho] = max(min(LpNormD(Remain,Select),[],2));
        Select(Remain(rho)) = true;
    end
    Select = Select(Len+1:end);
end

author: BIMK

repo: PlatEMO

save_path: github-repos/MATLAB/BIMK-PlatEMO


text:
function [N] = lbf2N(lbf)
% Convert force from pounds to newtons.
% Chad A. Greene 2012
N = lbf*4.4482216152605;

author: Sable

repo: mcbench-benchmarks


text:
%% Logistic logistic regression for multiclass classification
clear
k = 3;
n = 1000;
[X,t] = kmeansRnd(2,k,n);
[model, llh] = logitMn(X,t);
y = logitMnPred(model,X);
plotClass(X,y)

author: PRML

repo: PRMLT

save_path: github-repos/MATLAB/PRML-PRMLT


text:
function [mu, muw]= ViscosityCO2Water(T, p, x_CO2)
% This function calculates the viscosity of the pure water, brine, and
% water-CO2, and brine-CO2 systems
% p is in Pa, T in K, x_CO2 is the mole fraction of CO2, mu in Pa.s
% Written by Ali A. Eftekhari at TU Delft
% See the license file
Mwater = 0.0180153; % kg/mol
% MCO2 = 0.04401; % kg/mol

P = p/1e6; % convert Pa to MPa
m_CO2 = DuanSun(T, p, 0); % CO2 equilibrium molality (mol CO2/kg water)
x_CO2_eq = m_CO2/(m_CO2+1/Mwater);
muw = IAPWS_IF97('mu_pT', P, T); % pure water viscosity in Pa.s
vis_ratio = 1+(-4.069e-3*(T-273.15)+0.2531)*x_CO2/x_CO2_eq;
mu = muw*vis_ratio;

author: simulkade

repo: FVTool

save_path: github-repos/MATLAB/simulkade-FVTool


text:
function [featureVector,selectedheaderskeys] = featurize_bigram(inputcellarray, nminFeatures, removeStopWords, doStem)
%
%
% featureVector = featurize_bigram(inputcellarray, nminFeatures, removeStopWords, doStem)
%
% Takes an input cell array in which each cell is a review or text and
% outputs a feature vector with a number of features (nminFeatures is the
% number of times that a feature should be presented in the corpus to be
% included in the feature vector.)
%
% The approach is basically a bag-of-word but we also add bigrams to the
% feature vector too
%
% removeStopWords is a flag, if it is true it will remove the stop words.
%
% Inputs:
%       inputcellarray: a cell array with texts as the content of each cell
%       nFeatures: the number of features that we like to see in the vetor
%       removeStopWords: if ==1 it will remove all the stop words
%       doStem: a flag if true porter stemmer will be used
% Outputs:
%       featureVector

% Test case:

% inputcellarray = {' MATLAB desktop keyboard shortcuts, such as Ctrl+S,  are now customizable.';' To customize keyboard shortcuts, use Preferences. From there, you can also  restore previous default settings by following the steps outlined in Help.'}
% nminFeatures = 1;
% removeStopWords = 0;
% doStem =1;
% featureVector = featurize_bigram(inputcellarray, nminFeatures, removeStopWords, doStem)

%
fid = fopen('english.stop');

stopwords = textscan(fid, '%s');
stopwords = stopwords{1,1};
fclose(fid);

% params
% n:minimum appearances of a stem
n=nminFeatures;

g = containers.Map();
bigs = containers.Map();
lastword = '';
for i = 1:size(inputcellarray,1)
    fprintf('%d/%d\n', i, size(inputcellarray,1));
    comment = inputcellarray{i};
    comment = SanitizeComment(comment);
    comment = lower(comment);
    r=regexp(comment,' ','split');
    for j =1:size(r,2)

        if doStem
            word = porterStemmer(cell2mat(r(j)));
        else
            word = (cell2mat(r(j)));
        end

        if removeStopWords % if the function caller wants to remove stopwords
            tfflag = ~isStopWord(word, stopwords);
        else
            tfflag = 1;
        end

        if isKey(g, word) && tfflag
            g(word) = g(word)+1;
        elseif tfflag & (~strcmp(word,' ')) & (~strcmp(word,''))
            g(word) = 1;
        end
        if  (~strcmp(lastword,'')) & (~strcmp(word,' ')) & (~strcmp(lastword,' ')) & (~strcmp(word,''))
            bigram = char(strcat(lastword, {' '}, word));
            % for debugging :) fprintf('DDD%sDDDDDD%sDDD%sDDD\n',bigram,lastword,word);
            if isKey(g, bigram)
                g(bigram) = g(bigram)+1;
            else
                bigs(bigram) =1;
                g(bigram) = 1;

            end
        end
        lastword = word;

    end

end

selectedheaders =containers.Map();
gkeys = keys(g);

for i=1:size(gkeys,2)
    if g(gkeys{i})>=n

        selectedheaders(gkeys{i})=1;

    end

end
headers = keys(selectedheaders);

outputMatrix = zeros(size(inputcellarray,1),length(headers));
for i = 1:size(inputcellarray,1)
    fprintf('%d/%d ', i, size(inputcellarray,1));
    comment = inputcellarray{i};
    comment = SanitizeComment(comment);
    comment = lower(comment);

    r=regexp(comment,' ','split');
    comment = [];
    for j =1:size(r,2)
        if doStem
            word = porterStemmer(cell2mat(r(j)));
        else
            word = (cell2mat(r(j)));
        end

        comment = [comment,' ',word];
    end
    outputMatrix(i,:) = term_count(comment, headers);

    if mod(i,300)==0
        a = sprintf('%d', i);
        disp(a)
    end

end

featureVector = outputMatrix;
selectedheaderskeys = keys(selectedheaders);
%csvwrite('forWeka_featuresonly.csv', outputMatrix);
disp('');

end

author: faridani

repo: MatlabNLP

save_path: github-repos/MATLAB/faridani-MatlabNLP


text:
function mesh = genMesh3D(domain, nx, ny, nz)

%% Usage: mesh structure of a uniform cubic partition
%
% INPUTS:
% domain --- cubic domain = [xmin, xmax, ymin, ymax, zmin, zmax].
% nx --- the number of uniform partition in x direction.
% ny --- the number of uniform partition in y direction.
% nz --- the number of uniform partition in z direction.
%
% option.meshinfo --- basic (default): generate only p and t.
%                     all: enriched mesh information
%
% OUTPUTS:
% mesh --- a struct data contains mesh information.
%
% Last Modified: 08/07/2020 by Xu Zhang
%
%                A8-------------------A7        The Cube is divided into
%                /|                   /|        six congruent tetrahedrons
%               / |                  / |
%              /  |                 /  |
%             /   |                /   |        (1) A1-A2-A3-A7
%            /    |               /    |        (2) A1-A6-A2-A7
%          A5-----+-------------A6     |        (3) A1-A5-A6-A7
%           |     |             |      |        (4) A1-A8-A5-A7
%           |     |             |      |        (5) A1-A4-A8-A7
%           |     |             |      |        (6) A1-A3-A4-A7
%           |     A4------------+------A3
%           |     /             |      /
%           |    /              |     /
%           |   /               |    /
%           |  /                |   /
%           | /                 |  /
%           |/                  | /
%          A1-------------------A2
%
%% 1. Generate basic mesh info: p t
[p,T] = genMesh3DRectPT(domain, nx, ny, nz);
c1 = [1,2,3,7];
c2 = [1,6,2,7];
c3 = [1,5,6,7];
c4 = [1,8,5,7];
c5 = [1,4,8,7];
c6 = [1,3,4,7];
t = [T(:,c1); T(:,c2); T(:,c3); T(:,c4); T(:,c5); T(:,c6)];
mesh = struct('p',p,'t',t,'T',T);

repo: ifem


text:
function OUT = StarVariable(DATA)
% Computes foreign variables as in the GVAR with equal weights
% =======================================================================
% [OUT, wgts] = xstar(DATA)
% -----------------------------------------------------------------------
% INPUTS
%   - DATA: an (T x N) matrix of time series
%--------------------------------------------------------------------------
% OUTPUT
%    - OUT: (T x N) matrix of xstar time series
% =======================================================================
% Ambrogio Cesa Bianchi, March 2015
% ambrogio.cesabianchi@gmail.com

[nobs, nvar] = size(DATA);
OUT = nan(nobs, nvar);
w = ones(1,nvar)./nvar;

for ii=1:nvar
    % Create a logical vector to exclude variable ii
    select = ones(1,nvar); select(ii)=0; select = logical(select);
    % Remove variable ii from DATA
    data = DATA(:,select);
    % Remove weight/flow ii from w
    weights = w(select);
    % Compute xstar variable for cuntry ii
    OUT(:,ii) = CrossWeightAverage(data,weights);
end

author: ambropo

repo: VAR-Toolbox

save_path: github-repos/MATLAB/ambropo-VAR-Toolbox


text:
function o = spget(options,name,default,flag)
% SPGET   Get sparse interpolation OPTIONS parameters
%   VAL = SPGET(OPTIONS,'NAME') extracts the value of the named property
%   from the sparse grid options structure OPTIONS, returning an
%   empty matrix if the property value is not specified in
%   OPTIONS. It is sufficient to type only the leading characters
%   that uniquely identify the property. Case is ignored for
%   property names. [] is a valid OPTIONS argument.
%
%   VAL = SPGET(OPTIONS,'NAME',DEFAULT) extracts the named property as
%   above, but returns VAL = DEFAULT if the named property is not
%   specified in OPTIONS. For example
%
%       val = spget(opts,'RelTol',1e-2);
%
%   returns val = 1e-2 if the RelTol property is not specified in opts.
%
%   See also SPSET, SPINTERP, SPVALS.

% Author : Andreas Klimke, Universitaet Stuttgart
% Version: 1.5
% Date   : March 1, 2006

% Change log:
% V1.0   : September 23, 2003
%          Initial version
% V1.1   : January 27, 2004
%          Added "KeepGrid" and "KeepFunctionValues" options for
%          later use.
% V1.2   : Apr 22, 2004
%          Added options for dimension-adaptive sparse grids.
% V1.3   : Feb 22, 2005
%          Added handling of sparse sets of indices (useful for
%          high-dimensional interpolation)
% V1.4   : Feb 2, 2006
%          Added droptol property
% V1.5   : March 1, 2006
%          Added EnableFFT feature

% Note: SPGET is similar in syntax and code to the options
% handling with ODEGET and ODESET of the MATLAB ODE suite by Marc
% Reichelt and Lawrence Shampine.

% ------------------------------------------------------------
% Sparse Grid Interpolation Toolbox
% Copyright (c) 2006 W. Andreas Klimke, Universitaet Stuttgart
% Copyright (c) 2007-2008 W. A. Klimke. All Rights Reserved.
% See LICENSE.txt for license.
% email: klimkeas@ians.uni-stuttgart.de
% web  : http://www.ians.uni-stuttgart.de/spinterp
% ------------------------------------------------------------

% undocumented usage for fast access with no error checking
if (nargin == 4) && isequal(flag,'fast')
   o = getknownfield(options,name,default);
   return
end

if nargin < 2
  error('Not enough input arguments.');
end
if nargin < 3
  default = [];
end

if ~isempty(options) && ~isa(options,'struct')
  error('First argument must be an options structure created with SPSET.');
end

if isempty(options)
  o = default;
  return;
end

Names = {'GridType', 'RelTol', 'AbsTol', 'Vectorized', 'MinDepth', ...
				 'MaxDepth', 'VariablePositions', 'NumberOfOutputs', ...
				 'PrevResults', 'FunctionArgType', 'KeepFunctionValues', ...
				 'KeepGrid' 'DimensionAdaptive', 'MinPoints', ...
				 'MaxPoints', 'DimadaptDegree', 'DegreeStrategy', ...
				 'SparseIndices', 'DropTol', 'EnableDCT'};

m = length(Names);
names = cell(m,1);
for k = 1:m
	names{k} = lower(Names{k});
end

lowName = lower(name);
matched = strmatch(lowName,names);
if isempty(matched)               % if no matches
	msg = sprintf(['Unrecognized property name ''%s''.  ' ...
    'See SPSET for possibilities.'], name);
  error(msg);
elseif length(matched) > 1
	msg = sprintf('Ambiguous property name ''%s'' ', name);
	msg = [msg '(' Names{matched}];
	for k = j(2:length(matched))'
		msg = [msg ', ' Names{matched}];
	end
	msg = sprintf('%s).', msg);
	error(msg);
end

if any(strcmp(fieldnames(options),Names{matched}))
  o = options.(Names{matched});
  if isempty(o)
    o = default;
  end
else
  o = default;
end

% -------------------------------------------------------------------
function v = getknownfield(s, f, d)
% GETKNOWNFIELD   Get field f from struct s, or else yield default
% d.

if isfield(s,f)   % s could be empty.
  v = subsref(s, struct('type','.','subs',f));
  if isempty(v)
    v = d;
  end
else
  v = d;
end

author: johannesgerer

repo: jburkardt-m


text:
function binary = sample_bernoulli(probabilities)
    global report_calls_to_sample_bernoulli
    if report_calls_to_sample_bernoulli,
        fprintf('sample_bernoulli() was called with a matrix of size %d by %d. ', size(probabilities, 1), size(probabilities, 2));
    end
    seed = sum(probabilities(:));
    binary = +(probabilities > a4_rand(size(probabilities), seed)); % the "+" is to avoid the "logical" data type, which just confuses things.
end

repo: neural-nets


text:
function data=create_data()

mydat=load('cih.dat');

data=ts('1948Q1',mydat,{'LG_C','LG_I','LG_H'});
end

author: jmaih


text:
function [VertFacesConn, FaceConn] = tess_faceconn(Faces)
% TESS_FACECONN: Computes faces connectivity.
%
% USAGE:  [VertFacesConn, FaceConn] = tess_faceconn(Faces);
%
% INPUT:
%     - Faces    : Nx3 double matrix
% OUTPUT:
%     - FacesConn : sparse matrix [nVertices x nFaces]

% @=============================================================================
% This function is part of the Brainstorm software:
% https://neuroimage.usc.edu/brainstorm
%
% Copyright (c) University of Southern California & McGill University
% This software is distributed under the terms of the GNU General Public License
% as published by the Free Software Foundation. Further details on the GPLv3
% license can be found at http://www.gnu.org/copyleft/gpl.html.
%
% FOR RESEARCH PURPOSES ONLY. THE SOFTWARE IS PROVIDED "AS IS," AND THE
% UNIVERSITY OF SOUTHERN CALIFORNIA AND ITS COLLABORATORS DO NOT MAKE ANY
% WARRANTY, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO WARRANTIES OF
% MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE, NOR DO THEY ASSUME ANY
% LIABILITY OR RESPONSIBILITY FOR THE USE OF THIS SOFTWARE.
%
% For more information type "brainstorm license" at command prompt.
% =============================================================================@
%
% Authors: Anand Joshi, Dimitrios Pantazis, November 2007
%          Francois Tadel, 2008-2010

% Check matrices orientation
if (size(Faces, 2) ~= 3)
    error('Faces must have 3 columns (X,Y,Z).');
end

% Build VertFacesConn
nFaces = size(Faces,1);
rowno = double([Faces(:,1); Faces(:,2); Faces(:,3)]);
colno = [1:nFaces, 1:nFaces, 1:nFaces]';
data  = ones(3*nFaces, 1);
VertFacesConn = sparse(rowno,colno,data);

% Build FacesConn
if (nargout > 1)
    FaceConn = (VertFacesConn' * VertFacesConn) > 0;
end

author: brainstorm-tools


text:
function test_bug2342

% MEM 2gb
% WALLTIME 00:10:00
% DEPENDENCY

% first create some data
%--------------------------------------------------------
% make 3 channels with no direct link between 1 and 2
cfg             = [];
cfg.ntrials     = 100;
cfg.triallength = 0.5;
cfg.fsample     = 200;
cfg.nsignal     = 3;
cfg.method      = 'ar';
cfg.params(:,:,1) = [ 0.8 0   0;
                      0   0.9 0.5;
                      0.4 0   0.5];
cfg.params(:,:,2) = [-0.5    0  0;
                        0 -0.8  0;
                        0    0 -0.2];
cfg.noisecov      = [0.3 0 0;
                       0 1 0;
                       0 0 0.2];

data            = ft_connectivitysimulation(cfg);
data2           = data;
data2.label     = {'signal001b';'signal002b';'signal003b'};
data            = ft_appenddata([], data, data2);

% according to Martin, padding should get rid of most the zigzags.

% freqanalysis
cfgf           = [];
cfgf.method    = 'mtmfft';
cfgf.output    = 'fourier';
cfgf.tapsmofrq = 4;
freq           = ft_freqanalysis(cfgf, data);
cfgf.padding   = 10;
freqpad        = ft_freqanalysis(cfgf, data);

% connectivityanalysis
cfg = [];
cfg.method           = 'granger';
cfg.granger.sfmethod = 'bivariate';
g    = ft_connectivityanalysis(cfg, freq);
gpad = ft_connectivityanalysis(cfg, freqpad);

%...but it doesn't

% now make it more realistic and add a bit of noise to signals 4-6 (to make
% it not perfectly collinear).
data2 = data;
for k = 1:numel(data.trial)
  data2.trial{k}(4:6,:) = data.trial{k}(4:6,:)+randn(3,100)*0.000001;
end

cfgf           = [];
cfgf.method    = 'mtmfft';
cfgf.output    = 'fourier';
cfgf.tapsmofrq = 2;
freq2           = ft_freqanalysis(cfgf, data2);
cfgf.padding   = 10;
freq2pad        = ft_freqanalysis(cfgf, data2);

% connectivityanalysis
cfg = [];
cfg.method           = 'granger';
cfg.granger.sfmethod = 'bivariate';
g2    = ft_connectivityanalysis(cfg, freq2);
g2pad = ft_connectivityanalysis(cfg, freq2pad);

% conclusion: adding a tiny bit of noise removes the zigzags.

author: fieldtrip

repo: fieldtrip


text:
function ebsd = affinetrans(ebsd, A, b)
% perform an affine transformation on spatial ebsd data
%
% Input
%  ebsd - @EBSD
%  A    - transformation matrix or homogeneous coordinates, e.g. [1 0;0 1]  or  [1 0 dy; 0 1 dx; 0 0 1 ]
%  b    - shift term
%
% Output
%  transformed ebsd - @EBSD

% set up transformation matrix
if isa(A,'images.geotrans.internal.GeometricTransformation')

  xy = A.transformPointsInverse([ebsd.prop.x,ebsd.prop.y]);
  ebsd.prop.x = xy(:,1);
  ebsd.prop.y = xy(:,2);

else
  if all(size(A) == [3 3])
    T = A;
  elseif nargin < 3
    T(1:2,1:2) = A;
    T(3,3) = 1;
  else
    if ~isempty(A)
      T(1:2,1:2) = A;
    else
      T(1:2,1:2) = eye(2);
    end
    T(1:2,3) = b(:);
    T(3,3) = 1;
  end

  % rotate the spatial data
  xy = [ebsd.prop.x(:), ebsd.prop.y(:), ones(length(ebsd),1)] * T';

  ebsd.prop.x = xy(:,1);
  ebsd.prop.y = xy(:,2);

  % rotate the unit cells
  T(1:2,3) = 0; % no shift!
  if ~isempty(ebsd.unitCell)
    xy = [ebsd.unitCell ones(length(ebsd.unitCell),1)] * T';
    ebsd.unitCell = xy(:,1:2);
  end
end

% after affine transformation we loose any grid
ebsd = EBSD(ebsd);

author: mtex-toolbox

repo: mtex

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


text:
function M = feval(L, n, flag)
%FEVAL    Deprecated function, provided for limited backward compatibility.

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

warning('CHEBFUN:LINOP:feval:deprecated',...
    ['This function is provided only for limited backward compatibility.',...
    ' Use MATRIX instead.'] );
warning('off', 'CHEBFUN:LINOP:feval:deprecated')

if ( prod(size(L)) > 1 ) %#ok<PSIZE>
    error('CHEBFUN:LINOP:feval:multivariable', ...
        'This syntax is not available for multivariable systems.')
end

if ( nargin < 3 )
    flag = 'bc';
end
% TODO: Should we always use a chebcolloc2 representation here?
% if ( nargin < 4 )
%     pref = cheboppref();
%     discType = pref.discretization;
% elseif ( isa(discType, 'cheboppref') )
%     discType = discType.discretization;
% elseif ( ischar(discType) )
%     discType = str2func(discType);
% end
discType = @chebcolloc2;

if ( isa(n, 'chebfun') )
    M = L*n;   % application to a function

else

    disc = discType(L);
    disc.dimension = n;

    if ( strcmp(flag, 'oldschool') )
        disc.dimAdjust = 0;
    end

    [PA, ignored, B, A] = matrix(disc);

    % Depending on the flag, we will do different things about boundary
    % conditions.
    switch(flag)
        case 'nobc'
            M = A;
        case 'bc'
            M = [B; PA];
        otherwise  % oldschool

            k = size(B,1);     % number of rows to drop
            k2 = ceil(k/2);    % about half for top
            try
                A = cell2mat(A);
            catch
                error('CHEBFUN:LINOP:feval:oldschool', ...
                    'oldschool does not support this problem');
            end
            if ( k2 > 0 )
                A(1:k2,:) = B(1:k2,:);
                krem = k - k2;     % number remaining
                A(end-krem+1:end, :) = B(k2+1:end, :);
            end
            M = A;
    end
end

end

author: chebfun

repo: chebfun

save_path: github-repos/MATLAB/chebfun-chebfun


text:
function bord = copyMakeBorder(imgOriginal, t)

[linhas, cols] = size(imgOriginal);

tmp = [fliplr(imgOriginal(:,2:t+1)) imgOriginal fliplr(imgOriginal(:,cols-t:cols-1))];

bord = [flipud(tmp(2:t+1,:)); tmp; flipud(tmp(linhas-t:linhas-1,:))];

end

author: zhangqianqianQQ

repo: MachineVisionAlgorithm


text:
function mrGrayConvertClassToNifti(leftClassFile, rightClassFile, vAnatFile, t1NiftiFile)
% Convert mrGray class files to the NIFTI format used by ITKGray
%
%   mrGrayConvertClassToNifti(leftClassFile, rightClassFile, vAnatFile, t1NiftiFile)
%
% The vAnatFile and the t1NiftiFile may not align.  The NIFTI data are
% always in AC-PC coordinate frame and the old vAnatomy.dat files were not
% placed in this frame.  Hence, we find the transformation between the
% vAnatomy.dat and the new T1-NIFTI and apply that to the classifiction
% data.
%
% The left and right classification data are stored in a single file.
%
% See also:
%   analyze2mrGray, mrGray2Analyze
%
% Example:
%   mrGrayConvertClassToNifti - User is queried for files
%
% HISTORY:
% 2007.09.19 RFD wrote it.

labels = mrGrayGetLabels;

% false here will do 6-parameter rigid body alignment, true will do
% 9-parameter aligninment to account for scale differences.
adjustScales = true;

%% Get input file names
if(~exist('leftClassFile','var'))
    opts = {'*.class;*.Class;*.CLASS', 'Class files (*.class)'; '*.*','All Files (*.*)'};
    [f, p]=uigetfile(opts, 'Pick a mrGray class file for the left hemisphere');
    if(isequal(f,0)|| isequal(p,0))
        leftClassFile = '';
    else
        if(isempty(p)), p = pwd; end
        leftClassFile = fullFile(p,f);
    end
end
if(~exist('rightClassFile','var'))
    default = fullfile(fileparts(leftClassFile),'right.class');
    opts = {'*.class;*.Class;*.CLASS', 'Class files (*.class)'; '*.*','All Files (*.*)'};
    [f, p]=uigetfile(opts, 'Pick a mrGray class file for the right hemisphere',default);
    if(isequal(f,0)|| isequal(p,0))
        rightClassFile = '';
    else
        if(isempty(p)), p = pwd; end
        rightClassFile = fullFile(p,f);
    end
end
if(~exist('vAnatFile','var')||isempty(vAnatFile))
    default = fullfile(fileparts(rightClassFile),'vAnatomy.dat');
    opts = {'*.dat;*.Dat;*.DAT', 'MrGray dat files (*.dat)'; '*.nii;*.nii.gz', 'NIFTI files (*.nii,*.nii.gz)'; '*.*','All Files (*.*)'};
    [f, p]=uigetfile(opts, 'Pick a mrGray vAnatomy file',default);
    if(isequal(f,0)|| isequal(p,0)), disp('user canceled.'); return; end
    vAnatFile = fullFile(p,f);
end

if(~exist('t1NiftiFile','var')||isempty(t1NiftiFile))
    default = fullfile(fileparts(vAnatFile),'vAnatomy.nii.gz');
    opts = {'*.nii;*.nii.gz', 'NIFTI files (*.nii,*.nii.gz)'; '*.*','All Files (*.*)'};
    [f, p]=uigetfile(opts, 'Pick a NIFTI anatomy file',default);
    if(isequal(f,0)|| isequal(p,0)), disp('user canceled.'); return; end
    t1NiftiFile = fullFile(p,f);
end

[p,f,e] = fileparts(t1NiftiFile);
if(isempty(p)), p = pwd; end
if(strcmpi(e,'.gz')), [junk,f] = fileparts(f); end
outFileBaseName = fullfile(p, f);
outFileName = [outFileBaseName '_class.nii.gz'];
if(exist(outFileName,'file'))
    [f,p] = uiputfile('*.nii.gz','Save new NIFTI class file as...',outFileName);
    if(isequal(f,0)|| isequal(p,0)), disp('user canceled.'); return; end
    outFileName = fullfile(p,f);
    [p,f,e] = fileparts(outFileName);
    if(strcmpi(e,'.gz')), [junk,f] = fileparts(f); end
    outFileBaseName = fullfile(p, f);
end
disp(['NIFTI classification will be saved in ' outFileName]);

%% Read the data
% Allow t1NiftiFile and mrGrayClassFile to be passed as either filenames or
% structs that have already been loaded from file.
if(~isstruct(t1NiftiFile))
    t1NiftiFile = niftiRead(t1NiftiFile);
end

[vAnat,vAnatMm] = readVolAnat(vAnatFile);
if(all(vAnatMm==0))
    warning('vAnat voxel size is indeterminate - using [1 1 1]...'); %#ok<WNTAG>
    vAnatMm = [1 1 1];
end

%% Find the xForm between the vAnatomy data and the T1 NIFTI file data
[xformVAnatToAcpc] = dtiXformVanatCompute(t1NiftiFile.data, t1NiftiFile.qto_xyz, vAnat, vAnatMm, [], adjustScales);

% Create a bounding box identical to the nifti file bb.
bb = mrAnatXformCoords(t1NiftiFile.qto_xyz, [1 1 1; t1NiftiFile.dim]);
xform    = inv(xformVAnatToAcpc);
mm       = t1NiftiFile.pixdim;
interp   = [1 1 1 0 0 0];
outXform = t1NiftiFile.qto_xyz;
sz       = t1NiftiFile.dim;
vAnatSize = size(vAnat); %#ok<NASGU>
save([outFileBaseName '_xformToVanat.mat'],'sz','bb','mm','xform','vAnatSize');

% If you are concerned and want to check the alignments at this point use:
% newAnat = mrAnatResliceSpm(vAnat, xform, bb, mm, interp);
% showMontage(mrAnatHistogramClip(double(newAnat),0.4,0.99)-mrAnatHistogramClip(double(t1NiftiFile.data),0.4,0.99));
clear t1NiftiFile vAnat;
newClass = zeros(sz,'uint8');

% Get the left class file
if(~isempty(leftClassFile))
    if(~isstruct(leftClassFile))
        leftClassFile = readClassFile(leftClassFile, 0, 0);
    end
    disp('Resampling the left classification volumes...');
    cd = permute(leftClassFile.data,[2 1 3]);
    clear leftClassFile;
    types = unique(cd(:));
    types(types==0) = [];
    newTypes = mrGrayConvertClassTypeToLabel(types, 'l', labels);
    for(ii=1:numel(types))
    	c = double(cd==types(ii));
        c = mrAnatResliceSpm(c, xform, bb, mm, interp, false);
        newClass(c>=0.5) = newTypes(ii);
    end
    clear cd type newTypes;
end

% Get the right class file
if(~isempty(rightClassFile))
    if(~isstruct(rightClassFile))
        rightClassFile = readClassFile(rightClassFile, 0, 0);
    end
    disp('Resampling the right classification volumes...');
    cd = permute(rightClassFile.data,[2 1 3]);
    clear rightClassFile;
    types = unique(cd(:));
    types(types==0) = [];
    newTypes = mrGrayConvertClassTypeToLabel(types, 'r', labels);
    for(ii=1:numel(types))
    	c = double(cd==types(ii));
        c = mrAnatResliceSpm(c, xform, bb, mm, interp, false);
        newClass(c>=0.5) = newTypes(ii);
    end
    clear cd type;
end

% Both the right and left classification data are saved in a single file.
% They are transformed into the T1 NIFTI format.
disp(['Saving ' outFileName '...']);
dtiWriteNiftiWrapper(newClass, outXform, outFileName, 1, 'mrGray2 segmentation');

return;

author: vistalab

repo: vistasoft

save_path: github-repos/MATLAB/vistalab-vistasoft


text:
function [net,stats] = cnn_train_dag_modified(net, imdb, getBatch, varargin)
%CNN_TRAIN_DAG Demonstrates training a CNN using the DagNN wrapper
%    CNN_TRAIN_DAG() is similar to CNN_TRAIN(), but works with
%    the DagNN wrapper instead of the SimpleNN wrapper.

% Copyright (C) 2014-16 Andrea Vedaldi.
% 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).

opts.expDir = fullfile('data','exp') ;
opts.continue = true ;
opts.batchSize = 256 ;
opts.numSubBatches = 1 ;
opts.train = [] ;
opts.val = [] ;
opts.gpus = [] ;
opts.prefetch = false ;
opts.numEpochs = 300 ;
opts.learningRate = 0.001 ;
opts.weightDecay = 0.0005 ;
opts.momentum = 0.9 ;
opts.saveMomentum = true ;
opts.nesterovUpdate = false ;
opts.randomSeed = 0 ;
opts.profile = false ;
opts.parameterServer.method = 'mmap' ;
opts.parameterServer.prefix = 'mcn' ;

opts.derOutputs = {'objective', 1} ;
opts.extractStatsFn = @extractStats ;
opts.plotStatistics = true;
opts = vl_argparse(opts, varargin) ;

if ~exist(opts.expDir, 'dir'), mkdir(opts.expDir) ; end
if isempty(opts.train), opts.train = find(imdb.images.set==1) ; end
if isempty(opts.val), opts.val = find(imdb.images.set==2) ; end
if isnan(opts.train), opts.train = [] ; end
if isnan(opts.val), opts.val = [] ; end

% -------------------------------------------------------------------------
%                                                            Initialization
% -------------------------------------------------------------------------

evaluateMode = isempty(opts.train) ;
if ~evaluateMode
  if isempty(opts.derOutputs)
    error('DEROUTPUTS must be specified when training.\n') ;
  end
end

% -------------------------------------------------------------------------
%                                                        Train and validate
% -------------------------------------------------------------------------

modelPath = @(ep) fullfile(opts.expDir, sprintf('net-epoch-%d.mat', ep));
modelFigPath = fullfile(opts.expDir, 'net-train.pdf') ;

start = opts.continue * findLastCheckpoint(opts.expDir) ;
if start >= 1
  fprintf('%s: resuming by loading epoch %d\n', mfilename, start) ;
  [net, state, stats] = loadState(modelPath(start)) ;
else
  state = [] ;
end

for epoch=start+1:opts.numEpochs

  % Set the random seed based on the epoch and opts.randomSeed.
  % This is important for reproducibility, including when training
  % is restarted from a checkpoint.

  rng(epoch + opts.randomSeed) ;
%   prepareGPUs(opts, epoch == start+1) ;

  % Train for one epoch.
  params = opts ;
  params.epoch = epoch ;
  params.learningRate = opts.learningRate(min(epoch, numel(opts.learningRate))) ;
  params.train = opts.train(randperm(numel(opts.train))) ; % shuffle
  params.val = opts.val(randperm(numel(opts.val))) ;
  params.imdb = imdb ;
  params.getBatch = getBatch ;

  if numel(opts.gpus) <= 1
    [net, state] = processEpoch(net, state, params, 'train') ;
    [net, state] = processEpoch(net, state, params, 'val') ;

    if size(imdb.images.labels,1) > 1
        [net.meta.svm_model1,net.meta.svm_err1] = training_SVR(net,params);
    end

    if ~evaluateMode
      saveState(modelPath(epoch), net, state) ;
    end
    lastStats = state.stats ;
  else
    spmd
      [net, state] = processEpoch(net, state, params, 'train') ;
      [net, state] = processEpoch(net, state, params, 'val') ;
      if labindex == 1 && ~evaluateMode
        saveState(modelPath(epoch), net, state) ;
      end
      lastStats = state.stats ;
    end
    lastStats = accumulateStats(lastStats) ;
  end

  stats.train(epoch) = lastStats.train ;
  stats.val(epoch) = lastStats.val ;
  clear lastStats ;
  saveStats(modelPath(epoch), stats) ;

%   if opts.plotStatistics
%     switchFigure(1) ; clf ;
%     plots = setdiff(...
%       cat(2,...
%       fieldnames(stats.train)', ...
%       fieldnames(stats.val)'), {'num', 'time'}) ;
%     for p = plots
%       p = char(p) ;
%       values = zeros(0, epoch) ;
%       leg = {} ;
%       for f = {'train', 'val'}
%         f = char(f) ;
%         if isfield(stats.(f), p)
%           tmp = [stats.(f).(p)] ;
%           values(end+1,:) = tmp(1,:)' ;
%           leg{end+1} = f ;
%         end
%       end
%       subplot(1,numel(plots),find(strcmp(p,plots))) ;
%       plot(1:epoch, values','o-') ;
%       xlabel('epoch') ;
%       title(p) ;
%       legend(leg{:}) ;
%       grid on ;
%     end
%     drawnow ;
%     print(1, modelFigPath, '-dpdf') ;
%   end
end

% With multiple GPUs, return one copy
if isa(net, 'Composite'), net = net{1} ; end

% -------------------------------------------------------------------------
function [net, state] = processEpoch(net, state, params, mode)
% -------------------------------------------------------------------------
% Note that net is not strictly needed as an output argument as net
% is a handle class. However, this fixes some aliasing issue in the
% spmd caller.

% initialize with momentum 0
if isempty(state) || isempty(state.momentum)
  state.momentum = num2cell(zeros(1, numel(net.params))) ;
end

% move CNN  to GPU as needed
numGpus = numel(params.gpus) ;
if numGpus >= 1
  net.move('gpu') ;
  state.momentum = cellfun(@gpuArray, state.momentum, 'uniformoutput', false) ;
end
if numGpus > 1
  parserv = ParameterServer(params.parameterServer) ;
  net.setParameterServer(parserv) ;
else
  parserv = [] ;
end

% profile
if params.profile
  if numGpus <= 1
    profile clear ;
    profile on ;
  else
    mpiprofile reset ;
    mpiprofile on ;
  end
end

num = 0 ;
epoch = params.epoch ;
subset = params.(mode) ;
adjustTime = 0 ;

stats.num = 0 ; % return something even if subset = []
stats.time = 0 ;

start = tic ;
for t=1:params.batchSize:numel(subset)
  fprintf('%s: epoch %02d: %3d/%3d:', mode, epoch, ...
          fix((t-1)/params.batchSize)+1, ceil(numel(subset)/params.batchSize)) ;
  batchSize = min(params.batchSize, numel(subset) - t + 1) ;

  for s=1:params.numSubBatches
    % get this image batch and prefetch the next
    batchStart = t + (labindex-1) + (s-1) * numlabs ;
    batchEnd = min(t+params.batchSize-1, numel(subset)) ;
    batch = subset(batchStart : params.numSubBatches * numlabs : batchEnd) ;
    num = num + numel(batch) ;
    if numel(batch) == 0, continue ; end

    inputs = params.getBatch(params.imdb, batch) ;

    if params.prefetch
      if s == params.numSubBatches
        batchStart = t + (labindex-1) + params.batchSize ;
        batchEnd = min(t+2*params.batchSize-1, numel(subset)) ;
      else
        batchStart = batchStart + numlabs ;
      end
      nextBatch = subset(batchStart : params.numSubBatches * numlabs : batchEnd) ;
      params.getBatch(params.imdb, nextBatch) ;
    end

    if strcmp(mode, 'train')
      net.mode = 'normal' ;
      net.accumulateParamDers = (s ~= 1) ;
      net.eval(inputs, params.derOutputs, 'holdOn', s < params.numSubBatches) ;
    else
      net.mode = 'test' ;
      net.eval(inputs) ;
    end
  end

  % Accumulate gradient.
  if strcmp(mode, 'train')
    if ~isempty(parserv), parserv.sync() ; end
    state = accumulateGradients(net, state, params, batchSize, parserv) ;
  end

  % Get statistics.
  time = toc(start) + adjustTime ;
  batchTime = time - stats.time ;
  stats.num = num ;
  stats.time = time ;
  stats = params.extractStatsFn(stats,net) ;
  currentSpeed = batchSize / batchTime ;
  averageSpeed = (t + batchSize - 1) / time ;
  if t == 3*params.batchSize + 1
    % compensate for the first three iterations, which are outliers
    adjustTime = 4*batchTime - time ;
    stats.time = time + adjustTime ;
  end

  fprintf(' %.1f (%.1f) Hz', averageSpeed, currentSpeed) ;
  for f = setdiff(fieldnames(stats)', {'num', 'time'})
    f = char(f) ;
    fprintf(' %s: %.4f', f, stats.(f)) ;
  end
  fprintf('\n') ;
end

% Save back to state.
state.stats.(mode) = stats ;
if params.profile
  if numGpus <= 1
    state.prof.(mode) = profile('info') ;
    profile off ;
  else
    state.prof.(mode) = mpiprofile('info');
    mpiprofile off ;
  end
end
if ~params.saveMomentum
  state.momentum = [] ;
else
  state.momentum = cellfun(@gather, state.momentum, 'uniformoutput', false) ;
end

net.reset() ;
net.move('cpu') ;

% -------------------------------------------------------------------------
function state = accumulateGradients(net, state, params, batchSize, parserv)
% -------------------------------------------------------------------------
numGpus = numel(params.gpus) ;
otherGpus = setdiff(1:numGpus, labindex) ;

for p=1:numel(net.params)

  if ~isempty(parserv)
    parDer = parserv.pullWithIndex(p) ;
  else
    parDer = net.params(p).der ;
  end

  switch net.params(p).trainMethod

    case 'average' % mainly for batch normalization
      thisLR = net.params(p).learningRate ;
      net.params(p).value = vl_taccum(...
          1 - thisLR, net.params(p).value, ...
          (thisLR/batchSize/net.params(p).fanout),  parDer) ;

    case 'gradient'
      thisDecay = params.weightDecay * net.params(p).weightDecay ;
      thisLR = params.learningRate * net.params(p).learningRate ;

      if thisLR>0 || thisDecay>0
        % Normalize gradient and incorporate weight decay.
        parDer = vl_taccum(1/batchSize, parDer, ...
                           thisDecay, net.params(p).value) ;

        % Update momentum.
        state.momentum{p} = vl_taccum(...
          params.momentum, state.momentum{p}, ...
          -1, parDer) ;

        % Nesterov update (aka one step ahead).
        if params.nesterovUpdate
          delta = vl_taccum(...
            params.momentum, state.momentum{p}, ...
            -1, parDer) ;
        else
          delta = state.momentum{p} ;
        end

        % Update parameters.
        net.params(p).value = vl_taccum(...
          1,  net.params(p).value, thisLR, delta) ;
      end
    otherwise
      error('Unknown training method ''%s'' for parameter ''%s''.', ...
        net.params(p).trainMethod, ...
        net.params(p).name) ;
  end
end

% -------------------------------------------------------------------------
function stats = accumulateStats(stats_)
% -------------------------------------------------------------------------

for s = {'train', 'val'}
  s = char(s) ;
  total = 0 ;

  % initialize stats stucture with same fields and same order as
  % stats_{1}
  stats__ = stats_{1} ;
  names = fieldnames(stats__.(s))' ;
  values = zeros(1, numel(names)) ;
  fields = cat(1, names, num2cell(values)) ;
  stats.(s) = struct(fields{:}) ;

  for g = 1:numel(stats_)
    stats__ = stats_{g} ;
    num__ = stats__.(s).num ;
    total = total + num__ ;

    for f = setdiff(fieldnames(stats__.(s))', 'num')
      f = char(f) ;
      stats.(s).(f) = stats.(s).(f) + stats__.(s).(f) * num__ ;

      if g == numel(stats_)
        stats.(s).(f) = stats.(s).(f) / total ;
      end
    end
  end
  stats.(s).num = total ;
end

% -------------------------------------------------------------------------
function stats = extractStats(stats, net)
% -------------------------------------------------------------------------
sel = find(cellfun(@(x) isa(x,'dagnn.Loss'), {net.layers.block})) ;
for i = 1:numel(sel)
  stats.(net.layers(sel(i)).outputs{1}) = net.layers(sel(i)).block.average ;
end

% -------------------------------------------------------------------------
function saveState(fileName, net_, state)
% -------------------------------------------------------------------------
net = net_.saveobj() ;
% save(fileName, 'net', 'state') ;
save(fileName, 'net') ;

% -------------------------------------------------------------------------
function saveStats(fileName, stats)
% -------------------------------------------------------------------------
if exist(fileName)
  save(fileName, 'stats', '-append') ;
else
  save(fileName, 'stats') ;
end

% -------------------------------------------------------------------------
function [net, state, stats] = loadState(fileName)
% -------------------------------------------------------------------------
load(fileName, 'net', 'state', 'stats') ;
net = dagnn.DagNN.loadobj(net) ;
if isempty(whos('stats'))
  error('Epoch ''%s'' was only partially saved. Delete this file and try again.', ...
        fileName) ;
end

% -------------------------------------------------------------------------
function epoch = findLastCheckpoint(modelDir)
% -------------------------------------------------------------------------
list = dir(fullfile(modelDir, 'net-epoch-*.mat')) ;
tokens = regexp({list.name}, 'net-epoch-([\d]+).mat', 'tokens') ;
epoch = cellfun(@(x) sscanf(x{1}{1}, '%d'), tokens) ;
epoch = max([epoch 0]) ;

% -------------------------------------------------------------------------
function switchFigure(n)
% -------------------------------------------------------------------------
if get(0,'CurrentFigure') ~= n
  try
    set(0,'CurrentFigure',n) ;
  catch
    figure(n) ;
  end
end

% -------------------------------------------------------------------------
function clearMex()
% -------------------------------------------------------------------------
clear vl_tmove vl_imreadjpeg ;

% -------------------------------------------------------------------------
function prepareGPUs(opts, cold)
% -------------------------------------------------------------------------
numGpus = numel(opts.gpus) ;
if numGpus > 1
  % check parallel pool integrity as it could have timed out
  pool = gcp('nocreate') ;
  if ~isempty(pool) && pool.NumWorkers ~= numGpus
    delete(pool) ;
  end
  pool = gcp('nocreate') ;
  if isempty(pool)
    parpool('local', numGpus) ;
    cold = true ;
  end

end
if numGpus >= 1 && cold
  fprintf('%s: resetting GPU\n', mfilename)
  clearMex() ;
  if numGpus == 1
    gpuDevice(opts.gpus)
  else
    spmd
      clearMex() ;
      gpuDevice(opts.gpus(labindex))
    end
  end
end

function [svm_model,svm_err] = training_SVR(net,params)

move(net, 'gpu')
net.mode = 'test';
batch = 64;
sampling = floor(length(params.train) / 5000);
scores_train = params.imdb.images.score(1:sampling:end);
data = params.imdb.images.data(:,:,:,1:sampling:end);
num_batch = floor(length(scores_train) / batch);
index = find(arrayfun(@(a) strcmp(a.name, 'myLoss'), net.layers)==1);
tlayerName = net.layers(index).inputs{1};
for i = 1:num_batch
    images = gpuArray(data(:,:,:,(i-1)*batch+1:i*batch));
    images = bsxfun(@minus, single(images), net.meta.normalization.averageImage) ;
    net.eval({'data', images});
    if strcmp(net.layers(end).block.lossType,'CE')
        feats(:,(i-1)*batch+1:i*batch) = squeeze(gather(vl_nnsoftmax(net.vars(net.getVarIndex(tlayerName)).value))) ;
    else
        feats(:,(i-1)*batch+1:i*batch) = squeeze(gather(net.vars(net.getVarIndex(tlayerName)).value)) ;
    end
end
if num_batch*batch < length(scores_train)
    images = gpuArray(data(:,:,:,num_batch*batch+1:length(scores_train)));
    images = bsxfun(@minus, single(images), net.meta.normalization.averageImage) ;
    net.eval({'data', images}) ;
    if strcmp(net.layers(end).block.lossType,'CE')
        feats(:,num_batch*batch+1:length(scores_train)) = squeeze(gather(vl_nnsoftmax(net.vars(net.getVarIndex(tlayerName)).value)));
    else
        feats(:,num_batch*batch+1:length(scores_train)) = squeeze(gather(net.vars(net.getVarIndex(tlayerName)).value));
    end
end
svm_model = svmtrain(scores_train', double(feats'), '-s 4 -t 0');
[scores_pred, ~,~] = svmpredict(scores_train', double(feats'), svm_model,'-q');
svm_err = (sum(abs(scores_pred - scores_train'))) / length(scores_train);

author: HuiZeng


text:
function n = vectorNorm3d(v)
%VECTORNORM3D Norm of a 3D vector or of set of 3D vectors
%
%   N = vectorNorm3d(V);
%   Returns the norm of vector V.
%
%   When V is a N-by-3 array, compute norm for each vector of the array.
%   Vector are given as rows. Result is then a N-by-1 array.
%
%   NOTE: compute only euclidean norm.
%
%   See Also
%   vectors3d, normalizeVector3d, vectorAngle3d, hypot3
%

%   ---------
%   author : David Legland
%   INRA - TPV URPOI - BIA IMASTE
%   created the 21/02/2005.

%   HISTORY
%   19/06/2009 rename as vectorNorm3d

n = sqrt(sum(v.*v, 2));

author: rpng

repo: lips

save_path: github-repos/MATLAB/rpng-lips


text:
function varargout=loadbnm(fname)
% a compact binary matrix format.
%
% [M1,M2,M3,...]=loadbnm(fname);
%
% see savebnm
%
% Aslak Grinsted April 2004

% -------------------------------------------------------------------------
%   Copyright (C) 2002-2004, Aslak Grinsted
%   This software may be used, copied, or redistributed as long as it is not
%   sold and this copyright notice is reproduced on each copy made.  This
%   routine is provided as is without any express or implied warranties
%   whatsoever.

if (length(strfind(fname,'.'))==0)
    fname=[fname '.bnm'];
end

skip=(nargout==0); % should it skip loading data and just show info for the file?

header='BNM Matrix File';

[fid,message]=fopen(fname,'r','ieee-le');
if (fid<0)
    message=strrep(message,'Sorry. No help in figuring out the problem . . .','');
    error(['Could not open the specified file. ' message ' file=' fname])
end

hdr=[fread(fid,length(header),'uchar')' ''];
if ~strcmp(hdr,header)
    error([fname 'is not a ' header ' file!']);
end

version=fread(fid,1,'int16');
if version>2
    warning('File format version > 2')
end
if (version==1)
    nout=1;
else
    nout=fread(fid,1,'int16');
end

nn=min(nargout,nout);
if skip
    disp(sprintf('File: "%s" , Version: %i , MatrixCount : %i\n',fname,version,nout));
    nn=nout;
end

for ii=1:nn
    infos=fread(fid,2,'int16');
    sz=fread(fid,infos(1),'uint32')';
    rl=infos(2); %isreal???
    if skip
        disp(sprintf('Matrix %i,  size: [%s], isimag: %i',ii,num2str(sz),~rl))
        fseek(fid,prod(sz)*(1+(~rl))*4,0);
    else
        varargout{ii}=fread(fid,prod(sz),'float');
        if ~rl %if imaginary numbers...
            varargout{ii}=varargout{ii}+fread(fid,prod(sz),'float')*i;
        end
        varargout{ii}=reshape(varargout{ii},sz);
        %varargout{ii}=M;
    end
end

fclose(fid);

function s=loadstring(fid)
ll=fread(fid,1,'uint32');
s=[fread(fid,ll,'uchar')' ''];

author: buzsakilab

repo: buzcode

save_path: github-repos/MATLAB/buzsakilab-buzcode


text:
function grad = B_repmat(future_layers, curr_layer)
future_grad = GetFutureGrad(future_layers, curr_layer);

sourceDims = curr_layer.sourceDims;
% targetDims = curr_layer.targetDims;

if length(sourceDims)==1
    grad = sum(future_grad,2);
else
    % to be implemented
end

end

author: singaxiong

repo: SignalGraph