text:
% This small example illustrates how to use the remote API
% synchronous mode. The synchronous mode needs to be
% pre-enabled on the server side. You would do this by
% starting the server (e.g. in a child script) with:
%
% simRemoteApi.start(19999,1300,false,true)
%
% But in this example we try to connect on port
% 19997 where there should be a continuous remote API
% server service already running and pre-enabled for
% synchronous mode.
%
% IMPORTANT: for each successful call to simxStart, there
% should be a corresponding call to simxFinish at the end!

function irb140_abb_path_planning()
    global configuration

    % init basic data. Robots and number of collision objects.
    coppelia = [];
    coppelia.n_collision_objects = 0;
    coppelia.robots{1}.n_joints = 6;
    coppelia.robots{1}.end_effector.n_joints = 2;
    coppelia.dt = 50/1000; % default 50 ms
    %should match ARTE delta_time
    configuration.delta_time = 50/1000;
    %configuration.delta_time = 0.01;
    %coppelia = coppelia_start(coppelia);

    pick_and_place(coppelia)
    %coppelia_stop(coppelia);
end

function pick_and_place(coppelia)
%load in arte
global robot
robot = load_robot('ABB', 'IRB140');
robot_number = 1;
% initial position
Q = [0 0 0 0 0 0]';
robot.q = Q;
robot.qd = [0 0 0 0 0 0]';
lin_vel = 0.5; % m/s
ang_vel = 0.8; % rad/s

% %Initial point
% T1 = [-1 0 0 0.52;
%        0 1 0  0;
%        0 0 -1 0.4;
%        0 0 0 1];
% % Final point
% T2 = [1 0 0 0.5;
%        0 1 0  0.5;
%        0 0 1 0.5;
%        0 0 0 1];
% % Final point
% T3 = [-1 0 0 0.5;
%        0 1 0  -0.5;
%        0 0 -1 0.5;
%        0 0 0 1];
q2 = [0.5 0.5 0.5 0.5 0.5 0.5];
qd2 = [0 0 0 0 0 0];
% end speed is 20 % of max speed for all joints
[q, qdd]=AbsJPath(robot, q2, qd2, 20);
move_robot(coppelia, 1, q, qdd)

% go to pre pick point, then open, then to pick point
%MoveQ(coppelia, robot_number, T1)
%open_gripper(coppelia, robot_number);
%coppelia_wait(coppelia, 10)
%MoveL(coppelia, robot_number, T2, lin_vel, ang_vel)
%MoveL(coppelia, robot_number, T3, lin_vel, ang_vel)
%MoveL(coppelia, robot_number, T1, lin_vel, ang_vel)

% % now close gripper
% close_gripper(coppelia, robot_number);
% MoveQ(coppelia, robot_number, T1)
%
%
% %release
% MoveQ(coppelia, robot_number, T3)
% open_gripper(coppelia, robot_number);
% % over the release point
% MoveQ(coppelia, robot_number, T4)
% %start point again
% MoveAbsQ(coppelia, robot_number, Q)
end

% % Linear interpolation move to Q
% function MoveAbsQ(coppelia, robot_number, Q)
%     global robot
%     q0 = robot.q;
%     %q1 = [0.0 0.0 0.0 0.0 0.0 0.0]';
%     q_path = linear_q_path(q0, Q);
%     setjointtargetpositions(coppelia,robot_number, q_path);
%     %coppelia_wait(coppelia, 10)
%     robot.q = Q;
% end

% % MoveL
% % MoveQ
% function MoveQ(coppelia, robot_number, T)
% %drawrobot3d(robot, q0)
%     global robot
%     %start joint coordinates, as saved before
%     q0 = robot.q;
%
%     qinv = inversekinematic(robot, T);
%     q1 = qinv(:,1);
%     q_path = linear_q_path(q0, q1);
%     setjointtargetpositions(coppelia, robot_number, q_path)
%     %coppelia_wait(coppelia, 10)
%     robot.q = q1;
% end

function MoveL(coppelia, robot_number, Tf, lin_vel, ang_vel)

    global robot
    %start joint coordinates, as saved before
    q0 = robot.q;
    T0 = directkinematic(robot, q0);
    p0 = T0(1:3, 4);
    pf = Tf(1:3, 4);

    qinv = inversekinematic(robot, Tf);
    qf = qinv(:,1);
    Q0 = T2quaternion(T0);
    Qf = T2quaternion(Tf);

    v0 = abs(lin_vel)*compute_speed(p0, pf);
    w0 = abs(ang_vel)*compute_ang_speed(Q0, Qf);
    Vref = [v0' w0']';
    qs = [robot.q];
    qds = [robot.qd];
    q = q0;
    while 1
        Ti = directkinematic(robot, q);
        pi = Ti(1:3, 4);
        e = norm(pf - pi);
        if e < 0.05
            break
        end
        J = manipulator_jacobian(robot, q);
        qd = inv(J)*Vref;
        q = q + qd*coppelia.dt;
        qs = [qs q];
        qds = [qds qd];
    end
    %q_path = linear_q_path(q0, q1);
    move_robot(coppelia, robot_number, qs, qds);
    %coppelia_wait(coppelia, 10)
    robot.q = q;
    robot.qd = qd;
end

function q_path = linear_q_path(q0, q1)
q_path=[];
% find max difference
dq = abs(q0-q1);
[y, i] = max(dq);
n = floor(y/0.1 + 1);
for i=1:length(q0)
    q_path = [q_path; linspace(q0(i), q1(i), n)];
end
%q_path
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%En base a la posici�n y orientaci�n final, calcular cu�les deben ser las
%velocidades...
% Esto es diferente a calcular la velocidades cuando ya hay contacto y se
% trata de un problema de control... pero es parecido
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [v] = compute_speed(Pi, Pf)
%compute a constant linear speed till target
v = (Pf-Pi);
v = v(:)/norm(v);

end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%En base a la posici�n y orientaci�n final, calcular cu�les deben ser las
%velocidades...
% Esto es diferente a calcular la velocidades cuando ya hay contacto y se
% trata de un problema de control... pero es parecido
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [w] = compute_ang_speed(Qi, Qf)
%compute a constant angular speed till target
%asume the movement is performed in 1 second
w = angular_w_between_quaternions(Qi, Qf, 1);
w = w(:);
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Compute angular speed w that moves Q0 into Q1 in time total_time.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function w = angular_w_between_quaternions(Q0, Q1, total_time)
%global robot
%below this number, the axis is considered as [1 0 0]
%this is to avoid numerical errors
%this is the actual error allowed for w
%epsilon_len = robot.parameters.epsilonQ;
epsilon_len = 0.0001;
%Let's first find quaternion q so q*q0=q1 it is q=q1/q0
%For unit length quaternions, you can use q=q1*Conj(q0)
Q = qprod(Q1, qconj(Q0));

%To find rotation velocity that turns by q during time Dt you need to
%convert quaternion to axis angle using something like this:
len=sqrt(Q(2)^2 + Q(3)^2 + Q(4)^2);

if len > epsilon_len
    angle=2*atan2(len, Q(1));
    axis=[Q(2) Q(3) Q(4)]./len;
else
    angle=0;
    axis=[1 0 0];
end
w=axis*angle/total_time;
end

repo: ARTE


text:
function y = guard_interval(Ng,Nfft,NgType,ofdmSym)
if NgType==1
  y=[ofdmSym(Nfft-Ng+1:Nfft) ofdmSym(1:Nfft)];
elseif NgType==2
  y=[zeros(1,Ng) ofdmSym(1:Nfft)];
end

author: LyricYang


text:
function [surface_handles, pcl, ncl] = surface(r, varargin)
% Render image data on brain surfaces; options for cutaways and canonical surfaces
% Surface method for region object - renders blobs on multiple types of 3-D
% surfaces
%
% :Usage:
% ::
%
%    [surface_handles, pcl, ncl] = surface(r, [optional inputs])
%
%    - accepts combinations of existing surface handles and one or more
%    keywords accepted by addbrain.m to build surfaces
%
% ..
%     Author and copyright information:
%
%     Copyright (C) 2013 Tor Wager
%
%     This program is free software: you can redistribute it and/or modify
%     it under the terms of the GNU General Public License as published by
%     the Free Software Foundation, either version 3 of the License, or
%     (at your option) any later version.
%
%     This program is distributed in the hope that it will be useful,
%     but WITHOUT ANY WARRANTY; without even the implied warranty of
%     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
%     GNU General Public License for more details.
%
%     You should have received a copy of the GNU General Public License
%     along with this program.  If not, see <http://www.gnu.org/licenses/>.
% ..
%
% :Inputs:
%
%   **r:**
%        A region object
%
%   **'surface_handles':**
%        Followed by existing surface handles to render on
%
%   **addbrain.m keywords:**
%        Keywords for surfaces accepted by addbrain.m
%        A list is below:
%
%         % % Regions
%         % -----------------------------------------------------------------------
%         'vmpfc' 'nacc' 'BST' 'cau' 'caudate' 'put' 'GP' 'GPe' 'GPi' 'VeP' ...
%         'cm' 'md' 'stn' 'habenula' 'mammillary' 'hypothalamus','hy','hythal' ...
%         'midbrain' 'pag' 'PBP' 'sn' 'SNc' 'SNr' 'VTA' 'rn' ...
%         'pbn' 'lc' 'rvm' 'rvm_old' 'nts' 'sc' 'ic' 'drn' 'mrn' ...
%         'thalamus' 'thal' 'LGN' 'lgn' 'MGN' 'mgn' 'VPthal', 'VPLthal', 'VPL', 'intralaminar_thal', ...
%         'medullary_raphe' 'spinal_trigeminal' 'nuc_ambiguus' 'dmnx_nts' 'ncs_B6_B8' 'nrp_B5' 'pbn' 'ncf' 'vep' 'PBP'
%         'amygdala' 'amygdala hires' 'hippocampus', 'hipp' 'hippocampus hires'
%
%         % Cortical Surfaces
%         % -----------------------------------------------------------------------
%         'left' 'hires left' 'surface left' 'hires surface left' ...
%         'right' 'hires right' 'surface right' 'hires surface right' ...
%         'transparent_surface' 'foursurfaces' 'flat left'  'flat right' ...
%
%         % Macro subcortical surfaces
%         % -----------------------------------------------------------------------
%         'CIT168' 'cerebellum','cblm' 'brainstem' 'suit brainstem'
%
%         % Cutaways
%         % -----------------------------------------------------------------------
%         'brainbottom' 'cutaway', 'left_cutaway' 'right_cutaway' ...
%         'left_insula_slab' 'right_insula_slab' 'accumbens_slab' 'coronal_slabs' 'coronal_slabs_4' 'coronal_slabs_5' ...
%
%         % Groups
%         % -----------------------------------------------------------------------
%         'bg', 'basal ganglia' 'midbrain_group' 'limbic' 'limbic hires' 'brainstem_group' 'thalamus_group'
%
%   **render_on_surface keywords:**
%   allowable = {'clim' 'color' 'colormap' 'colormapname' 'axis_handle' 'pos_colormap' 'neg_colormap'};
%   'color'         followed by [r g b] triplet for solid color
%   'clim'          Limits for image values mapped to most extreme colors
%   'colormap'      Followed by the name of a Matlab colormap to use (
%   'colormapname'
%   'axis_handle'
%   'pos_colormap'
%   'neg_colormap'
%
% :Outputs:
%
%   **all_surf_handles:**
%        surface patch handles
%
%   **pcl:**
%        region object with positive-only clusters
%
%   **ncl:**
%        region object with negative-only clusters
%
%
% :Examples:
% ::
%
%    % Render on combined cortical cutaway and subcortical surfaces
%    Use surface(r), with optional arguments taken by surface_cutaway:
%    poscm = colormap_tor([1 .3 0], [1 1 0]); % orange to yellow
%    [all_surf_handles, pcl, ncl] = surface(r, 'cutaway', 'ycut_mm', -30, 'pos_colormap', poscm, 'existingfig');
%
%    % Render on four cortical surfaces
%    [all_surf_handles2, pcl, ncl] = surface(r, 'foursurfaces', 'pos_colormap', poscm, 'neg_colormap', negcm);
%    [all_surf_handles2, pcl, ncl] = surface(r, 'foursurfaces', 'existingfig', 'color_upperboundpercentile', 95, 'color_lowerboundpercentile', 5, 'neg_colormap', colormap_tor([0 0 1], [.4 0 .7]));
%
%    % use mediation_brain_surface_figs and re-make colors [OLD METHOD! SLOW!]
%    all_surf_handles = mediation_brain_surface_figs([]);
%    surface(r, 'cutaway', 'surface_handles', all_surf_handles, 'color_upperboundpercentile', 95, 'color_lowerboundpercentile', 5, 'neg_colormap', colormap_tor([0 0 1], [.2 0 .5]));
%
%    % Make a region of interest surface and render on that
%    create_figure; p = addbrain('thalamus'); lightRestoreSingle;
%    [all_surf_handles2, pcl, ncl] = surface(region(t_age), 'color_upperboundpercentile', 95, 'color_lowerboundpercentile', 5, 'neg_colormap', colormap_tor([0 0 1], [.4 0 .7]), 'surface_handles', p);
%
% figure; hh = addbrain('foursurfaces');
% t.surface('surface_handles', hh);
%
% Note:
% To erase rendered blobs, use:
% sh = addbrain('eraseblobs', sh);
%
% :See also:*
%
% surface_cutaway, cluster_surf, mediation_brain_surface_figs

% ..
%    DEFAULTS AND INPUTS
% ..

existingfig = false;
surface_handles = [];

addbrain_allowable_args = {'vmpfc' 'nacc' 'BST' 'cau' 'caudate' 'put' 'GP' 'GPe' 'GPi' 'VeP' ...
    'cm' 'md' 'stn' 'habenula' 'mammillary' 'hypothalamus','hy','hythal' ...
    'midbrain' 'pag' 'PBP' 'sn' 'SNc' 'SNr' 'VTA' 'rn' ...
    'pbn' 'lc' 'rvm' 'rvm_old' 'nts' 'sc' 'ic' 'drn' 'mrn' ...
    'thalamus' 'thal' 'LGN' 'lgn' 'MGN' 'mgn' 'VPthal', 'VPLthal', 'VPL', 'intralaminar_thal', ...
    'medullary_raphe' 'spinal_trigeminal' 'nuc_ambiguus' 'dmnx_nts' 'ncs_B6_B8' 'nrp_B5' 'pbn' 'ncf' 'vep' 'PBP' ...
    'left' 'hires left' 'surface left' 'hires surface left' ...
    'right' 'hires right' 'surface right' 'hires surface right' ...
    'transparent_surface' 'foursurfaces' 'flat left'  'flat right' ...
    'brainbottom' 'cutaway', 'left_cutaway' 'right_cutaway' ...
    'left_insula_slab' 'right_insula_slab' 'accumbens_slab' 'coronal_slabs' 'coronal_slabs_4' 'coronal_slabs_5' ...
    'brainstem' 'suit brainstem' 'amygdala' 'amygdala hires' 'hippocampus', 'hipp' 'hippocampus hires' 'cerebellum','cblm' 'CIT168' ...
    'bg', 'basal ganglia' 'midbrain_group' 'limbic' 'limbic hires' 'brainstem_group' 'thalamus_group'};

render_on_surface_allowable_args = {'clim' 'cmaprange' 'color' 'colormap' 'colormapname' 'axis_handle' 'pos_colormap' 'neg_colormap' ...
    'summer' 'cool' 'hot' 'bone' 'copper' 'prism' 'hsv' 'winter' 'nolegend'};

% optional inputs with default values
for i = 1:length(varargin)
    if ischar(varargin{i})
        switch varargin{i}

            case 'surface_handles', surface_handles = varargin{i+1}; varargin{i} = []; varargin{i+1} = [];
            case {'existingfig', 'nofigure'}, existingfig = true; varargin{i} = [];

            case 'noverbose'

            case addbrain_allowable_args
                % do nothing, handle later

            case render_on_surface_allowable_args
                % no nothing, will pass in to render_on_surface

            otherwise, warning(['Unknown input string option:' varargin{i}]);
        end
    end
end

% BUILD BRAIN SURFACE BASED ON ALLOWABLE INPUTS TO ADDBRAIN
% -------------------------------------------------------------------------

for i = 1:length(varargin)

    if ischar(varargin{i}) && any(strcmp(addbrain_allowable_args, varargin{i}))

        surface_handles = [surface_handles addbrain(varargin{i})];

    end

end

if isempty(surface_handles)

    surface_handles = addbrain('left_cutaway');

end

% END DEFAULTS AND INPUTS
% -------------------------------------------------------------------------

% Legacy
if nargout > 1
    [pcl, ncl] = posneg_separate(r, 'Z');
end

% Run color change
% Convert to image vector first
% -------------------------------------------------------------------------
obj = region2imagevec(r);

if isempty(obj)
    disp('No data to display. Returning')
    return
end

if isempty(surface_handles)
    disp('No surface handles to display on. Returning')
    return
end

render_on_surface(obj, surface_handles, varargin{:});

end % function

author: canlab

repo: CanlabCore

save_path: github-repos/MATLAB/canlab-CanlabCore


text:
function obj = runICPSaveCorr(obj, p, g, CP)

% Pre ICP (= first correspondences)
if g.nItICP == 1

    obj.CP.PreICP = CP;

% Post ICP (= final correspondences)
elseif (g.nItICP == p.MaxNoIt) || (obj.D.stats{g.nItICP}.normdx < p.StopConditionNormdx)

    obj.CP.PostICP = CP;

end

end

author: pglira


text:
function [i_min,i_max,j_min,j_max,lew] = get_lew(img,smap_b)
% cut a local window from the original img, and also return
% lew's coordinates in the original img

    [all_i,all_j] = find(smap_b~=0);
    rescale = 0.2;
    i_min = round(min(all_i)*(1-rescale));
    i_max = round(max(all_i)*(1+rescale));
    j_min = round(min(all_j)*(1-rescale));
    j_max = round(max(all_j)*(1+rescale));

    img_size = size(img(:,:,1));
    [i_min,i_max,j_min,j_max] = handle_cross_boundary(i_min,i_max,j_min,j_max,img_size);
    lew = img(i_min:i_max,j_min:j_max,:);
end

author: zhangqianqianQQ

repo: MachineVisionAlgorithm


text:
% PrintSegyFormatRevsisionNumber: Prints/returns the SEGY format revision number
%
% Call:
%     Revision=PrintSegyFormatRevsisionNumber(SegyHeader);
%     Revision=PrintSegyFormatRevsisionNumber(SegyHeader,vlevel);
%

function Revision=PrintSegyFormatRevsisionNumber(SegyHeader,vlevel);
if nargin<2
    vlevel=0;
end
Revision=str2num(dec2base(SegyHeader.SegyFormatRevisionNumber,16))/100;
SegymatVerbose(sprintf('%s: Using SEG Y revision %3.2f',mfilename,Revision),vlevel)

author: cultpenguin

repo: segymat

save_path: github-repos/MATLAB/cultpenguin-segymat


text:
function jed = transition_to_jed_common ( )

%*****************************************************************************80
%
%% TRANSITION_TO_JED_COMMON returns the Common calendar transition as a JED.
%
%  Discussion:
%
%    In the Common calendar, the last moment of the Julian calendar was
%      11:59 pm, 4 October 1582 Julian/CE,
%      11:59 pm, 14 October 1582 Gregorian.
%    The first minute of the Gregorian calendar ended at
%      12:01 am, 5 October 1582 Julian,
%      12:01 am, 15 October 1582 Gregorian/CE.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    21 December 2012
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Output, real JED, the Julian Ephemeris Date of the date.
%
  jed = 2299160.5;

  return
end

author: johannesgerer

repo: jburkardt-m


text:
function data = crescentfullmoon(N, r1, r2, r3)

if nargin < 1
    N = 1000;
end
if mod(N,4) ~= 0
    N = round(N/4) * 4;
end
if nargin < 2
    r1 = 5;
end
if nargin < 3
    r2 = 10;
end
if nargin < 4
    r3 = 15;
end

N1 = N/4;
N2 = N-N1;

phi1 = rand(N1,1) * 2 * pi;
R1 = sqrt(rand(N1, 1));
moon = [cos(phi1) .* R1 * r1 sin(phi1) .* R1 * r1 zeros(N1,1)];

d = r3 - r2;
phi2 = pi + rand(N2,1) * pi;
R2 = sqrt(rand(N2,1));
crescent = [cos(phi2) .* (r2 + R2 * d) sin(phi2) .* (r2 + R2 * d) ones(N2,1)];

data = [moon; crescent];

repo: gspbox


text:
function [inter, inside]= intersectRayPolygon3d(ray, poly)
%INTERSECTRAYPOLYGON3D Intersection point of a 3D ray and a 3D polygon
%
%   INTER = intersectRayPolygon3d(RAY, POLY)
%   Compute coordinates of intersection point between the 3D ray RAY and
%   the 3D polygon POLY. RAY is a 1-by-6 row vector containing origin and
%   direction vector of the ray, POLY is a Np-by-3 array containing
%   coordinates of 3D polygon vertices.
%   INTER is a 1-by-3 row vector containing coordinates of intersection
%   point, or [NaN NaN NaN] if ray and polygon do not intersect.
%
%   INTERS = intersectRayPolygon3d(RAYS, POLY)
%   If RAYS is a N-by-6 array representing several rays, the result
%   INTERS is a N-by-3 array containing coordinates of intersection of each
%   ray with the polygon.
%
%   [INTER INSIDE] = intersectRayPolygon3d(RAY, POLY)
%   Also return a N-by-1 boolean array containing TRUE if both the polygon
%   and the corresponding ray contain the intersection point.
%
%   Example
%     % Compute intersection between a 3D ray and a 3D triangle
%     pts3d = [3 0 0; 0 6 0;0 0 9];
%     ray1 = [0 0 0 3 6 9];
%     inter = intersectRayPolygon3d(ray1, pts3d)
%     inter =
%           1   2   3
%
%     % keep only valid intersections with several rays
%     pts3d = [3 0 0; 0 6 0;0 0 9];
%     rays = [0 0 0 3 6 9;10 0 0 1 2 3;3 6 9 3 6 9];
%     [inter inside] = intersectRayPolygon3d(rays, pts3d);
%     inter(inside, :)
%     ans =
%           1   2   3
%
%   See Also
%   intersectRayPolygon, intersectLinePolygon3d, intersectLineTriangle3d
%

% ------
% Author: David Legland
% e-mail: david.legland@inra.fr
% Created: 2011-05-22,    using Matlab 7.9.0.529 (R2009b)
% Copyright 2011 INRA - Cepia Software Platform.

% supporting plane of polygon vertices
plane   = createPlane(poly(1:3, :));

% intersection of 3D ray with the plane
inter   = intersectLinePlane(ray, plane);

% project all points on reference plane
pts2d   = planePosition(projPointOnPlane(poly, plane), plane);
pInt2d  = planePosition(projPointOnPlane(inter, plane), plane);

% need to check polygon orientation
inPoly  = xor(isPointInPolygon(pInt2d, pts2d), polygonArea(pts2d) < 0);
onRay   = linePosition3d(inter, ray) >= 0;
inside  = inPoly & onRay;

% intersection points outside the polygon are set to NaN
inter(~inside, :) = NaN;

author: rpng

repo: lips

save_path: github-repos/MATLAB/rpng-lips


text:
%TRPLOT Plot a 3D coordinate frame
%
% TRPLOT(T, OPTIONS) draws a 3D coordinate frame represented by the SE(3) homogeneous
% transform T (4x4).
%
% H = TRPLOT(T, OPTIONS) as above but returns a handle.
%
% TRPLOT(R, OPTIONS) as above but the coordinate frame is rotated about the
% origin according to the orthonormal rotation matrix R (3x3).
%
% H = TRPLOT(R, OPTIONS) as above but returns a handle.
%
% H = TRPLOT() creates a default frame EYE(3,3) at the origin and returns a
% handle.
%
% Animation::
%
% Firstly, create a plot and keep the the handle as per above.
%
% TRPLOT(T, 'handle', H) moves the coordinate frame described by the handle H to
% the pose T (4x4).
%
% Options::
% 'handle',h             Update the specified handle
% 'axhandle',A           Draw in the MATLAB axes specified by the axis handle A
%
% 'color',C              The color to draw the axes, MATLAB ColorSpec
% 'axes'                 Show the MATLAB axes, box and ticks (default true)
% 'axis',A               Set dimensions of the MATLAB axes to A=[xmin xmax ymin ymax zmin zmax]
% 'frame',F              The coordinate frame is named {F} and the subscript on the axis labels is F.
% 'framelabel',F         The coordinate frame is named {F}, axes have no subscripts.
% 'framelabeloffset',O   Offset O=[DX DY] frame labels in units of text box height
% 'text_opts', opt       A cell array of MATLAB text properties
% 'length',s             Length of the coordinate frame arms (default 1)
% 'thick',t              Thickness of lines (default 0.5)
% 'text'                 Enable display of X,Y,Z labels on the frame (default true)
% 'labels',L             Label the X,Y,Z axes with the 1st, 2nd, 3rd character of the string L
% 'rgb'                  Display X,Y,Z axes in colors red, green, blue respectively
% 'rviz'                 Display chunky rviz style axes%
% 'arrow'                Use arrows rather than line segments for the axes
% 'width', w             Width of arrow tips (default 1)
%
% 'perspective'          Display the axes with perspective projection (default off)
% '3d'                   Plot in 3D using anaglyph graphics
% 'anaglyph',A           Specify anaglyph colors for '3d' as 2 characters for
%                        left and right (default colors 'rc'): chosen from
%                        r)ed, g)reen, b)lue, c)yan, m)agenta.
% 'dispar',D             Disparity for 3d display (default 0.1)
% 'view',V               Set plot view parameters V=[az el] angles, or 'auto'
%                        for view toward origin of coordinate frame
% 'lefty'                Draw left-handed frame (dangerous)
%
% Examples::
%
%       trplot(T, 'frame', 'A')
%       trplot(T, 'frame', 'A', 'color', 'b')
%       trplot(T1, 'frame', 'A', 'text_opts', {'FontSize', 10, 'FontWeight', 'bold'})
%       trplot(T1, 'labels', 'NOA');
%
%       h = trplot(T, 'frame', 'A', 'color', 'b');
%       trplot(h, T2);
%
% 3D anaglyph plot
%       trplot(T, '3d');
%
% Notes::
% - Multiple frames can be added using the HOLD command
% - When animating a coordinate frame it is best to set the axis bounds initially.
% - The 'rviz' option is equivalent to 'rgb', 'notext', 'noarrow',
%   'thick', 5.
% - The 'arrow' option requires https://www.mathworks.com/matlabcentral/fileexchange/14056-arrow3

%## 3d homogeneous graphics

% Copyright (C) 1993-2019 Peter I. Corke
%
% This file is part of The Spatial Math Toolbox for MATLAB (SMTB).
%
% Permission is hereby granted, free of charge, to any person obtaining a copy
% of this software and associated documentation files (the "Software"), to deal
% in the Software without restriction, including without limitation the rights
% to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies
% of the Software, and to permit persons to whom the Software is furnished to do
% so, subject to the following conditions:
%
% The above copyright notice and this permission notice shall be included in all
% copies or substantial portions of the Software.
%
% THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
% IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS
% FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR
% COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER
% IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN
% CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
%
% https://github.com/petercorke/spatial-math

%TODO:
% 'rviz', chunky RGB lines, no arrows

% Copyright (C) 1993-2017, by Peter I. Corke
%
% This file is part of The Robotics Toolbox for MATLAB (RTB).
%
% RTB 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.
%
% RTB 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 RTB.  If not, see <http://www.gnu.org/licenses/>.
%
% http://www.petercorke.com

% TODO
%  need to decide how to handle scaling
%  what does hold on mean?  don't touch scaling?

function hout = trplot(T, varargin)

    if nargin == 0
        T = eye(3,3);
    end

    if size(T,3) > 1
        error('SMTB:trplot:badarg', 'trplot cannot operate on a sequence');
    end
    if isa(T, 'SE3')
        T = T.T;
    elseif ~ishomog(T) && ~isrot(T)
        error('SMTB:trplot:badarg', 'trplot operates only on transform (4x4) or rotation matrix (3x3)');
    end

    opt.color = 'b';
    opt.textcolor = [];
    opt.rgb = false;
    opt.axes = true;
    opt.axis = [];
    opt.frame = [];
    opt.framelabel = [];
    opt.text_opts = [];
    opt.view = [];
    opt.width = 1;
    opt.arrow = false;
    opt.labels = 'XYZ';
    opt.axhandle = [];
    opt.anaglyph = 'rc';
    opt.d_3d = false;
    opt.dispar = 0.1;
    opt.thick = 0.5;
    opt.length = 1;
    opt.text = true;
    opt.lefty = false;
    opt.rviz = false;
    opt.framelabeloffset = [0 0];
    opt.handle = [];
    opt.perspective = false;

    [opt,args] = tb_optparse(opt, varargin);

    if opt.arrow && ~exist('arrow3')
        opt.arrow = false;
        warning('SMTB:trplot:badarg', 'arrow option requires arrow3 from FileExchange');
    end

    if isscalar(T) && ishandle(T)
        warning('SMTB:trplot:deprecated', 'Use ''handle'' option');
        % trplot(H, T)
        opt.handle = T; T = args{1};
    end

    % ensure it's SE(3)
    if isrot(T)
        T = r2t(T);
    end

    if ~isempty(opt.handle)
        set(opt.handle, 'Matrix', T);
        % for the 3D case retrieve the right hgtransform and set it
        hg2 = get(opt.handle, 'UserData');
        if ~isempty(hg2)
            set(hg2, 'Matrix', T);
        end
        if nargout > 0
            hout = opt.handle;
        end
        return;
    end

    if opt.rviz
        opt.thick = 5;
        opt.arrow = false;
        opt.rgb = true;
        opt.text = false;
    end

    if opt.rgb && opt.d_3d
        error('SMTB:trplot:badarg', 'cannot specify ''rgb'' and ''3d'', use ''anaglyph'' option');
    end
    if isempty(opt.textcolor)
        opt.textcolor = opt.color;
    end
    if isempty(opt.text_opts)
        opt.text_opts = {};
    end

    if opt.d_3d
        opt.color = ag_color(opt.anaglyph(1));
    end

    % figure the dimensions of the axes, if not given
    if isempty(opt.axis)
        % determine some default axis dimensions

        % get the origin of the frame
        if isrot(T)
            c = [0 0 0];  % at zero for a rotation matrix
        else
            c = transl(T);
        end

        d = 1.2;
        opt.axis = [c(1)-d c(1)+d c(2)-d c(2)+d c(3)-d c(3)+d];

    end

    if ~isempty(opt.axhandle)
        hax = opt.axhandle;
        hold(hax);
    else
        ih = ishold;
        if ~ih
            % if hold is not on, then clear the axes and set scaling
            cla
            if ~isempty(opt.axis)
                axis(opt.axis);
            end
            daspect([1 1 1]);

            if opt.axes
                xlabel( 'X');
                ylabel( 'Y');
                zlabel( 'Z');
                rotate3d on
            end
            new_plot = true;
        end
        hax = gca;
        hold on
    end
    % hax is the handle for the axis we will work with, either new or
    % passed by option 'handle'

    opt.text_opts = [opt.text_opts, 'Color', opt.color];

    if opt.perspective
        hax.Projection = 'perspective';
    end
    hg = hgtransform('Parent', hax);

    % trplot( Q.R, fmt, color);
    if isrot(T)
        T = r2t(T);
    end

    % create unit vectors
    o =  [0 0 0]';
    x1 = opt.length*[1 0 0]';
    y1 = opt.length*[0 1 0]';
    if opt.lefty
        z1 = opt.length*[0 0 -1]';
    else
        z1 = opt.length*[0 0 1]';
    end

    % draw the axes

    mstart = [o o o]';
    mend = [x1 y1 z1]';

    if opt.rgb
        axcolors = {'r', 'g', 'b'};
    else
        axcolors = { opt.color, opt.color, opt.color};
    end

    if opt.arrow
%         % draw the 3 arrows
%         S = [opt.color num2str(opt.width)];
%         ha = arrow3(mstart, mend, S);
%         for h=ha'
%             set(h, 'Parent', hg);
%         end
          daspect([1,1,1])
          for i=1:3
              ha = arrow3(mstart(i,1:3), mend(i,1:3), axcolors{i}, opt.width);
              set(ha, 'Parent', hg);
          end
    else
        for i=1:3
            plot2([mstart(i,1:3); mend(i,1:3)], 'Color', axcolors{i}, ...
                'LineWidth', opt.thick, ...
                'Parent', hg);
        end
    end

    % label the axes
    if isempty(opt.frame)
        fmt = '%c';
    else
        fmt = sprintf('%%c_{%s}', opt.frame);
    end

    if opt.text
        % add the labels to each axis
        h = text(x1(1), x1(2), x1(3), sprintf(fmt, opt.labels(1)), 'Parent', hg);
        set(h, opt.text_opts{:});

        h = text(y1(1), y1(2), y1(3), sprintf(fmt, opt.labels(2)), 'Parent', hg);
        set(h, opt.text_opts{:});

        h = text(z1(1), z1(2), z1(3), sprintf(fmt, opt.labels(3)), 'Parent', hg);
        set(h, opt.text_opts{:});
    end

    if ~isempty(opt.framelabel)
        opt.frame = opt.framelabel;
    end
    % label the frame
    if ~isempty(opt.frame)
        h = text(o(1), o(2), o(3), ...
            ['\{' opt.frame '\}'], 'Parent', hg);
        set(h, 'VerticalAlignment', 'top', ...
            'HorizontalAlignment', 'center', opt.text_opts{:}, ...
            'FontUnits', 'normalized');
        e = get(h, 'Extent');
        d = e(4); % use height of text box as a scale factor
        e(1:2) = e(1:2) - opt.framelabeloffset * d;
        set(h, 'Position', e(1:2));

    end

    if ~opt.axes
        set(gca, 'visible', 'off');
    end
    if ischar(opt.view) && strcmp(opt.view, 'auto')
        cam = x1+y1+z1;
        view(cam(1:3));
    elseif ~isempty(opt.view)
        view(opt.view);
    end
    if isempty(opt.handle) && ~ih
        grid on
        hold off
    end

    % now place the frame in the desired pose
    set(hg, 'Matrix', T);

    if opt.d_3d
        % 3D display.  The original axes are for the left eye, and we add
        % another set of axes to the figure for the right eye view and
        % displace its camera to the right of that of that for the left eye.
        % Then we recursively call trplot() to create the right eye view.

        left = gca;
        right = axes;

        % compute the offset in world coordinates
        off = -t2r(view(left))'*[opt.dispar 0 0]';
        pos = get(left, 'CameraPosition');

        set(right, 'CameraPosition', pos+off');
        set(right, 'CameraViewAngle', get(left, 'CameraViewAngle'));
        set(right, 'CameraUpVector', get(left, 'CameraUpVector'));
        target = get(left, 'CameraTarget');
        set(right, 'CameraTarget', target+off');

        % set perspective projections
        set(left, 'Projection', 'perspective');
        set(right, 'Projection', 'perspective');

        % turn off axes for right view
        set(right, 'Visible', 'Off');

        % set color for right view
        hg2 = trplot(T, 'color', ag_color(opt.anaglyph(2)));

        % the hgtransform for the right view is user data for the left
        % view hgtransform, we need to change both when we rotate the
        % frame.
        set(hg, 'UserData', hg2);
    end

    % optionally return the handle, for later modification of pose
    if nargout > 0
        hout = hg;
    end
end

function out = ag_color(c)

% map color character to an color triple, same as anaglyph.m

    % map single letter color codes to image planes
    switch c
    case 'r'
        out = [1 0 0];        % red
    case 'g'
        out = [0 1 0];        % green
    case 'b'
        % blue
        out = [0 0 1];
    case 'c'
        out = [0 1 1];        % cyan
    case 'm'
        out = [1 0 1];        % magenta
    case 'o'
        out = [1 1 0];        % orange
    end
end

author: petercorke

repo: spatialmath-matlab


text:
function [perc,n,which] = misclass(Y,Yest)
% The rate of misclassifications.
%
% '[perc,n,which] = misclass(Y,Yest)'
%
%  'Y' contains the real class labels;
%  'Yest' contains the estimated class labels;
%  'perc' is the rate of misclassifications (between 0 and 1);
%  'n' is the number of misclassifications;
%  'which' contains the indices of the misclassificated instances
%     (the first column gives the row, the second the column index)
%
%
% see also:
%    validate, mse, linf, medae, mae

% Copyright (c) 2011,  KULeuven-ESAT-SCD, License & help @ http://www.esat.kuleuven.be/sista/lssvmlab

n = sum(sum(Y~=Yest));
perc = n/numel(Y);
[I,J] = find(Y~=Yest);
which = [J I];

author: peterhcharlton

repo: RRest


text:
function F = eq(X,Y)
% Internal class for constraint lists

superiorto('sdpvar');
superiorto('double');

if isa(X,'sdpvar')
    if is(X,'binary')
        F=iff(X,Y);
        return;
    end
end
if isa(Y,'sdpvar')
    if is(Y,'binary')
        F=iff(Y,X);
        return;
    end
end

error('Equalities can not be used in double-sided constraints')

author: yalmip

repo: YALMIP

save_path: github-repos/MATLAB/yalmip-YALMIP


text:
function z = logMn(x, p)
% Compute log pdf of a multinomial distribution.
% Input:
%   x: d x 1 integer vector
%   p: d x 1 probability
% Output:
%   z: probability density in logrithm scale z=log p(x)
% Written by Mo Chen (sth4nth@gmail.com).
z = gammaln(sum(x)+1)-sum(gammaln(x+1))+dot(x,log(p));

author: PRML

repo: PRMLT

save_path: github-repos/MATLAB/PRML-PRMLT


text:
classdef MergeMertens < handle
    %MERGEMERTENS  Merge exposure sequence to a single image
    %
    % Pixels are weighted using contrast, saturation and well-exposedness
    % measures, than images are combined using laplacian pyramids.
    %
    % The resulting image weight is constructed as weighted average of
    % contrast, saturation and well-exposedness measures.
    %
    % The resulting image doesn't require tonemapping and can be converted to
    % 8-bit image by multiplying by 255, but it's recommended to apply gamma
    % correction and/or linear tonemapping.
    %
    % For more information see [MK07].
    %
    % ## References
    % [MK07]:
    % > Tom Mertens, Jan Kautz, and Frank Van Reeth. "Exposure fusion".
    % > In Computer Graphics and Applications, 2007. PG'07. 15th Pacific
    % > Conference on, pages 382-390. IEEE, 2007.
    %
    % See also: cv.MergeDebevec, cv.MergeRobertson, makehdr
    %

    properties (SetAccess = private)
        % Object ID
        id
    end

    properties (Dependent)
        % contrast measure weight
        ContrastWeight
        % saturation measure weight
        SaturationWeight
        % well-exposedness measure weight
        ExposureWeight
    end

    %% MergeMertens
    methods
        function this = MergeMertens(varargin)
            %MERGEMERTENS  Creates MergeMertens object
            %
            %     obj = cv.MergeMertens()
            %     obj = cv.MergeMertens('OptionName',optionValue, ...)
            %
            % ## Options
            % * __ContrastWeight__ default 1.0
            % * __SaturationWeight__ default 1.0
            % * __ExposureWeight__ default 0.0
            %
            % See also: cv.MergeMertens.process
            %
            this.id = MergeMertens_(0, 'new', varargin{:});
        end

        function delete(this)
            %DELETE  Destructor
            %
            %     obj.delete()
            %
            % See also: cv.MergeMertens
            %
            if isempty(this.id), return; end
            MergeMertens_(this.id, 'delete');
        end
    end

    %% MergeExposures
    methods
        function dst = process(this, src)
            %PROCESS  Merges images
            %
            %     dst = obj.process(src)
            %
            % ## Input
            % * __src__ vector of input images (1- or 3-channels), all of the
            %   same size and `uint8` type.
            %
            % ## Output
            % * __dst__ result image, same size as `src{i}` and `single` type.
            %
            % See also: cv.MergeMertens.MergeMertens
            %
            dst = MergeMertens_(this.id, 'process', src);
        end
    end

    %% Algorithm
    methods (Hidden)
        function clear(this)
            %CLEAR  Clears the algorithm state
            %
            %     obj.clear()
            %
            % See also: cv.MergeMertens.empty, cv.MergeMertens.load
            %
            MergeMertens_(this.id, 'clear');
        end

        function b = empty(this)
            %EMPTY  Returns true if the algorithm is empty
            %
            %     b = obj.empty()
            %
            % ## Output
            % * __b__ Returns true if the object is empty (e.g in the
            %   very beginning or after unsuccessful read).
            %
            % See also: cv.MergeMertens.clear, cv.MergeMertens.load
            %
            b = MergeMertens_(this.id, 'empty');
        end

        function name = getDefaultName(this)
            %GETDEFAULTNAME  Returns the algorithm string identifier
            %
            %     name = obj.getDefaultName()
            %
            % ## Output
            % * __name__ This string is used as top level XML/YML node tag
            %   when the object is saved to a file or string.
            %
            % See also: cv.MergeMertens.save, cv.MergeMertens.load
            %
            name = MergeMertens_(this.id, 'getDefaultName');
        end

        function save(this, filename)
            %SAVE  Saves the algorithm parameters to a file
            %
            %     obj.save(filename)
            %
            % ## Input
            % * __filename__ Name of the file to save to.
            %
            % This method stores the algorithm parameters in the specified
            % XML or YAML file.
            %
            % See also: cv.MergeMertens.load
            %
            MergeMertens_(this.id, 'save', filename);
        end

        function load(this, fname_or_str, varargin)
            %LOAD  Loads algorithm from a file or a string
            %
            %     obj.load(fname)
            %     obj.load(str, 'FromString',true)
            %     obj.load(..., 'OptionName',optionValue, ...)
            %
            % ## Input
            % * __fname__ Name of the file to read.
            % * __str__ String containing the serialized model you want to
            %   load.
            %
            % ## Options
            % * __ObjName__ The optional name of the node to read (if empty,
            %   the first top-level node will be used). default empty
            % * __FromString__ Logical flag to indicate whether the input is a
            %   filename or a string containing the serialized model.
            %   default false
            %
            % This method reads algorithm parameters from the specified XML or
            % YAML file (either from disk or serialized string). The previous
            % algorithm state is discarded.
            %
            % See also: cv.MergeMertens.save
            %
            MergeMertens_(this.id, 'load', fname_or_str, varargin{:});
        end
    end

    %% Getters/Setters
    methods
        function value = get.ContrastWeight(this)
            value = AlignMTB_(this.id, 'get', 'ContrastWeight');
        end
        function set.ContrastWeight(this, value)
            AlignMTB_(this.id, 'set', 'ContrastWeight', value);
        end

        function value = get.SaturationWeight(this)
            value = AlignMTB_(this.id, 'get', 'SaturationWeight');
        end
        function set.SaturationWeight(this, value)
            AlignMTB_(this.id, 'set', 'SaturationWeight', value);
        end

        function value = get.ExposureWeight(this)
            value = AlignMTB_(this.id, 'get', 'ExposureWeight');
        end
        function set.ExposureWeight(this, value)
            AlignMTB_(this.id, 'set', 'ExposureWeight', value);
        end
    end

end

author: kyamagu

repo: mexopencv

save_path: github-repos/MATLAB/kyamagu-mexopencv


text:
%%
 %  Copyright (c) 2014, Facebook, Inc.
 %  All rights reserved.
 %
 %  This source code is licensed under the BSD-style license found in the
 %  LICENSE file in the root directory of this source tree. An additional grant
 %  of patent rights can be found in the PATENTS file in the same directory.
 %
 %%

function [ap,rec,prec] = get_precision_recall(scores,labels,num_truths)
if nargin==2
    num_truths=sum(labels==1);
end
[srt1,srtd]=sort(scores,'descend');
fp=cumsum(labels(srtd)==-1);
tp=cumsum(labels(srtd)==1);
rec=tp/num_truths;
prec=tp./(fp+tp);

%print_precrec(srt1,prec,rec);

mrec=[0 ; rec ; 1];
mpre=[0 ; prec ; 0];
for i=numel(mpre)-1:-1:1
    mpre(i)=max(mpre(i),mpre(i+1));
end
i=find(mrec(2:end)~=mrec(1:end-1))+1;
ap=sum((mrec(i)-mrec(i-1)).*mpre(i));

end

function print_precrec(score,prec,rec)

for rc=0.1:0.1:0.9
    t = find(rec>rc,1);
    fprintf('recall: %4.2f  precision: %4.2f  score: %4.2f\n',rec(t),prec(t),score(t));
end

end

author: facebookarchive

repo: pose-aligned-deep-networks


text:
function Y = lp_fuse(M1, M2, zt, ap, mp)
%Y = fuse_lap(M1, M2, zt, ap, mp) image fusion with laplacian pyramid
%
%    M1 - input image A
%    M2 - input image B
%    zt - maximum decomposition level
%    ap - coefficient selection highpass (see selc.m)
%    mp - coefficient selection base image (see selb.m)
%
%    Y  - fused image

%    (Oliver Rockinger 16.08.99)
% whos
% check inputs
[z1 s1] = size(M1);
[z2 s2] = size(M2);
M1=double(M1);
M2=double(M2);
if (z1 ~= z2) | (s1 ~= s2)
  error('Input images are not of same size');
end;

% define filter
w  = [1 4 6 4 1] / 16;

% cells for selected images
E = cell(1,zt);
% tic
% loop over decomposition depth -> analysis
for i1 = 1:zt
    tic
  % calculate and store actual image size
  [z s]  = size(M1);
  zl(i1) = z; sl(i1)  = s;

  % check if image expansion necessary
  if (floor(z/2) ~= z/2), ew(1) = 1; else, ew(1) = 0; end;
  if (floor(s/2) ~= s/2), ew(2) = 1; else, ew(2) = 0; end;

  % perform expansion if necessary
  if (any(ew))
  	M1 = adb(M1,ew);
  	M2 = adb(M2,ew);
  end;

  % perform filtering
  G1 = conv2(conv2(es2(M1,2), w, 'valid'),w', 'valid');
  G2 = conv2(conv2(es2(M2,2), w, 'valid'),w', 'valid');

  % decimate, undecimate and interpolate
  M1T = conv2(conv2(es2(undec2(dec2(G1)), 2), 2*w, 'valid'),2*w', 'valid');
  M2T = conv2(conv2(es2(undec2(dec2(G2)), 2), 2*w, 'valid'),2*w', 'valid');
%   toc76

% tic
  % select coefficients and store them
  E(i1) = {selc(M1-M1T, M2-M2T, ap)};
% toc
  % decimate
%  tic
  M1 = dec2(G1);
  M2 = dec2(G2);
%   toc
% toc
% feature ('memstats')

end;
% toc
% select base coefficients of last decompostion stage
% tic
M1 = selb(M1,M2,mp);
% toc
% whos
% feature ('memstats')
% toc
% loop over decomposition depth -> synthesis
% tic
for i1 = zt:-1:1
  % undecimate and interpolate
%   tic
  M1T = conv2(conv2(es2(undec2(M1), 2), 2*w, 'valid'), 2*w', 'valid');
  % add coefficients
  M1  = M1T + E{i1};
%   toc
  % select valid image region
  M1 	= M1(1:zl(i1),1:sl(i1));
%   feature ('memstats')
% whos
end;

% feature ('memstats')
% copy image
Y = M1;
% toc

author: xingchenzhang

repo: VIFB

save_path: github-repos/MATLAB/xingchenzhang-VIFB


text:
% =========================================================================
function parser = matlab2tikzInputParser ()
  % Initialize the structure.
  parser = struct ();
  % Public Properties
  parser.Results = {};
  % Enabel/disable parameters case sensitivity.
  parser.CaseSensitive = false;
  % Enable/disable error for parameters not defined by the constructor.
  parser.KeepUnmatched = false;
  % Enable/disable passing arguments in a structure.
  parser.StructExpand = true;
  % Names of parameters defined in input parser constructor.
  parser.Parameters = {};
  % Names of parameters not defined in the constructor.
  parser.Unmatched = struct ();
  % Names of parameters using default values.
  parser.UsingDefaults = {};

  % Handles for functions that act on the object.
  parser.addRequired   = @addRequired;
  parser.addOptional   = @addOptional;
  parser.addParamValue = @addParamValue;
  parser.parse         = @parse;

  % Initialize the parser plan
  parser.plan = {};
end
% =========================================================================
function p = parser_plan (q, arg_type, name, default, validator)
  p = q;
  plan = p.plan;
  if (isempty (plan))
    plan = struct ();
    n = 1;
  else
    n = numel (plan) + 1;
  end
  plan(n).type      = arg_type;
  plan(n).name      = name;
  plan(n).default   = default;
  plan(n).validator = validator;
  p.plan = plan;
end
% =========================================================================
function p = addRequired (p, name, validator)
  p = parser_plan (p, 'required', name, [], validator);
end
% =========================================================================
function p = addOptional (p, name, default, validator)
  p = parser_plan (p, 'optional', name, default, validator);
end
% =========================================================================
function p = addParamValue (p, name, default, validator)
  p = parser_plan (p, 'paramvalue', name, default, validator);
end
% =========================================================================
function p = parse (p, varargin)
  plan = p.plan;
  results = p.Results;
  using_defaults = {};
  if (p.CaseSensitive)
    name_cmp = @strcmp;
  else
    name_cmp = @strcmpi;
  end
  if (p.StructExpand)
    k = find (cellfun (@isstruct, varargin));
    for m = numel(k):-1:1
      n = k(m);
      s = varargin{n};
      c = [fieldnames(s).'; struct2cell(s).'];
      c = c(:).';
      if (n > 1 && n < numel (varargin))
        varargin = horzcat (varargin(1:n-1), c, varargin(n+1:end));
      elseif (numel (varargin) == 1)
        varargin = c;
      elseif (n == 1);
        varargin = horzcat (c, varargin(n+1:end));
      else % n == numel (varargin)
        varargin = horzcat (varargin(1:n-1), c);
      end
    end
  end
  if (isempty (results))
    results = struct ();
  end
  type = {plan.type};
  n = find( strcmp( type, 'paramvalue' ) );
  m = setdiff (1:numel( plan ), n );
  plan = plan ([n,m]);
  for n = 1 : numel (plan)
    found = false;
    results.(plan(n).name) = plan(n).default;
    if (~ isempty (varargin))
      switch plan(n).type
      case 'required'
        found = true;
        if (strcmpi (varargin{1}, plan(n).name))
          varargin(1) = [];
        end
        value = varargin{1};
        varargin(1) = [];
      case 'optional'
        m = find (cellfun (@ischar, varargin));
        k = find (name_cmp (plan(n).name, varargin(m)));
        if (isempty (k) && validate_arg (plan(n).validator, varargin{1}))
          found = true;
          value = varargin{1};
          varargin(1) = [];
        elseif (~ isempty (k))
          m = m(k);
          found = true;
          value = varargin{max(m)+1};
          varargin(union(m,m+1)) = [];
        end
      case 'paramvalue'
        m = find( cellfun (@ischar, varargin) );
        k = find (name_cmp (plan(n).name, varargin(m)));
        if (~ isempty (k))
          found = true;
          m = m(k);
          value = varargin{max(m)+1};
          varargin(union(m,m+1)) = [];
        end
      otherwise
        error( sprintf ('%s:parse', mfilename), ...
               'parse (%s): Invalid argument type.', mfilename ...
             )
      end
    end
    if (found)
      if (validate_arg (plan(n).validator, value))
        results.(plan(n).name) = value;
      else
        error( sprintf ('%s:invalidinput', mfilename), ...
               '%s: Input argument ''%s'' has invalid value.\n', mfilename, plan(n).name ...
             );
      end
      p.Parameters = union (p.Parameters, {plan(n).name});
    elseif (strcmp (plan(n).type, 'required'))
      error( sprintf ('%s:missinginput', mfilename), ...
             '%s: input ''%s'' is missing.\n', mfilename, plan(n).name ...
           );
    else
      using_defaults = union (using_defaults, {plan(n).name});
    end
  end

  if (~ isempty (varargin))
    if (p.KeepUnmatched)
    % Include properties that do not match specified properties
      for n = 1:2:numel(varargin)
        if (ischar (varargin{n}))
          results.(varargin{n}) = varargin{n+1};
        else
          error (sprintf ('%s:invalidinput', mfilename), ...
                 '%s: invalid input', mfilename)
        end
      end
    end
    p.Unmatched = varargin{1:2:end};
  end

  % Store the results of the parsing
  p.Results = results;
  p.UsingDefaults = using_defaults;

end
% =========================================================================
function result = validate_arg (validator, arg)
  try
    result = validator (arg);
  catch
    result = false;
  end
end
% =========================================================================

author: epfl-lasa


text:
classdef myLoss < dagnn.Loss

  properties
    lossType = 'CE'
  end

  methods
    function outputs = forward(obj, inputs, params)
        X = inputs{1};
        c = inputs{2};
        c = reshape(c,size(X));
        switch obj.lossType
            case 'CE'
                X = vl_nnsoftmax(X);
                Y = squeeze(- c .* log(X./(c  + eps(1))));
            case 'MSE'
                Y = squeeze((X - c).^2);
        end
        outputs{1} = sum(Y(:));
        n = obj.numAveraged ;
        m = n + size(inputs{1},4);
        obj.average = (n * obj.average + gather(outputs{1})) / m ;
        obj.numAveraged = m ;
    end

    function [derInputs, derParams] = backward(obj, inputs, params, derOutputs)
        X = inputs{1};
        c = inputs{2};
        c = reshape(c,size(X));
        switch obj.lossType
            case 'CE'
                X = vl_nnsoftmax(X);
                Y = X - c;
            case 'MSE'
                Y = X - c;
        end
        derInputs = {Y, []};
        derParams = {};
    end

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

author: HuiZeng


text:
%SETFEATSIZE Set featsize of datafile
%
%   A = SETFEATSIZE(A,FEATSIZE)
%
% INPUT
%   A        - Datafile
%   FEATSIZE - Desired feature size,
%              Default: compute from first object in A
%
% OUTPUT
%   A       - Datafile
%
% DESCRIPTION
% If A is a dataset, its featsize is set. The featsize of a datafile
% is unclear as objects may have different sizes. For that reason
% the feature sizes of datafiles might better be set to 0, unless the
% objects are images of the same size.
%
% SEE ALSO (PRTools Guide)
% DATAFILES

author: marianux

repo: ecg-kit

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


text:
function A = group_sum_matrix(G,k)
  % GROUP_SUM_MATRIX Builds a matrix A such that A*V computes the sum of
  % vertices in each group specified by G
  %
  % A = group_sum_matrix(G,k)
  %
  % Inputs:
  %   G  #V list of group indices (1 to k) for each vertex, such that vertex i
  %     is assigned to group G(i)
  %   k  #groups, default is max(G)
  % Outputs:
  %   A  #groups by #V sparse matrix such that A*V = group_sums
  %
  % See also: centroid_matrix

  % number of vertices
  n = numel(G);
  % number of groups
  if ~exist('k','var')
    k = max(G);
  end

  indices =  1:n;

  % builds A such that A(i,j) = 1 where i corresponds to group i and j
  % corresponds to vertex j
  A = sparse(G,indices,1,k,n);

end

author: alecjacobson

repo: gptoolbox


text:
function varargout = ctranspose(varargin)
% Transposing not allowed
% _______________________________________________________________________
% Copyright (C) 2008 Wellcome Trust Centre for Neuroimaging

%
% $Id$

error('file_array objects can not be transposed.');

author: fieldtrip

repo: fieldtrip

save_path: github-repos/MATLAB/fieldtrip-fieldtrip


text:
function FPt(option, up)
%FPT detects Fiducial Points from PPG peak and trough annotations
% as specified in PC's literature review.
%	            FPt(option, up)
%
%	Inputs:
%		option          the option which has led to this function being used
%       up              universal parameters structure
%
%	Outputs:
%       ...
%

fprintf('\n--- Identifying Fiducial Points ');
log_int_respSig = 1;             % Has value 1 unless this is a final respiratory signal

for subj = up.paramSet.subj_list

    %% Skip if this processing has been done previously
    iden_resp_sig_file_ending
    savepath = [up.paths.data_save_folder, num2str(subj), ending];
    filecontents = whos('-file', savepath);
    var_names = extractfield(filecontents, 'name');
    rel_log = zeros(size(var_names));
    if strcmp(option(1:3), 'ekg')
        rel_string = [option(1:3), up.paths.filenames.qrss];
    else
        rel_string = [option(1:3), up.paths.filenames.pulse_peaks];
    end
    for s = 1 : length(var_names)
        if strfind(var_names{s}, rel_string)
            rel_log(s) = 1;
        end
    end
    rel_var_names = var_names(logical(rel_log));
    for rel_var_name_no = 1 : length(rel_var_names)
        eval(['save_name = ''' option(1:3), up.paths.filenames.fid_pts, rel_var_names{rel_var_name_no}(4:end) ''';']);
        exist_log = check_exists(savepath, save_name);
        if exist_log
            continue
        end

        %% Load relevant data
        loadpath = [up.paths.data_save_folder, num2str(subj), up.paths.filenames.int_respSigs];

        if strcmp(option(1:3), 'ppg')
            rel_name1 = rel_var_names{rel_var_name_no};
            rel_name2 = [option(1:3) up.paths.filenames.elim_vhf];
            load(loadpath, rel_name1, rel_name2);
            eval(['rel_data.beats = ' rel_name1 ';']);
            eval(['rel_data.s = ' rel_name2 ';']);
            eval(['rel_data.fs = ' rel_name1 '.fs;']);
            eval(['rel_data.timings = ' rel_name1 '.timings;']);
        else
            %% ECG
            rel_name1 = rel_var_names{rel_var_name_no};
            rel_name2 = [option(1:3) up.paths.filenames.elim_vhf];
            load(loadpath, rel_name1, rel_name2);
            eval(['rel_data.beats = ' rel_name1 ';']);
            eval(['rel_data.s = ' rel_name2 ';']);
            eval(['rel_data.fs = ' rel_name1 '.fs;']);
            eval(['rel_data.timings = ' rel_name1 '.timings;']);
        end

        %% Identify Fiducial Points

        if strcmp(option(1:3), 'ppg')
            %% PPG Peaks
            % Find peaks as max between detected onsets
            temp.p_max.t = nan(length(rel_data.beats.tr.t)-1,1);
            temp.p_max.v = nan(length(rel_data.beats.tr.t)-1,1);
            for beat_no = 1 : (length(rel_data.beats.tr.t)-1)
                rel_range = find(rel_data.s.t >= rel_data.beats.tr.t(beat_no) & rel_data.s.t < rel_data.beats.tr.t(beat_no+1));
                [~, rel_el] = max(rel_data.s.v(rel_range));
                temp.p_max.t(beat_no) = rel_data.s.t(rel_range(rel_el));
                temp.p_max.v(beat_no) = rel_data.s.v(rel_range(rel_el));
            end

            %% PPG Troughs

            % Find troughs as min between detected peaks
            temp.tr_min.t = nan(length(temp.p_max.t)-1,1);
            temp.tr_min.v = nan(length(temp.p_max.t)-1,1);
            for beat_no = 1 : (length(temp.p_max.t)-1)
                rel_range = find(rel_data.s.t >= temp.p_max.t(beat_no) & rel_data.s.t < temp.p_max.t(beat_no+1));
                [~, rel_el] = min(rel_data.s.v(rel_range));
                temp.tr_min.t(beat_no) = rel_data.s.t(rel_range(rel_el));
                temp.tr_min.v(beat_no) = rel_data.s.v(rel_range(rel_el));
            end

        else
            %% ECG Peaks
            % Used to find peaks as max between detected onsets, but it now
            % appears that this doesn't work if it detects the onsets as a
            % mixture of Q and S waves.

            temp.p_max = rel_data.beats.p;

            %% ECG Troughs
            % Find troughs as min within search range 0.1s before peaks
            temp.tr_min.t = nan(length(temp.p_max.t),1);
            temp.tr_min.v = nan(length(temp.p_max.t),1);
            thresh = 0.1;
            for beat_no = 1 : (length(temp.p_max.t))
                rel_range = find(rel_data.s.t >= (temp.p_max.t(beat_no) - thresh) & rel_data.s.t < temp.p_max.t(beat_no));
                % used to be:     rel_range = find(rel_data.s.t >= (temp.p_max.t(beat_no) - thresh) & rel_data.s.t < (temp.p_max.t(beat_no) + thresh));
                [~, rel_el] = min(rel_data.s.v(rel_range));
                if ~isempty(rel_el)   % it is empty if the peak is at the very first element of the signal
                    temp.tr_min.t(beat_no) = rel_data.s.t(rel_range(rel_el));
                    temp.tr_min.v(beat_no) = rel_data.s.v(rel_range(rel_el));
                end
            end
            % get rid of any nans (arise if the peak is at the very first element of the signal)
            bad_els = isnan(temp.tr_min.t);
            temp.tr_min.t = temp.tr_min.t(~bad_els);
            temp.tr_min.v = temp.tr_min.v(~bad_els);

            % very ocassionally it picks out the same trough or peak twice (if two consecutive peaks are ridiculously close together so share some of the same search range)
            [~, rel_els, ~] = unique(temp.tr_min.t);
            temp.tr_min.t = temp.tr_min.t(rel_els);
            temp.tr_min.v = temp.tr_min.v(rel_els);

        end

        % Carry forward detected peaks and onsets
        temp.det_p.t = rel_data.beats.p.t;
        temp.det_p.v = rel_data.beats.p.v;
        temp.det_tr.t = rel_data.beats.tr.t;
        temp.det_tr.v = rel_data.beats.tr.v;

        % carry forward fs and timings
        temp.fs = rel_data.fs;
        temp.timings = rel_data.timings;

        %% Save processed data
        eval([save_name ' = temp;']);
        save_or_append_data

    end

end

end


text:
function plot_path(map, path)
% PLOT_PATH Visualize a path through an environment
%   PLOT_PATH(map, path) creates a figure showing a path through the
%   environment.  path is an N-by-3 matrix where each row corresponds to the
%   (x, y, z) coordinates of one point along the path.

path = points_to_idx(map, path);

[x,y,z] = meshgrid(1:size(map.map3d,1),1:size(map.map3d,2),1:size(map.map3d,3));
x = x(:); y = y(:); z = z(:);

xyz = [];
for i = 1:numel(x)
    if map.map3d(x(i),y(i),z(i),1)~= 255 | map.map3d(x(i),y(i),z(i),2)~= 255 | map.map3d(x(i),y(i),z(i),2)~= 255
        xyz = [xyz; x(i),y(i),z(i)];
    end
end
rgb = zeros(size(xyz,1),3);
for i = 1:size(xyz,1)
    rgb(i,:) = map.map3d(xyz(i,1),xyz(i,2),xyz(i,3),:);
end
if size(xyz,1) > 0
pcshow(xyz, rgb, 'MarkerSize', 200);
end
hold on;
if size(path,1) > 0
pcshow(path, [0,0,0],'MarkerSize', 20);
end
hold off;
end

author: yrlu

repo: quadrotor

save_path: github-repos/MATLAB/yrlu-quadrotor


text:
classdef TimeDepthPlotter < TimeSomethingPlotter
    % Used to create time-depth plots

    methods
        function obj=TimeDepthPlotter()
            obj@TimeSomethingPlotter('Depth','km');
        end

        function pl=plot(obj,ax,catalog,bigcat)
            % plot plot a time-depth series for this catalog, with symbol sizes representing
            % event size
            % pl = plot(catalog)
            %

            if ~exist('ax','var') || isempty(ax)
                f=figure('Name','Time Depth','NumberTitle','off', ...
                    'Tag',obj.Tags.Figure);
                addAboutMenuItem();
                ax=axes(f);
            end
            obj.ax=ax;

            if isempty(ax.Tag)
                ax.Tag=obj.Tags.Axes;
            end

            ax.Visible = 'off';

            if obj.hasLotsOfEvents(catalog)
                obj.Marker='.';
            else
                obj.Marker='s';
            end
            pl=obj.scatter(catalog);

            obj.prepare_axes(catalog,'YDir','reverse');

            f=ancestor(ax,'figure');
            c=findobj(f,'Tag','TimeDepthContext');
            if isempty(c)
                c=uicontextmenu('Tag','TimeDepthContext');
                uimenu(c,'Label','Use Log Scale','MenuSelectedFcn',@(s,~)logtoggle(s,'Y'));
            end
            ax.YLabel.UIContextMenu=c;

            if exist('bigcat','var')
                obj.overlayBigEvents(bigcat);
            end
            ax.Visible = 'on';

        end

    end
end

author: CelsoReyes


text:
function SO3F = dot(SO3VF1, SO3VF2, varargin)
% pointwise inner product
%
% Syntax
%   SO3F = dot(SO3VF1, SO3VF2)
%   SO3F = dot(v, SO3VF2)
%   SO3F = dot(SO3VF1, v)
%
% Input
%   SO3VF1, SO3VF2 - @SO3VectorField
%   v - @vector3d
%
% Output
%   SO3F - @SO3Fun
%

if isa(SO3VF2, 'vector3d')
  SO3F = SO3VF1.SO3F;
  v = SO3VF2;
  SO3F = reshape(sum( SO3F .* reshape(v.xyz.',[3,size(v)]),1),size(v));
  if SO3VF2.antipodal
    SO3F = abs(SO3F);
  end
  return
end

if isa(SO3VF1, 'vector3d')
  SO3F = dot(SO3VF2,SO3VF1);
  return
end

ensureCompatibleSymmetries(SO3VF1,SO3VF2)
f = SO3FunHandle(@(rot) dot(SO3VF1.eval(rot),SO3VF2.eval(rot)),SO3VF1.CS,SO3VF1.SS);
SO3F = SO3FunHarmonic(f, varargin{:});

end

author: mtex-toolbox

repo: mtex

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


text:
function test_suite=test_mask_dim_intersect
% tests for test_mask_dim_intersect
%
% #   For CoSMoMVPA's copyright information and license terms,   #
% #   see the COPYING file distributed with CoSMoMVPA.           #
    try % assignment of 'localfunctions' is necessary in Matlab >= 2016
        test_functions=localfunctions();
    catch % no problem; early Matlab versions can use initTestSuite fine
    end
    initTestSuite;

function test_mask_dim_intersect_feature_dim()
    helper_test_mask_dim_intersect(2)

function test_mask_dim_intersect_sample_dim()
    helper_test_mask_dim_intersect(1)

function helper_test_mask_dim_intersect(dim)
    other_dim=3-dim;

    n_ds=ceil(rand()*2+2);
    keep_ratio=.8;

    ds_cell=cell(n_ds,1);
    for k=1:n_ds
        if dim==1
            % optimalization due to somewhat slow dim_transpose
            sz='small';
        else
            sz='big';
        end

        ds=cosmo_synthetic_dataset('size',sz,'seed',0,...
                                    'nchunks',1,'ntargets',1);
        n_features=size(ds.samples,2);
        ds.samples=1:n_features;

        rp=randperm(n_features);
        keep_indices=round(1:(keep_ratio*n_features));

        ds_keep=cosmo_slice(ds,rp(keep_indices),2);

        if dim==1
            ds_keep=cosmo_dim_transpose(ds_keep,ds.a.fdim.labels,1);
        end

        ds_indices=ds_keep;

        % fill with random data so that the function to test cannot
        % use the indices
        ds_data=ds_keep;
        ds_data.samples=randn(size(ds_keep.samples));

        % the first row has random data, and is used as input for %
        % cosmo_mask_dim_intersect; the second row has indices, and is
        % used to verify the proper contents of the datasets
        ds_cell{k}=cosmo_stack({ds_data,ds_indices},other_dim);
    end

    ds_data_cell=cellfun(@(x)cosmo_slice(x,1,other_dim),ds_cell,...
                            'UniformOutput',false);

    [indices,ds_intersect_cell]=cosmo_mask_dim_intersect(ds_data_cell,dim);

    if dim==2
        % must have same output using default second argument
        [indices2,ds_intersect_cell2]=cosmo_mask_dim_intersect(...
                                            ds_data_cell);

        assertEqual(indices,indices2);
        assertEqual(ds_intersect_cell,ds_intersect_cell2);
    end

    % verify ds_intersect_cell based on indices
    assert(iscell(ds_intersect_cell));
    assertEqual(size(ds_intersect_cell),size(ds_data_cell));
    for k=1:n_ds
        idx=indices{k};
        assert(all(idx>0));
        assert(all(isfinite(idx)));
        assertEqual(cosmo_slice(ds_data_cell{k},idx,dim),...
                                ds_intersect_cell{k})
    end

    ds_indices_cell=cellfun(@(x)cosmo_slice(x,2,other_dim),ds_cell,...
                            'UniformOutput',false);

    % see which indices are common across all datasets
    ds_keep_indices=1:n_features;
    for k=1:n_ds
        ds_idx=ds_indices_cell{k}.samples;
        assertEqual(sort(ds_idx),unique(ds_idx));
        ds_keep_indices=intersect(ds_keep_indices, ds_idx);
    end

    % verify that indices match across all datasets
    for k=1:n_ds
        % select indices
        idx=ds_indices_cell{k}.samples(indices{k});

        % indices must be unique
        if isempty(idx)
            assert(isempty(ds_keep_indices));
        else
            assertEqual(sort(idx(:)),sort(ds_keep_indices(:)));
        end

        % must return indices in the same order
        ds_sel=cosmo_slice(ds_indices_cell{k},indices{k},dim);
        if k==1
            first_ds_sel=ds_sel;
        else
            assertEqual(first_ds_sel,ds_sel)
        end
    end

function test_mask_dim_intersect_identity
    % after permuting the features, dataset should be identical following
    % unpermuting them

    types={'fmri','surface','source','timelock','timefreq'};

    for k=1:numel(types)
        type=types{k};
        ds=cosmo_synthetic_dataset('size','big','type',type,...
                                        'ntargets',1,'nchunks',1);

        switch type
            case 'source'
                args={'matrix_labels',{'pos'}};
            otherwise
                args={};
        end
        n_features=size(ds.samples,2);
        ds.samples(1,:)=1:n_features;

        rp=randperm(n_features);
        ds=cosmo_slice(ds,rp,2);

        [indices_cell,ds_cell]=cosmo_mask_dim_intersect({ds},2,args{:});
        assert(numel(indices_cell)==1);
        assert(numel(ds_cell)==1);

        idx=indices_cell{1};
        ds_perm=ds_cell{1};

        assertEqual(sort(idx),1:n_features);
        ds_reordered=cosmo_slice(ds,idx,2);
        assertEqual(ds_perm,ds_reordered);

        % indices must be sorted
        assertEqual(ds_perm.samples,1:n_features);

    end

function test_mask_dim_intersect_exceptions()
    aet=@(varargin)assertExceptionThrown(@()...
                        cosmo_mask_dim_intersect(varargin{:}),'');
    % cannot deal with empty
    aet(cosmo_synthetic_dataset());

    % input must be a cell with datasets
    aet(struct);
    aet('foo');
    aet({struct,struct});

    % repeated features not supported
    ds=cosmo_synthetic_dataset('size','big');
    ds2=cosmo_stack({ds,ds},2);
    aet({ds2});

    % not even a single feature duplicated is allowed
    col_index=ceil(rand()*size(ds.samples,2));
    ds_extra_col=cosmo_stack({ds,cosmo_slice(ds,col_index,2)},2);
    aet(ds_extra_col);

    % dim must be 1 or 2
    aet({ds},3);
    aet({ds},struct);
    aet({ds},[1 2]);

function test_mask_dim_intersect_missing_dim()
    ds1=cosmo_synthetic_dataset();
    ds2=ds1;
    ds2.fa=rmfield(ds2.fa,'k');

    assertEqual(ds1.samples,ds2.samples);
    assertExceptionThrown(@()cosmo_mask_dim_intersect({ds1,ds2}),'');

function test_mask_dim_intersect_nonmatching_dim()
    ds1=cosmo_synthetic_dataset();
    ds2=ds1;
    ds2=cosmo_dim_remove(ds2,'k');

    assertEqual(ds1.samples,ds2.samples);
    assertExceptionThrown(@()cosmo_mask_dim_intersect({ds1,ds2}),'');

function test_mask_dim_intersect_renamed_dim()
    ds1=cosmo_synthetic_dataset();
    ds2=ds1;
    ds2=cosmo_dim_rename(ds2,'k','kk');

    assertEqual(ds1.samples,ds2.samples);
    assertExceptionThrown(@()cosmo_mask_dim_intersect({ds1,ds2}),'');

author: CoSMoMVPA

repo: CoSMoMVPA

save_path: github-repos/MATLAB/CoSMoMVPA-CoSMoMVPA


text:
function [D,W,Lw,Lb]=BuildPatch(fea,gnd,k1,k2,zeta)

% Build two patches for MD-NMF and NDLA by using KNN + sparse binary model
%   data: nSample x nFeature matrix,
%   zeta: perturbution coefficient.
% Written by Naiyang Guan (ny.guan@gmail.com)
% Copyright by Naiyang Guan and Dacheng Tao

% Initialization: assume data is MxN matrix
N=size(fea, 1);
Sw=zeros(N);
Sb=zeros(N);

% Extreme position
if k2<=0
    fprintf('error! k2 cannot be smaller than 1\n');
    return;
end

% Euclidean distance between sample points in the PCA space
dist=EuDist2(fea,fea,0);

% Construct sparse Laplacian matrix on (k1,k2)-patch
for i=1:N
    % k1 nearest neighbors within same class
    sam_idx=find(gnd==gnd(i));
    sam_dist=dist(i,sam_idx');
    [junk,index]=sort(sam_dist);
    sam_idx=sam_idx(index);
    if k1>=1
        Fi=sam_idx(2:(k1+1));
    else
        Fi=i;
    end
    Sw(i,Fi)=1;
    Sw(Fi,i)=1;

    % k2 nearest neighbors between classes
    diff_idx=find(gnd~=gnd(i));
    diff_dist=dist(i,diff_idx');
    [junk,index]=sort(diff_dist);
    diff_idx=diff_idx(index);
    Fi=diff_idx(1:k2);
    Sb(i,Fi)=1;
    Sb(Fi,i)=1;
end

% Laplacian regularization and normalization
Dw=diag(sum(Sw));
Db=diag(sum(Sb));
Lw=Dw-Sw;
Lb=Db-Sb+eye(N)*sum(sum(Sb))*zeta;

[U,S,V]=svd(Lb);
S=diag(S);
S=diag(1./sqrt(S));
D=U*S*V'*Dw*V*S*U';
W=U*S*V'*Sw*V*S*U';

return;

author: andrewssobral

repo: lrslibrary


text:
function calpak_test002 ( )

%*****************************************************************************80
%
%% CALPAK_TEST002 tests EASTER_JULIAN and EASTER_JULIAN2.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    19 June 2012
%
%  Author:
%
%    John Burkardt
%
  n_test = 10;
  d_test = [  27,    19,   11,   30,   15,    5,   27,   11,    1,   23 ];
  m_test = [   4,     4,    4,    4,    4,    5,    4,    4,    5,    4 ];
  y_test = [ 1997, 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006 ];

  fprintf ( 1, '\n' );
  fprintf ( 1, 'CALPAK_TEST002\n' );
  fprintf ( 1, '  For the Julian calendar,\n' );
  fprintf ( 1, '  for a given year, compute the day and month of Easter.\n' );
  fprintf ( 1, '  EASTER_JULIAN uses Richard''s algorithm.\n' );
  fprintf ( 1, '  EASTER_JULIAN2 uses Richards''s algorithm.\n' );

  for i = 1 : n_test

    y = y_test(i);
    m = m_test(i);
    d = d_test(i);
    f = 0.5;

    fprintf ( 1, '\n' );
    s = ymd_to_s_gregorian ( y, m, d );
    fprintf ( 1, '  CORRECT (Gregorian): %s\n', s )

    jed = ymdf_to_jed_gregorian ( y, m, d, f );
    [ y, m, d, f ] = jed_to_ymdf_julian ( jed );

    s = ymd_to_s_julian ( y, m, d );
    fprintf ( 1, '  CORRECT (Julian):    %s\n', s )

    [ m, d ] = easter_julian ( y );
    s = ymd_to_s_julian ( y, m, d );
    fprintf ( 1, '  EASTER_JULIAN:       %s\n', s )

    [ m, d ] = easter_julian2 ( y );
    s = ymd_to_s_julian ( y, m, d );
    fprintf ( 1, '  EASTER_JULIAN2:      %s\n', s )

  end

  return
end

author: johannesgerer

repo: jburkardt-m


text:
% function [enc]=jp2write(img,filename,indexing,arguments)
%
% It takes a matlab image and stores it in a JPEG2000 compliant file.
% Important notice: the JVM matlab uses is important in order to
% see if jars provided in classpath are compatible or not. If not,
% the only message thrown is 'class not found'.
% Here I've built the jars using -source 1.4 -target 1.4 options in javac
% while it is possible to change the JVM that Matlab is running
% by setting the MATLAB_JAVA system variable. In order to know which
% version of JVM Matlab is running the command is version -java
%
% Using Matlab 7.0.1.24704 (R14) Service Pack 1
% the corresponding JVM is
% Java 1.4.2_05 with Sun Microsystems Inc. Java HotSpot(TM) Client VM
%    (mixed mode)
%
% To resemble the use of the standard command line environment, the
% arguments variable is intented to be a big string with all the
% jj2000 parameters passes as if they were to be driving the command
% line implementation of the encoder.

function enc=jp2write(img,filename,indexing,arguments)

    import middleware.*;
    import jj2000.j2k.util.*;
    import jj2000.j2k.encoder.*;

    % la sistemo sopra l'immagine in modo che qui arrivi sempre RGB

    [image]=im2java(img);

    %imgr=ImgReaderMTLB(image);

    defpl=ParameterList();

    params=Encoder.getAllParameters();

    for i=1:size(params,1)
        if ne(params(i,4),'')
            defpl.setProperty(params(i,1),params(i,4));
            disp(sprintf('%s=%s',char(params(i,1)),char(params(i,4))));
        end
    end

    pl=ParameterList(defpl);

    pl.setProperty('v','on');
    pl.setProperty('debug','on');
    pl.setProperty('Qguard_bits','2'); % default values
    pl.setProperty('file_format','on'); % this to use XML
    pl.setProperty('o',filename); % to resemble matlab

    na=1;
    [t(na).s,r]=strtok(arguments);
    while(~isempty(r))
        na=na+1;
        [t(na).s,r]=strtok(r)
    end
    argv=javaArray('java.lang.String',na);
    for i=1:na
        argv(i)=java.lang.String(t(i).s);

    end
% [PENDING] Gestione degli spazi nei nomi dei file delle ROI
% Serve agire sullo strato di middleware
%     argv(na)=java.lang.String(t(na).s);
%     argv(na)=javaMethod('replaceAll',argv(na),'\*','\\ ')

    pl.parseArgs(argv);

    enc=MyEncoder(ParameterList(pl),image,indexing);

    javaMethod('run',enc);

return

author: cerr

repo: CERR

save_path: github-repos/MATLAB/cerr-CERR