text:
%ZMQ socket communication via ZMQ
%
% ZMQ is a wrapper for libzmq, which exposes a simple ZMQ REP
% endpoint. You can receive messages from the socket, and send messages
% to the socket.
%
% ZMQ Methods:
%   receive - read a message from the socket
%   send    - write a message to the socket

% Copyright (c) 2016, Bastian Bechtold
% This code is published under the terms of the BSD 3-clause license

classdef ZMQ < handle
    properties (Access=private)
        context
        socket
    end
    methods
        function obj = ZMQ(libname, address)
            ZMQ_REP = 4;
            try
                if not(libisloaded('libzmq'))
                    [notfound, warnings] = ...
                        loadlibrary(libname, 'transplantzmq.h', ...
                                    'alias', 'libzmq');
                    % the library did not contain the functions we need:
                    assert(isempty(notfound), 'Could not load ZMQ library')
                end
                obj.context = calllib('libzmq', 'zmq_ctx_new');
                assert(not(obj.context.isNull), ...
                       'zmq_ctx_new failed: Could not create context');
                obj.socket = calllib('libzmq', 'zmq_socket', obj.context, ZMQ_REP);
                assert(not(obj.socket.isNull), obj.errortext('zmq_socket'));
                err = calllib('libzmq', 'zmq_connect', obj.socket, address);
                assert(err == 0, obj.errortext('zmq_connect'));
            catch exception
                % print exception, since we probably don't have a working connection
                % to the transplant master yet for reporting errors properly:
                disp(['Error loading libzmq: ' exception.message]);
                throw(exception);
            end
        end

        function str = receive(obj)
            msg = libstruct('zmq_msg_t', struct('hidden', zeros(1, 64, 'uint8')));
            calllib('libzmq', 'zmq_msg_init', msg); % always returns 0
            msglen = calllib('libzmq', 'zmq_msg_recv', msg, obj.socket, 0);
            assert(msglen >= 0, obj.errortext('zmq_msg_recv'));
            msgptr = calllib('libzmq', 'zmq_msg_data', msg);
            if not(msgptr.isNull)
                setdatatype(msgptr, 'uint8Ptr', 1, msglen);
                str = uint8(msgptr.Value);
            else
                str = uint8([]);
            end
            err = calllib('libzmq', 'zmq_msg_close', msg);
            assert(err == 0, obj.errortext('zmq_msg_close'));
        end

        function send(obj, data)
            dataptr = libpointer('uint8Ptr', data);
            msglen = calllib('libzmq', 'zmq_send', obj.socket, dataptr, numel(data), 0);
            assert(msglen >= 0, obj.errortext('zmq_send'));
        end

        function delete(obj)
            % if we crashed in the constructor:
            if ~libisloaded('libzmq')
                return
            end
            err = calllib('libzmq', 'zmq_close', obj.socket);
            assert(err == 0, obj.errortext('zmq_close'));
            err = calllib('libzmq', 'zmq_ctx_term', obj.context);
            assert(err == 0, obj.errortext('zmq_ctx_term'));
            unloadlibrary('libzmq');
        end
    end

    methods (Hidden=true)
        function str = errortext(obj, instruction)
            base = 156384712;
            errno = calllib('libzmq', 'zmq_errno');
            switch errno
                case 4 % EINTR
                    if strcmp(instruction, 'zmq_ctx_term')
                        str = 'Termination was interrupted by a signal. It can be restarted if needed.';
                    elseif strcmp(instruction, 'zmq_send')
                        str = 'The operation was interrupted by delivery of a signal before the message was sent.';
                    elseif strcmp(instruction, 'zmq_msg_recv')
                        str = 'The operation was interrupted by delivery of a signal before a message was available.';
                    else
                        str = 'Interrupted system call';
                    end
                case 11 % EAGAIN
                    if strcmp(instruction, 'zmq_send')
                        str = 'Non-blocking mode was requested and the message cannot be sent at the moment.';
                    elseif strcmp(instruction, 'zmq_msg_recv')
                        str = 'Non-blocking mode was requested and no message are available at the moment.';
                    else
                        str = 'Try again';
                    end
                case 22 % EINVAL
                    if strcmp(instruction, 'zmq_socket')
                        str = 'The requested socket *type* is invalid.';
                    elseif strcmp(instruction, 'zmq_connect')
                        str = 'The endpoint supplied is invalid.';
                    else
                        str = 'Invalid argument';
                    end
                case 14 % EFAULT
                    if strcmp(instruction, 'zmq_socket') || strcmp(instruction, 'zmq_ctx_term')
                        str = 'The provided *context* is invalid.';
                    elseif strcmp(instruction, 'zmq_msg_close')
                        str = 'Invalid message';
                    elseif strcmp(instruction, 'zmq_msg_recv')
                        str = 'The message passed to the function was invalid.';
                    else
                        str = 'Invalid argument';
                    end
                case 24 % EMFILE
                    if strcmp(instruction, 'zmq_socket')
                        str = 'The limit on the total number of open ZMQ sockets has been reached.';
                    else
                        str = 'Too many open files';
                    end
                case base + 1 % ENOTSUP
                    str = 'The *zmq_send()* operation is not supported by this socket type.';
                case base + 2 % EPROTONOSUPPORT
                    str = 'The requested *transport* protocol is not supported.';
                case base + 9 % ENOTSOCK
                    str = 'The provided *socket* was invalid.';
                case base + 17 % EHOSTUNREACH
                    str = 'The message cannot be routed.';
                case base + 51 % EFSM
                    if strcmp(instruction, 'zmq_send')
                        str = 'The *zmq_send()* operation cannot be performed on this socket at the moment due to the socket not being in the appropriate state.';
                    else
                        str = 'socket not in appropriate state.';
                    end
                case base + 52 % ENOCOMPATPROTO
                    str = 'The requested *transport* protocol is not compatible with the socket type.';
                case base + 53 % ETERM
                    str = 'The ZMQ *context* associated with the specified *socket* was terminated.';
                case base + 54 % EMTHREAD
                    str = 'No I/O thread is available to accomplish the task.';
                otherwise
                    str = sprintf('Unknown error %d', errno);
            end
            str = sprintf('Error in %s: %s', instruction, str);
        end
    end
end

author: bastibe

repo: transplant

save_path: github-repos/MATLAB/bastibe-transplant


text:
function varargout = KUKA_GUI(varargin)
% Begin initialization code - DO NOT EDIT
gui_Singleton = 1;
gui_State = struct('gui_Name',       mfilename, ...
                   'gui_Singleton',  gui_Singleton, ...
                   'gui_OpeningFcn', @KUKA_GUI_OpeningFcn, ...
                   'gui_OutputFcn',  @KUKA_GUI_OutputFcn, ...
                   'gui_LayoutFcn',  [] , ...
                   'gui_Callback',   []);
if nargin && ischar(varargin{1})
    gui_State.gui_Callback = str2func(varargin{1});
end

if nargout
    [varargout{1:nargout}] = gui_mainfcn(gui_State, varargin{:});
else
    gui_mainfcn(gui_State, varargin{:});
end
% End initialization code - DO NOT EDIT

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%                                                                        %%
%%%                  Written by Rohit Mishra                               %%
%%%                  A.C. Patil College of Engg.                                            %%
%%%                  March 2013                                            %%
%%%                  Simulation of six axis kuka kr6.                      %%
%%%                                                                        %%
%%%      The purpose of this program is to simulate a six axis kuka kr6    %%
%%%      robot located in the DRHR DEPT. at BARC MUMBAI .  This            %%
%%%      program simulates Forward Kinematics of this robot.               %%
%%%      To operate this program:                                          %%
%%%           1)     Make sure that KUKA.WRL, KUKA_GUI.fig, KUKA_GUI.m     %%
%%%                  are in the appropriate path.                          %%
%%%           2)     Execute KUKA_GUI.m                                    %%
%%%           4)     Have fun!                                             %%
%%%      If you have question about this program email Rohit Mishra        %%
%%%      <mishra.rohit2410@gmail.com> .  Rohit Mishra is Final Year Student%%
%%%      from Mumbai University(ACPCE, TRONIX 2013)has developed this      %%
%%%      project with the help of his friends Tejas Mahadik, Kaushik       %%
%%%      Natrajan and  Siddhesh Chavan .                                   %%
%%%      We are thankful to our project guide and HOD Dr. V.N. Pawar for   %%
%%%      his unwavering support and guidance in helping us with this       %%
%%%      project and also to Mrs. Varsha Shirwalkar at Barc, Mumbai.       %%
%%%                                                                        %%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function KUKA_GUI_OpeningFcn(hObject, eventdata, handles, varargin)
handles.output = hObject;
guidata(hObject, handles);
KUKA = vrworld('KUKA.wrl');
open(KUKA)
view(KUKA)

    radian=0
    KUKA.shoulder.rotation = [0, 1, 0, radian];

    KUKA.elbow.rotation = [0, 0, 1, radian];

    KUKA.pitch.rotation = [0, 0, 1, radian];

    KUKA.roll.rotation = [1, 0, 0, radian];

    KUKA.yaw.rotation = [0, 0, 1, radian];

    KUKA.drill.rotation = [1, 0, 0, radian];

% --- Executes on slider movement.
function slider1_Callback(hObject, eventdata, handles)
KUKA = vrworld('KUKA.wrl');
open(KUKA)
T1=get(handles.slider1,'value');
set(handles.edit57,'String',num2str(T1));
guidata(hObject, handles);
    radian=T1*pi/180;
    q1=radian;
    KUKA.shoulder.rotation = [0, 1, 0, radian];

% --- Executes on slider movement.
function slider2_Callback(hObject, eventdata, handles)
KUKA = vrworld('KUKA.wrl');
open(KUKA)
T2=get(handles.slider2,'value');
set(handles.edit58,'String',num2str(T2));
guidata(hObject, handles);
    radian=T2*pi/180;
    q2=radian;
    KUKA.elbow.rotation = [0, 0, 1, radian];

% --- Executes on slider movement.
function slider3_Callback(hObject, eventdata, handles)
KUKA = vrworld('KUKA.wrl');
open(KUKA)
T3=get(handles.slider3,'value');
set(handles.edit59,'String',num2str(T3));
guidata(hObject, handles);
    radian=T3*pi/180;
    q3=radian;
    KUKA.pitch.rotation = [0, 0, 1, radian];

% --- Executes on slider movement.
function slider4_Callback(hObject, eventdata, handles)
KUKA = vrworld('KUKA.wrl');
open(KUKA)
T4=-get(handles.slider4,'value');
set(handles.edit60,'String',num2str(T4));
guidata(hObject, handles);
    radian=T4*pi/180;
    q4=radian;
    KUKA.roll.rotation = [1, 0, 0, radian];

% --- Executes on slider movement.
function slider5_Callback(hObject, eventdata, handles)
KUKA = vrworld('KUKA.wrl');
open(KUKA)
T5=-get(handles.slider5,'value');
    set(handles.edit61,'String',num2str(T5));
    guidata(hObject, handles);
    radian=T5*pi/180;
    q5=radian;
    KUKA.yaw.rotation = [0, 0, 1, radian];

    % --- Executes on slider movement.
function slider6_Callback(hObject, eventdata, handles)
KUKA = vrworld('KUKA.wrl');
open(KUKA)
T6=-get(handles.slider6,'value');
set(handles.edit62,'String',num2str(T6));
guidata(hObject, handles);
    radian=T6*pi/180;
    q6=radian;
    KUKA.drill.rotation = [1, 0, 0, radian];

% --- Executes on button press in pushbutton2.
function pushbutton2_Callback(hObject, eventdata, handles)

% % % GIVEN JOINT VARIABLES

d=[335   0   0   -295    0   -80];
a=[75    270 90  0   0   0];

% % TRANSFORMATION MATRIX 1.
q1=get(handles.slider1,'Value');
T1=[    cos(q1)     0       -sin(q1) a(1)*cos(q1)
        sin(q1)     0       -cos(q1) a(1)*sin(q1)
        0           -1          0       d(1)
        0           0           0        1      ];

% % TRANSFORMATION MATRIX 2.
q2=get(handles.slider2,'Value');
T2=[    cos(q2)     -sin(q2)    0   a(2)*cos(q2)
        sin(q2)     cos(q2)     0   a(2)*sin(q2)
        0               0       1       0
        0               0       0       1   ];

% % TRANSFORMATION MATRIX 3.
q3=get(handles.slider3,'Value');
T3=[    cos(q3)         0    sin(q3)   a(3)*cos(q3)
        sin(q3)         0    -cos(q3)   a(3)*sin(q3)
        0               0       0       0
        0               0       0       1   ];
% % TRANSFORMATION MATRIX 4.
q4=get(handles.slider4,'Value');
T4=[    cos(q4)         0   -sin(q4)    0
        sin(q4)         0    cos(q4)    0
        0              -1       0      d(4)
        0               0       0       1   ];

% % TRANSFORMATION MATRIX 5.
q5=get(handles.slider5,'Value');
T5=[    cos(q5)     0       sin(q5)           0
        sin(q5)     0       -cos(q5)          0
        0           1          0               0
        0           0          0               1       ];

% % TRANSFORMATION MATRIX 6.
q6=get(handles.slider6,'Value');
T6=[    cos(q6)     -sin(q6)        0           0
        sin(q6)     cos(q6)     0               0
        0           0           1           0
        0           0           0               1   ];

A=T1*T2*T3*T4*T5*T6;

set(handles.edit42,'String',num2str(A(13)));
set(handles.edit43,'String',num2str(A(14)));
set(handles.edit44,'String',num2str(A(15)));
set(handles.edit48,'String',num2str(A(1)));
set(handles.edit49,'String',num2str(A(5)));
set(handles.edit50,'String',num2str(A(9)));
set(handles.edit51,'String',num2str(A(2)));
set(handles.edit52,'String',num2str(A(6)));
set(handles.edit53,'String',num2str(A(10)));
set(handles.edit54,'String',num2str(A(3)));
set(handles.edit55,'String',num2str(A(7)));
set(handles.edit56,'String',num2str(A(11)));
guidata(hObject, handles);

author: Sable

repo: mcbench-benchmarks


text:
% load data

load initData

% set irradiance

irradiance = StationData(1).Irradiance;

author: mathworks

repo: MathWorks-Excellence-in-Innovation


text:
% Author: Francesco Ledda
%
% The Image Smoothing Filter implements a set of utilities for the digital treatment of images. The main purpose is the image cleaning from noise
% but there are also other utilities like noise addition, conversion in gray-scale, border emphasis, bidimensional FFT calculation and so on.
%
% This file is just for a brief description of the project and it aims to be the driver of the actual GUI.
%
% For more information about the code, you can send me an e-mail at leddafrancesco89@gmail.com
%
% If you use this code, please don't forget to mention the author. :)
%

% Start the GUI

Progetto


text:
function dt6 = dtiLoadTensorElements(teBasename)
% dt6 = dtiLoadTensorElements(teBasename)
%
%
% HISTORY:
% 2004.03.26 RFD (bobd@stanford.edu) Wrote it.
%
% TODO:
%
%

if(~exist('teBasename', 'var') | isempty(teBasename))
    [f,p] = uigetfile({'*.001';'*.*'},'Select a slice from the Tensor Elements file...');
    teBasename = fullfile(p,f(1:end-3));
end

if(teBasename(end) ~= '.')
    teBasename = [teBasename '.'];
end

nSlices = length(dir([teBasename,'*']));
d = dir([teBasename,'001']);
% OK- there are 6*4=24 bytes per voxel -> d.bytes/24 voxels. We assume that
% the bits are divided equally between X and Y.
xySize = sqrt(d.bytes./24);
imdim = [xySize,xySize,nSlices];
disp([mfilename ': I''m guessing that the image dimension is [',num2str(imdim),'] (based on file size).']);

npix = imdim(1)*imdim(2);
dt6 = zeros([imdim,6]);
disp(['Reading ',num2str(nSlices),' slices...']);
for(ii=1:imdim(3))
    fname = sprintf('%s%03d', teBasename, ii);
    fid = fopen(fname, 'rb', 'ieee-le');
    dt6(:,:,ii,:) = reshape(fread(fid, inf, 'float32'),[imdim(1),imdim(2),1,6]);
    fclose(fid);
end
return;

author: vistalab

repo: vistasoft

save_path: github-repos/MATLAB/vistalab-vistasoft


text:
% From examples/imagenet/cnn_imagenet_deploy.m
% -------------------------------------------------------------------------
function dagMergeBatchNorm(net, names)
% -------------------------------------------------------------------------
for name = names
  name = char(name) ;
  layer = net.layers(net.getLayerIndex(name)) ;

  % merge into previous conv layer
  playerName = dagFindLayersWithOutput(net, layer.inputs{1}) ;
  playerName = playerName{1} ;
  playerIndex = net.getLayerIndex(playerName) ;
  player = net.layers(playerIndex) ;
  if ~isa(player.block, 'dagnn.Conv')
    error('Batch normalization cannot be merged as it is not preceded by a conv layer.') ;
  end

  % if the convolution layer does not have a bias,
  % recreate it to have one
  if ~player.block.hasBias
    block = player.block ;
    block.hasBias = true ;
    net.renameLayer(playerName, 'tmp') ;
    net.addLayer(playerName, ...
                 block, ...
                 player.inputs, ...
                 player.outputs, ...
                 {player.params{1}, sprintf('%s_b',playerName)}) ;
    net.removeLayer('tmp') ;
    playerIndex = net.getLayerIndex(playerName) ;
    player = net.layers(playerIndex) ;
    biases = net.getParamIndex(player.params{2}) ;
    net.params(biases).value = zeros(block.size(4), 1, 'single') ;
  end

  filters = net.getParamIndex(player.params{1}) ;
  biases = net.getParamIndex(player.params{2}) ;
  multipliers = net.getParamIndex(layer.params{1}) ;
  offsets = net.getParamIndex(layer.params{2}) ;
  moments = net.getParamIndex(layer.params{3}) ;

  [filtersValue, biasesValue] = mergeBatchNorm(...
    net.params(filters).value, ...
    net.params(biases).value, ...
    net.params(multipliers).value, ...
    net.params(offsets).value, ...
    net.params(moments).value) ;

  net.params(filters).value = filtersValue ;
  net.params(biases).value = biasesValue ;
end

% -------------------------------------------------------------------------
function [filters, biases] = mergeBatchNorm(filters, biases, multipliers, offsets, moments)
% -------------------------------------------------------------------------
% wk / sqrt(sigmak^2 + eps)
% bk - wk muk / sqrt(sigmak^2 + eps)
a = multipliers(:) ./ moments(:,2) ;
b = offsets(:) - moments(:,1) .* a ;
biases(:) = biases(:) + b(:) ;
sz = size(filters) ;
numFilters = sz(4) ;
filters = reshape(bsxfun(@times, reshape(filters, [], numFilters), a'), sz) ;

% -------------------------------------------------------------------------
function layers = dagFindLayersWithOutput(net, outVarName)
% -------------------------------------------------------------------------
layers = {} ;
for l = 1:numel(net.layers)
  if any(strcmp(net.layers(l).outputs, outVarName))
    layers{1,end+1} = net.layers(l).name ;
  end
end

author: bertinetto

repo: siamese-fc


text:
classdef tree2

    properties
        % parent pointers
        pp = [];
        nodeNames;
        nodeFeatures_unnormalized =[];
        nodeFeatures = [];
        leafFeatures=[];
        % the parent pointers do not save which is the left and right child of each node, hence:
        % numNodes x 2 matrix of kids, [0 0] for leaf nodes
        kids = [];
        numkids = [];
        freq = [];
        nodeLabels=[];
        score=0;
        nodeScores=[];
        nodeZ = [];
        nodeDelta_out1 = [];
        nodeDelta_out2 = [];
        parentDelta = [];
        catDelta = [];
        catDelta_out = [];
        node_dW1 = [];
        node_dW2 = [];
        node_dW3 = [];
        node_dW4 = [];
        node_dL = [];
        node_y1c1 = [];
        node_y2c2 = [];

    end

    methods
        function id = getTopNode(obj)
            id = find(obj.pp==0);
        end

        function kids = getKids(obj,node)
            %kids = find(obj.pp==node);
            kids = obj.kids(node,:);
        end

        %TODO: maybe compute leaf-node-ness once and then just check for it
        function l = isLeaf(obj,node)
            l = ~any(obj.pp==node);
        end

function plotTree(obj)
            %TREEPLOT Plot picture of tree.
            %   TREEPLOT(p) plots a picture of a tree given a row vector of
            %   parent pointers, with p(i) == 0 for a root and labels on each node.
            %
            %   Example:
            %      myTreeplot([2 4 2 0 6 4 6],{'i' 'like' 'labels' 'on' 'pretty' 'trees' '.'})
            %   returns a binary tree with labels.
            %
            %   Copyright 1984-2004 The MathWorks, Inc.
            %   $Revision: 5.12.4.2 $  $Date: 2004/06/25 18:52:28 $
            %   Modified by Richard @ Socher . org to display text labels

            p = obj.pp';
            [x,y,h]=treelayout(p);
            f = find(p~=0);
            pp = p(f);
            X = [x(f); x(pp); NaN(size(f))];
            Y = [y(f); y(pp); NaN(size(f))];
            X = X(:);
            Y = Y(:);

            n = length(p);
            if n < 500,
                plot (x, y, 'wo', X, Y, 'b-');
            else
                plot (X, Y, 'r-');
            end;
            xlabel(['height = ' int2str(h)]);
            axis([0 1 0 1]);

            if ~isempty(obj.nodeNames)
                for l=1:length(obj.nodeNames)
                        if isnumeric(obj.nodeNames(l))
                            text(x(l),y(l),num2str(obj.nodeNames(l)),'Interpreter','none',...
                                'HorizontalAlignment','center','FontSize',8,'BackgroundColor',[1 1 .6])
                        else
                            text(x(l),y(l),obj.nodeNames{l},'Interpreter','none',...
                                'HorizontalAlignment','center','FontSize',8,'BackgroundColor',[1 1 .6])
                        end
%                     if ~isempty(obj.nodeLabels)
%                         if iscell(obj.nodeNames)
%                             text(x(l),y(l),[labels{l} '(' obj.nodeLabels{l} ')'],'Interpreter','none',...
%                                 'HorizontalAlignment','center','FontSize',8,'BackgroundColor',[1 1 .6])
%                         else
%                             % for numbers
%                             if isnumeric(obj.nodeLabels(l))
% %                                 if isinteger(obj.nodeLabels(l))
%                                     allL = obj.nodeLabels(:,l);
%                                     allL = find(allL);
%                                     if isempty(allL)
%                                         text(x(l),y(l),[num2str(obj.nodeNames(l))],'Interpreter','none',...
%                                         'HorizontalAlignment','center','FontSize',8,'BackgroundColor',[1 1 .6])
%                                     else
%                                         text(x(l),y(l),[num2str(obj.nodeNames(l)) ' (' mat2str(allL) ')'],'Interpreter','none',...
%                                             'HorizontalAlignment','center','FontSize',8,'BackgroundColor',[1 1 .6])
%                                     end
%
% %                                 else
% %                                     text(x(l),y(l),[obj.nodeLabels(l) ' ' num2str(obj.nodeLabels(l),'%.1f') ],'Interpreter','none',...
% %                                         'HorizontalAlignment','center','FontSize',8,'BackgroundColor',[1 1 .6])
% %                                 end
%                                 % change to font size 6 for nicer tree prints
%                             else
%                                 text(x(l),y(l),[obj.nodeNames{l}],'Interpreter','none',...
%                                     'HorizontalAlignment','center','FontSize',8,'BackgroundColor',[1 1 .6])
%                             end
%                         end
%                     end
                end
            end

        end

    end
end

author: jacoxu


text:
%  Return time frame of a NIFTI dataset. Support both *.nii and
%  *.hdr/*.img file extension. If file extension is not provided,
%  *.hdr/*.img will be used as default.
%
%  It is a lightweighted "load_nii_hdr", and is equivalent to
%  hdr.dime.dim(5)
%
%  Usage: [ total_scan ] = get_nii_frame(filename)
%
%  filename - NIFTI file name.
%
%  Returned values:
%
%  total_scan - total number of image scans for the time frame
%
%  NIFTI data format can be found on: http://nifti.nimh.nih.gov
%
%  - Jimmy Shen (jimmy@rotman-baycrest.on.ca)
%
function [ total_scan ] = get_nii_frame(fileprefix)

   if ~exist('fileprefix','var'),
      error('Usage: [ hdr, fileprefix, machine ] = get_nii_frame(filename)');
   end

   if ~exist('machine','var'), machine = 'ieee-le'; end

   new_ext = 0;

   if findstr('.nii',fileprefix)
      new_ext = 1;
      fileprefix = strrep(fileprefix,'.nii','');
   end

   if findstr('.hdr',fileprefix)
      fileprefix = strrep(fileprefix,'.hdr','');
   end

   if findstr('.img',fileprefix)
      fileprefix = strrep(fileprefix,'.img','');
   end

   if new_ext
      fn = sprintf('%s.nii',fileprefix);

      if ~exist(fn)
         msg = sprintf('Cannot find file "%s.nii".', fileprefix);
         error(msg);
      end
   else
      fn = sprintf('%s.hdr',fileprefix);

      if ~exist(fn)
         msg = sprintf('Cannot find file "%s.hdr".', fileprefix);
         error(msg);
      end
   end

   fid = fopen(fn,'r',machine);

   if fid < 0,
      msg = sprintf('Cannot open file %s.',fn);
      error(msg);
   else
      hdr = read_header(fid);
      fclose(fid);
   end

   if hdr.sizeof_hdr ~= 348
      % first try reading the opposite endian to 'machine'
      switch machine,
      case 'ieee-le', machine = 'ieee-be';
      case 'ieee-be', machine = 'ieee-le';
      end

      fid = fopen(fn,'r',machine);

      if fid < 0,
         msg = sprintf('Cannot open file %s.',fn);
         error(msg);
      else
         hdr = read_header(fid);
         fclose(fid);
      end
   end

   if hdr.sizeof_hdr ~= 348
      % Now throw an error
      msg = sprintf('File "%s" is corrupted.',fn);
      error(msg);
   end

   total_scan = hdr.dim(5);

   return;					% get_nii_frame

%---------------------------------------------------------------------
function [ dsr ] = read_header(fid)

    fseek(fid,0,'bof');
    dsr.sizeof_hdr    = fread(fid,1,'int32')';  % should be 348!

    fseek(fid,40,'bof');
    dsr.dim           = fread(fid,8,'int16')';

    return;					% read_header


text:
function display(obj)
% Display a NIFTI-1 object
%__________________________________________________________________________
% Copyright (C) 2005-2017 Wellcome Trust Centre for Neuroimaging

%
% $Id: display.m 7147 2017-08-03 14:07:01Z spm $

disp(' ');
disp([inputname(1),' = '])
disp(' ');
disp(obj)
disp(' ')

author: spm


text:
function proj2d = projection(data3d,param, angle_rad)

proj2d = zeros(param.nu,param.nv,'single');

[uu,vv] = meshgrid(param.us,param.vs);

[xx,yy] = meshgrid(param.xs,param.ys);

rx = xx.*cos(angle_rad) - yy.*sin(angle_rad);
ry = xx.*sin(angle_rad) + yy.*cos(angle_rad);

for iz = 1:param.nz

    data3d(:,:,iz) = interp2(xx,yy ,data3d(:,:,iz),rx,ry,'linear');

end

data3d(isnan(data3d))=0;

data3d = permute(data3d,[1 3 2]);

[xx,zz] = meshgrid(param.xs,param.zs);

for iy = 1:param.ny

    Ratio = (param.ys(iy)+param.DSO)/(param.DSD);

    pu = uu*Ratio;
    pv = vv*Ratio;

    pu = (pu - xx(1,1))/(param.dx)+1;
    pv = (pv - zz(1,1))/(param.dz)+1;

    tmp = interp2(data3d(:,:,iy),pv,pu,'linear');

    tmp(isnan(tmp))=0;

    proj2d = proj2d + tmp';
end

dist = sqrt((param.DSD)^2 + uu.^2 + vv.^2)./(param.DSD)*param.dy;

proj2d = proj2d .* dist';


text:
function varargout = cluster_tool_subdivide(varargin)
% CLUSTER_TOOL_SUBDIVIDE M-file for cluster_tool_subdivide.fig
%      CLUSTER_TOOL_SUBDIVIDE, by itself, creates a new CLUSTER_TOOL_SUBDIVIDE or raises the existing
%      singleton*.
%
%      H = CLUSTER_TOOL_SUBDIVIDE returns the handle to a new CLUSTER_TOOL_SUBDIVIDE or the handle to
%      the existing singleton*.
%
%      CLUSTER_TOOL_SUBDIVIDE('CALLBACK',hObject,eventData,handles,...) calls the local
%      function named CALLBACK in CLUSTER_TOOL_SUBDIVIDE.M with the given input arguments.
%
%      CLUSTER_TOOL_SUBDIVIDE('Property','Value',...) creates a new CLUSTER_TOOL_SUBDIVIDE or raises the
%      existing singleton*.  Starting from the left, property value pairs are
%      applied to the GUI before cluster_tool_subdivide_OpeningFunction gets called.  An
%      unrecognized property name or invalid value makes property application
%      stop.  All inputs are passed to cluster_tool_subdivide_OpeningFcn via varargin.
%
%      *See GUI Options on GUIDE's Tools menu.  Choose "GUI allows only one
%      instance to run (singleton)".
%
% See also: GUIDE, GUIDATA, GUIHANDLES

% Edit the above text to modify the response to help cluster_tool_subdivide

% Last Modified by GUIDE v2.5 04-Sep-2006 11:37:00

% Begin initialization code - DO NOT EDIT
gui_Singleton = 1;
gui_State = struct('gui_Name',       mfilename, ...
                   'gui_Singleton',  gui_Singleton, ...
                   'gui_OpeningFcn', @cluster_tool_subdivide_OpeningFcn, ...
                   'gui_OutputFcn',  @cluster_tool_subdivide_OutputFcn, ...
                   'gui_LayoutFcn',  [] , ...
                   'gui_Callback',   []);
if nargin && ischar(varargin{1})
    gui_State.gui_Callback = str2func(varargin{1});
end

if nargout
    [varargout{1:nargout}] = gui_mainfcn(gui_State, varargin{:});
else
    gui_mainfcn(gui_State, varargin{:});
end
% End initialization code - DO NOT EDIT

% --- Executes just before cluster_tool_subdivide is made visible.
function cluster_tool_subdivide_OpeningFcn(hObject, eventdata, handles, varargin)
movegui(hObject,'center');
set(handles.title,'string',[cluster_tool_getver]);
handles.pwd=pwd;handles.fdir=pwd;
for k=1:size(varargin(:),1)
    if strcmp(varargin{k},'cl'),handles.cl=varargin{k+1};
    elseif strcmp(varargin{k},'fdir'),handles.fdir=varargin{k+1};
    end
end
handles.criterion=0;handles.k=1;handles.min_adj=0;
% This function has no output args, see OutputFcn.
% hObject    handle to figure
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)
% varargin   command line arguments to cluster_tool_subdivide (see VARARGIN)

% Choose default command line output for cluster_tool_subdivide
handles.output = hObject;

% Update handles structure
guidata(hObject, handles);

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

% --- Outputs from this function are returned to the command line.
function varargout = cluster_tool_subdivide_OutputFcn(hObject, eventdata, handles)
% varargout  cell array for returning output args (see VARARGOUT);
% hObject    handle to figure
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

% Get default command line output from handles structure
varargout{1} = handles.output;

% --- Executes on button press in quit.
function quit_Callback(hObject, eventdata, handles)
delete(handles.figure1);

% --- Executes on button press in ret.
function ret_Callback(hObject, eventdata, handles)
if isfield(handles,'cl'),cluster_tool_tools('cl',handles.cl,'fdir',handles.fdir);
else cluster_tool_tools('fdir',handles.fdir);end
delete(handles.figure1);

% --- Executes on button press in save.
function save_Callback(hObject, eventdata, handles)
cluster_tool_savefunc

% --- Executes on selection change in adjacency.
function adjacency_Callback(hObject, eventdata, handles)

% --- Executes during object creation, after setting all properties.
function adjacency_CreateFcn(hObject, eventdata, handles)
set(hObject,'String','Surface|Edge|Corner');
set(hObject,'Value',2);
disp('**If a warning prints, you can safely ignore it.**')
% Hint: popupmenu controls usually have a white background on Windows.
%       See ISPC and COMPUTER.
if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
    set(hObject,'BackgroundColor','white');
end

% --- Executes on button press in done.
function done_Callback(hObject, eventdata, handles)
if get(handles.adjacency,'Value')==1
    adjacency='surface';
elseif get(handles.adjacency,'Value')==2
    adjacency='edge';
elseif get(handles.adjacency,'Value')==3
    adjacency='corner';
end
cl=cl_subdivide(handles.cl,handles.criterion,handles.min_adj,adjacency,handles.k);
if ~isempty(cl)
    handles.cl=cl;
else
    warndlg('Resulting ''cl'' structure was empty, the initial structure is being kept in memory.');
end
guidata(hObject,handles);

function criterion_Callback(hObject, eventdata, handles)
try handles.criterion=str2double(get(hObject,'string'));catch set(hObject,'string','0');handles.criterion=0;end
guidata(hObject,handles);

% --- Executes during object creation, after setting all properties.
function criterion_CreateFcn(hObject, eventdata, handles)

% Hint: edit controls usually have a white background on Windows.
%       See ISPC and COMPUTER.
if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
    set(hObject,'BackgroundColor','white');
end

function k_Callback(hObject, eventdata, handles)
try handles.k=str2double(get(hObject,'string'));catch set(hObject,'string','1');handles.k=1;end
guidata(hObject,handles);

% --- Executes during object creation, after setting all properties.
function k_CreateFcn(hObject, eventdata, handles)

% Hint: edit controls usually have a white background on Windows.
%       See ISPC and COMPUTER.
if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
    set(hObject,'BackgroundColor','white');
end

function minadj_Callback(hObject, eventdata, handles)
try handles.min_adj=str2double(get(hObject,'string'));catch set(hObject,'string','0');handles.min_adj=0;end
guidata(hObject,handles);

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

% Hint: edit controls usually have a white background on Windows.
%       See ISPC and COMPUTER.
if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
    set(hObject,'BackgroundColor','white');
end

author: canlab

repo: CanlabCore

save_path: github-repos/MATLAB/canlab-CanlabCore


text:
function [ypred]=knn(xapp,yapp,valY,X,k)

%
% knn implementation
%
% USE : [ypred]=knn(xapp,yapp,valY,X,k)
%
% Vincent Guigue 08/01/03

% check nargin

if nargin<4
  error('too few argumemnts');
elseif nargin<5
  k=3;
else
  if mod(k,2)==0
    error('k must be odd');
  end
end

if size(xapp,2)~=size(X,2)
  error('dimension incompatibility');
end

ndim = size(xapp,2);
nptxapp = size(xapp,1);
nptX = size(X,1);

% distance de X a xapp :
mat1 =  repmat(xapp, nptX,1);
%mat21 = reshape(X',1,nptX*ndim)
mat22 = repmat(X,1,nptxapp)';
mat2 = reshape(mat22 ,ndim, nptxapp*nptX)';
distance = mat1 - mat2 ;

distance = sum(distance.^2,2);
distance = reshape(distance,nptxapp,nptX);
[val kppv] = sort(distance,1);

% bilan sur les k premieres lignes
kppv = reshape(kppv(1:k,:),k*nptX,1);
Ykppv = yapp(kppv,1);
Ykppv = reshape(Ykppv,k,nptX);

% trouver le plus de reponses identique par colonne

vote = [];
for i=1:nptX
  for j=1:length(valY)
    vote(j,i)=size(find(Ykppv(:,i)==valY(j)),1);
  end
end

[val ind]=max(vote,[],1);
ypred = valY(ind);

author: bayesnet

repo: bnt

save_path: github-repos/MATLAB/bayesnet-bnt


text:
function c = r8ut_mxm ( n, a, b )

%*****************************************************************************80
%
%% R8UT_MXM multiplies two R8UT matrices.
%
%  Discussion:
%
%    The R8UT storage format is used for an M by N upper triangular
%    matrix.  The format stores all M*N entries, even those which are zero.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    06 February 2004
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, integer N, the order of the matrices.
%    N must be positive.
%
%    Input, real A(N,N), B(N,N), the R8UT factor matrices.
%
%    Output, real C(N,N), the R8UT product matrix.
%
  c = zeros ( n, n );

  for i = 1 : n
   for j = i : n
      for k = i : j
        c(i,j) = c(i,j) + a(i,k) * b(k,j);
      end
    end
  end

  return
end

author: johannesgerer

repo: jburkardt-m


text:
%% Plotting Spherical Functions
%
%%
% In this chapter various ways of plotting spherical functions are
% explained. We start by defining some example functions.

% the smiley
sF1 = S2Fun.smiley;

% some osilatory function
f = @(v) 0.1*(v.theta+sin(8*v.x).*sin(8*v.y));
sF2 = S2FunHarmonic.quadrature(f, 'bandwidth', 150);

%% Smooth Plot
% The default <S2Fun.plot.html |plot|> command generates a colored plot
% without contours

plot(sF1)

%%
% * |plot(sF1)| is the same as |pcolor(sF1)|

%% Contour Plot
% nonfilled contour plot plots only the contour lines
contour(sF1, 'LineWidth', 2);

%% Filled Contour Plot
% filled contour plot plots the contour lines
contourf(sF1, 'LineWidth', 2);

%% 3D Plot
% 3D plot of a sphere colored accordingly to the function values.
plot3d(sF1);

%% Surface Plot
% 3D plot where the radius of the sphere is transformed according to the function values
surf(sF2);

%% Section Plot
% Plot the intersection of the surf plot with a plane defined by a normal vector |v|

plotSection(sF2, zvector,'color','interp','linewidth',10)
colormap spring
mtexTitle('Flowerpower!')

%% Spectral Plot
% plotting the Fourier coefficients

close all
plotSpektra(sF1,'FontSize',15);

%%
% The more specific plot options are covered in the respective classes.

author: mtex-toolbox

repo: mtex

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


text:
function pass = test_chebcoeffs2( pref )
% Test chebpoly2

if ( nargin == 0)
    pref = chebfunpref;
end

tol = 1000*pref.cheb2Prefs.chebfun2eps;
j = 1;

% Rank-2 function
n = 10;
m = 8;
Tn = chebpoly(n);
Tm = chebpoly(m);
f = Tn * Tn' + Tm * Tn';
X = chebcoeffs2( f );
Exact = zeros(n+1); Exact(n+1,n+1) = 1; Exact(m+1,n+1) = 1;
pass(j) = norm( X - Exact ) < tol; j = j + 1;

%f = Tm * Tn';
% check inverses
Z = chebfun2.vals2coeffs( chebpolyval2( f ) );
pass(j) = norm( Z - Exact ) < tol; j = j + 1;

f = chebfun2( @(x,y) x + y );
pass(j) = isequal( size(coeffs2(f,2,1)), [2 1] ); j = j + 1;

end

author: chebfun

repo: chebfun

save_path: github-repos/MATLAB/chebfun-chebfun


text:
function grid_test ( code )

%*****************************************************************************80
%
%% GRID_TEST tests the grid routines.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    18 February 2006
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, string CODE, the code for the element.
%    Legal values include 'Q4', 'Q8', 'Q9', 'Q12', 'Q16', 'QL',
%    'T3', 'T4', 'T6' and 'T10'.
%

%
%  NODE is defined as a vector rather than a two dimensional array,
%  so that we can handle the various cases using a single array.
%
  fprintf ( 1, '\n' );
  fprintf ( 1, '  GRID_TEST: Test the grid routine for element "%s"\n', code );

  nelemx = 3;
  nelemy = 2;

  fprintf ( 1, '  Number of elements in X direction = %d\n', nelemx );
  fprintf ( 1, '  Number of elements in Y direction = %d\n', nelemy );

  element_order = order_code ( code );

  node_num = grid_node_num ( code, nelemx, nelemy );

  fprintf ( 1, '  Element order =       %d\n', element_order );
  fprintf ( 1, '  Nodes in grid =       %d\n', node_num );

  element_num = grid_element_num ( code, nelemx, nelemy );

  element_node = grid_element ( code, element_order, nelemx, nelemy );

  grid_print ( element_order, element_num, element_node );

  width = grid_width ( element_order, element_num, element_node );

  fprintf ( 1, '\n' );
  fprintf ( 1, '  The grid width is %d\n', width );

  return
end


text:
close all;
clc;

num_samples_per_digit = 10;
rng('default');
digit_set = 0:9;
K = length(digit_set);
% Number of subspaces
trial.K = K;
cluster_sizes = num_samples_per_digit*ones(1, K);
trial.cluster_sizes = cluster_sizes;
% maximum dimension for each subspace
D = 10;
trial.D = D;
S = sum(cluster_sizes);
trial.S = S;
% identify sample indices for each digit
sample_list = [];
% experiment number
r = 1;
for k=1:K
    digit = digit_set(k);
    digit_indices = md.digit_indices(digit);
    num_digit_samples = length(digit_indices);
    % initialize the random number generator for repeatability
    rng( (r-1) * K + k);
    choices = randperm(num_digit_samples, cluster_sizes(k));
    selected_indices = digit_indices(choices);
    % fprintf('%d ', selected_indices);
    % fprintf('\n');
    sample_list = [sample_list selected_indices];
end
[Y, true_labels] = md.selected_samples(sample_list);
% Perform PCA to reduce dimensionality
Y = spx.la.pca.low_rank_approx(Y, 100);
Yn = spx.norm.normalize_l2(Y);
if true
    fprintf('\n\n\n Point Statistics:\n\n');
    [M, S]  = size(Yn);
    angle_result = spx.cluster.subspace.nearest_same_subspace_neighbors_by_inner_product(Yn, cluster_sizes);
    spx.cluster.subspace.print_nearest_neighbor_result(angle_result);
end

if false
    mf = spx.graphics.Figures;
    % Ambient space dimension and number of data points
    [trial.M, trial.S] = size(Y);
    % Number of subspaces
    trial.K = num_samples_per_digit;
    % maximum dimension for each subspace
    trial.D = 5;
    trial.cluster_sizes = cluster_sizes;
    % Solve the sparse subspace clustering problem
    tstart = tic;
    try
        clustering_result = solver(Y, trial.D, trial.K);
    catch ME
        % we will move on to next one
        fprintf('Problem in processing this example. %s: %s\n', ME.identifier, ME.message);
        % move on to next example
        rethrow(ME);
        error('cannot continue.');
    end
    trial.elapsed_time = toc (tstart);
    % graph connectivity
    trial.connectivity = clustering_result.connectivity;
    % estimated number of clusters
    trial.estimated_num_subspaces = clustering_result.num_clusters;
    % Time to compare the clustering
    cluster_labels = clustering_result.labels;
    comparison_result = spx.cluster.clustering_error(cluster_labels, true_labels, trial.K);
    trial.clustering_error_perc = comparison_result.error_perc;
    trial.clustering_acc_perc = 100 - comparison_result.error_perc;
    % Compute the statistics related to subspace preservation
    spr_stats = spx.cluster.subspace.subspace_preservation_stats(clustering_result.Z, trial.cluster_sizes);
    trial.spr_error = spr_stats.spr_error;
    trial.spr_flag = spr_stats.spr_flag;
    trial.spr_perc = spr_stats.spr_perc;
    fprintf('\nclustering error: %0.2f %% , clustering accuracy: %0.2f %%, mean spr error: %0.2f preserving : %0.2f %%, connectivity: %0.2f, elapsed time: %0.2f sec', trial.clustering_error_perc, trial.clustering_acc_perc, spr_stats.spr_error, spr_stats.spr_perc, trial.connectivity, trial.elapsed_time);
    fprintf('\n\n');
    Z = full(abs(clustering_result.Z));
    s = spx.stats.format_descriptive_statistics(Z(:));
    fprintf(s);
    fprintf('\n');
    fprintf('Missed points: \n');
    fprintf('%d ', find(comparison_result.misses));
    fprintf('\n');
    mf.new_figure('Representations');
    imshow(Z);
end

author: indigits

repo: sparse-plex


text:
function s = mfileparse(mfile, mdirs, names, options)
%MFILEPARSE Parsing of an M-file to obtain synopsis, help and references
%  S = MFILEPARSE(MFILE, MDIRS, NAMES, OPTIONS) parses the M-file MFILE looking
%  for synopsis (function), H1 line, subroutines and todo tags (if any).
%  It also fills in a boolean array indicating whether MFILE calls M-files
%  defined by MDIRS (M-files directories) AND NAMES (M-file names).
%  The input OPTIONS comes from M2HTML: fields used are 'verbose', 'global'
%  and 'todo'.
%  Output S is a structure whose fields are:
%     o synopsis: char array (empty if MFILE is a script)
%     o h1line: short one-line description into the first help line
%     o subroutine: cell array of char containing subroutines synopsis
%     o hrefs: boolean array with hrefs(i) = 1 if MFILE calls mdirs{i}/names{i}
%     o todo: structure containing information about potential todo tags
%
%  See also M2HTML

%  Copyright (C) 2003 Guillaume Flandin <Guillaume@artefact.tk>
%  $Revision: 1.0 $Date: 2003/29/04 17:33:43 $

narginchk(3,4);
if nargin == 3,
	options = struct('verbose',1, 'globalHypertextLinks',0, 'todo',0);
end

%- Delimiters used in strtok: some of them may be useless (% " .), removed '.'
strtok_delim = sprintf(' \t\n\r(){}[]<>+-*~!|\\@&/,:;="''%%');

%- Open for reading the M-file
fid = openfile(mfile,'r');
it = 0; % line number

%- Initialize Output
s = struct('synopsis',   '', ...
		   'h1line',     '', ...
		   'subroutine', {{}}, ...
		   'hrefs',      sparse(1,length(names)), ...
		   'todo',       struct('line',[],'comment',{{}}), ...
		   'ismex',      zeros(size(mexexts)));

%- Initialize flag for synopsis cont ('...')
flagsynopcont = 0;
%- Look for synopsis and H1 line
%  Help is the first set of contiguous comment lines in an m-file
%  The H1 line is a short one-line description into the first help line
while 1
	tline = fgetl(fid);
	if ~ischar(tline), break, end
	it = it + 1;
	tline = deblank(fliplr(deblank(fliplr(tline))));
	%- Synopsis line
	if ~isempty(strmatch('function',tline))
		s.synopsis = tline;
		if ~isempty(strmatch('...',fliplr(tline)))
			flagsynopcont = 1;
			s.synopsis = deblank(s.synopsis(1:end-3));
		end
	%- H1 Line
	elseif ~isempty(strmatch('%',tline))
		% allow for the help lines to be before the synopsis
		if isempty(s.h1line)
			s.h1line = fliplr(deblank(tline(end:-1:2)));
		end
		if ~isempty(s.synopsis), break, end
	%- Go through empty lines
	elseif isempty(tline)

	%- Code found. Stop.
	else
		if flagsynopcont
			if isempty(strmatch('...',fliplr(tline)))
				s.synopsis = [s.synopsis tline];
				flagsynopcont = 0;
			else
				s.synopsis = [s.synopsis deblank(tline(1:end-3))];
			end
		else
			break;
		end
	end
end

%- Global Hypertext Links option
%  If false, hypertext links are done only among functions in the same
%  directory.
if options.globalHypertextLinks
	hrefnames = names;
else
	indhref = find(strcmp(fileparts(mfile),mdirs));
	hrefnames = {names{indhref}};
end

%- Compute cross-references and extract subroutines
%  hrefs(i) is 1 if mfile calls mfiles{i} and 0 otherwise
while ischar(tline)
	% Remove blanks at both ends
	tline = deblank(fliplr(deblank(fliplr(tline))));

	% Split code into meaningful chunks
	splitc = splitcode(tline);
	for j=1:length(splitc)
		if isempty(splitc{j}) | ...
			splitc{j}(1) == '''' | ...
			~isempty(strmatch('...',splitc{j}))
			% Forget about empty lines, char strings or conts
		elseif splitc{j}(1) == '%'
			% Cross-references are not taken into account in comments
			% Just look for potential TODO line
			if options.todo
				if ~isempty(strmatch('% TODO %',splitc{j}))
					s.todo.line   = [s.todo.line it];
					s.todo.comment{end+1} = splitc{j}(9:end);
				end
			end
		else
			% detect if this line is a declaration of a subroutine
			if ~isempty(strmatch('function',splitc{j}))
				s.subroutine{end+1} = splitc{j};
			else
				% get list of variables and functions
				symbol = {};
				while 1
					[t,splitc{j}] = strtok(splitc{j},strtok_delim);
					if isempty(t), break, end;
					symbol{end+1} = t;
				end
				if options.globalHypertextLinks
					s.hrefs = s.hrefs + ismember(hrefnames,symbol);
				else
					s.hrefs(indhref) = s.hrefs(1,indhref) + ...
									   ismember(hrefnames,symbol);
				end
			end
		end
	end
	tline = fgetl(fid);
	it = it + 1;
end

fclose(fid);

%- Look for Mex files
[pathstr,name] = fileparts(mfile);
samename = dir(fullfile(pathstr,[name	'.*']));
samename = {samename.name};
ext = {};
for i=1:length(samename)
	[dummy, dummy, ext{i}] = fileparts(samename{i});
end
s.ismex = ismember(mexexts,ext);

author: pdollar

repo: toolbox

save_path: github-repos/MATLAB/pdollar-toolbox


text:
filename='Bridge_quad_coarse';
 ptype = 'MACRO';
method = 'SIMPALL';
materialType = 'ISOTROPIC';
initial_case = 'full';
cost = {'compliance';'perimeter'};
weights = [1, 0.1];
constraint = {'volumeConstraint'};
optimizer = 'PROJECTED GRADIENT'; incrementFactor = 1;
filterType = 'P1';

nsteps = 5;
Vfrac_final = 0.2;
Perimeter_target = 1;
optimality_final = 1e-3;
constr_final =1e-3;

Vfrac_initial = 1;
optimality_initial = 1e-2;
constr_initial = 1e-3;
TOL.rho_plus = 1;
TOL.rho_minus = 0;
TOL.E_plus = 1;
TOL.E_minus = 1e-3;
TOL.nu_plus = 1/3;
TOL.nu_minus = 1/3;

author: SwanLab

repo: Swan

save_path: github-repos/MATLAB/SwanLab-Swan


text:
classdef MatNatProject < MatNatBase
    % MatNatProject A Matlab class representing an XNAT project
    %
    %     Licence
    %     -------
    %     Part of MatNat. https://github.com/tomdoel/matnat
    %     Author: Tom Doel, 2015.  www.tomdoel.com
    %     Distributed under the GNU GPL v3 licence. Please see website for details.
    %

    properties (SetAccess = protected)
        Name
        Id
        SecondaryId
        Description
    end

    properties (Access = private)
        RestClient
        SubjectMap
    end

    methods
        function obj = MatNatProject(restClient)
            obj.RestClient = restClient;
        end

        function subjectMap = getSubjectMap(obj)
            % Returns a map of subject identifiers to MatNatSubject objects

            obj.populateSubjectMapIfNecessary;
            subjectMap = obj.SubjectMap;
        end

        function subject = getSubject(obj, patient_id)
            % Returns the MatNatSubject object corresponding to the given
            % subject identifier

            if isempty(patient_id)
                subject = [];
                return;
            end

            obj.populateSubjectMapIfNecessary;
            if ~obj.SubjectMap.isKey(patient_id)
                subject = [];
            else
                subject = obj.SubjectMap(patient_id);
            end
        end

        function resource = GetResourceForSeriesUid(obj, patient_id, series_uid)
            % Gets the first resource for the series with this identifier

            subject = obj.getSubject(patient_id);
            if ~isempty(subject)
                scan = subject.FindScan(series_uid);
                if ~isempty(scan)
                    resources = scan.getResources;
                    if numel(resources) > 0
                        resource = resources(1);
                        return;
                    end
                end
            end
            resource = [];
        end
    end

    methods (Access = private)
        function populateSubjectMapIfNecessary(obj)
            if isempty(obj.SubjectMap)
                obj.SubjectMap = obj.RestClient.getSubjectMap(obj.Id);
            end
        end
    end

    methods (Static)
        function obj = createFromServerObject(restClient, serverObject)
            % Creates a MatNatProject based on the project information
            % structure returned from the XNAT server

            obj = MatNatProject(restClient);
            obj.Name = MatNatBase.getOptionalProperty(serverObject, 'name');
            obj.Id = MatNatBase.getOptionalProperty(serverObject, 'id');
            obj.SecondaryId = MatNatBase.getOptionalProperty(serverObject, 'secondary_id');
            obj.Description = MatNatBase.getOptionalProperty(serverObject, 'description');
        end
    end
end

author: tomdoel

repo: pulmonarytoolkit


function [local_ids, labels, scores2, THRESH]=paste_segments(chosenpercateg,scores, region_meta_info,MAX_AVG_N_SEGM, NMS_MAX_SEGMS, THRESH)
n_segms_per_img = inf;

NMS_MAX_OVER = 0; % spatial non-maximum supression (1 means nms is not performed)

SIMP_BIAS_STEP = 0;%0.02; % background threshold increase for each additional segment above 1

% max number of segments per image on average (used to set the background threshold)
% we set this value to the average number of objects in the
% training set. Of course, that is just a coincidence. ;-)
%MAX_AVG_N_SEGM = 2.2;
for i=1:numel(scores)
    scores{i}=scores{i}(region_meta_info.gt{i}==0,:);
    for j=1:20
        scores{i}(~ismember([1:size(scores{i},1)], chosenpercateg{j}{i}),j)=-inf;
    end
end

chosen=chosenpercateg{1};
for k=2:20
    for i=1:numel(chosen)
        chosen{i}=[chosen{i}(:); chosenpercateg{k}{i}(:)];
    end
end
for i=1:numel(chosen)
    chosen{i}=unique(chosen{i});
end
while(n_segms_per_img > MAX_AVG_N_SEGM)

    [local_ids, labels, scores2] = nms_inference_simplicity_bias(chosen, scores, NMS_MAX_OVER, NMS_MAX_SEGMS, SIMP_BIAS_STEP, THRESH);

    n_segms_per_img = numel(cell2mat(labels')) / numel(labels)

    THRESH = THRESH
    THRESH = THRESH+0.01;
end
THRESH=THRESH-0.01;


text:
function a = frank_inverse ( n )

%*****************************************************************************80
%
%% FRANK_INVERSE returns the inverse of the FRANK matrix.
%
%  Formula:
%
%    if ( I = J-1 )
%      A(I,J) = -1
%    elseif ( I = J )
%      if ( I = 1 )
%        A(I,J) = 1
%      else
%        A(I,J) = N + 2 - I
%    elseif ( J < I )
%      A(I,J) = - (N+1-I) * A(I-1,J)
%    else
%      A(I,J) = 0
%
%  Example:
%
%    N = 5
%
%     1  -1  0  0  0
%    -4   5 -1  0  0
%    12 -15  4 -1  0
%   -24  30 -8  3 -1
%    24 -30  8 -3  2
%
%  Properties:
%
%    A is generally not symmetric: A' /= A.
%
%    A is lower Hessenberg.
%
%    det ( A ) = 1.
%
%    A is unimodular.
%
%    A is integral: int ( A ) = A.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    15 October 2007
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, integer N, the order of A.
%
%    Output, real A(N,N), the matrix.
%
  a = zeros ( n, n );

  for i = 1 : n
    for j = 1 : n

      if ( i == j-1 )
        a(i,j) = - 1.0;
      elseif ( i == j )
        if ( i == 1 )
          a(i,j) = 1.0;
        else
          a(i,j) = ( n + 2 - i );
        end
      elseif ( j < i )
        a(i,j) = - ( n + 1 - i ) * a(i-1,j);
      else
        a(i,j) = 0.0;
      end

    end
  end

  return
end


text:
function [flag] = isHrfStable(X,kas,kaf,phi)

try
    phi;
catch
    phi = -0.5785;
end

% get vasodilatory signal and blood inflow
X1 = X(1);
X2 = X(2);

% if X(1) > 0
%     flag = 1;
%     return
% end

% first compute separatrix
x2 = [-5:1e-3:5]';
n = length(x2);
x1 = phi.*exp(kaf/2-abs(kas-kaf/2)).*exp(x2);
dx1 = diff(x1);
dx2 = diff(x2);

% get closest point on separatrix
d2 = sum([(x1-X1).^2,(x2-X2).^2],2);
[md,ind] = min(d2);

% see whether the point is inside the domain
vec = [x2(ind)-X2,x1(ind)-X1];
tan = [-dx1(ind),dx2(ind)];
nv = sqrt(sum(vec.^2));
nt = sqrt(sum(tan.^2));
angle = acos(sum(vec.*tan)./(nv*nt));
flag = -sign(cos(angle));

hf = figure;
ha = axes(...
    'parent',hf,...
    'nextplot','add');
plot(ha,x2,x1,'.')
plot(ha,x2(ind),x1(ind),'ro')
quiver(x2(1:end-1),...
    x1(1:end-1),...
    -diff(x1),...
    diff(x2),...
    2)
quiver(X2,...
    X1,...
    vec(1),...
    vec(2),...
	'color',[1 0 0])
grid(ha,'on')
axis(ha,'tight')

plot(ha,X2,X1,'r+')
plot(ha,X2,X1,'ro')

author: MBB-team

repo: VBA-toolbox


text:
function n = ndims( x )
n = length( size( x ) );

% Copyright 2010 Michael C. Grant and Stephen P. Boyd.
% See the file COPYING.txt for full copyright information.
% The command 'cvx_where' will show where this file is located.

repo: REST


text:
%% Gammatone-like spectrograms
% Gammatone filters are a popular linear approximation to the
% filtering performed by the ear.  This routine provides a simple
% wrapper for generating time-frequency surfaces based on a
% gammatone analysis, which can be used as a replacement for a
% conventional spectrogram.  It also provides a fast approximation
% to this surface based on weighting the output of a conventional
% FFT.

%% Introduction
% It is very natural to visualize sound as a time-varying
% distribution of energy in frequency - not least because this is
% one way of describing the information our brains get from our
% ears via the auditory nerve.  The spectrogram is the traditional
% time-frequency visualization, but it actually has some important
% differences from how sound is analyzed by the ear, most
% significantly that the ear's frequency subbands get wider for
% higher frequencies, whereas the spectrogram has a constant
% bandwidth across all frequency channels.
%
% There have been many signal-processing approximations proposed
% for the frequency analysis performed by the ear; one of the most
% popular is the Gammatone filterbank originally proposed by
% Roy Patterson and colleagues in 1992.  Gammatone filters were
% conceived as a simple fit to experimental observations of
% the mammalian cochlea, and have a repeated pole structure leading
% to an impulse response that is the product of a Gamma envelope
% g(t) = t^n e^{-t} and a sinusoid (tone).
%
% One reason for the popularity of this approach is the
% availability of an implementation by Malcolm Slaney, as
% described in:
%
% Malcolm Slaney (1998) "Auditory Toolbox Version 2",
% Technical Report #1998-010, Interval Research Corporation, 1998.
% http://cobweb.ecn.purdue.edu/~malcolm/interval/1998-010/
%
% Malcolm's toolbox includes routines to design a Gammatone
% filterbank and to process a signal by every filter in a bank,
% but in order to convert this into a time-frequency visualization
% it is necessary to sum up the energy within regular time bins.
% While this is not complicated, the function here provides a
% convenient wrapper to achieve this final step, for applications
% that are content to work with time-frequency magnitude
% distributions instead of going down to the waveform levels.  In
% this mode of operation, the routine uses Malcolm's MakeERBFilters
% and ERBFilterBank routines.
%
% This is, however, quite a computationally expensive approach, so
% we also provide an alternative algorithm that gives very similar
% results.  In this mode, the Gammatone-based spectrogram is
% constructed by first calculating a conventional, fixed-bandwidth
% spectrogram, then combining the fine frequency resolution of the
% FFT-based spectra into the coarser, smoother Gammatone responses
% via a weighting function.  This calculates the time-frequency
% distribution some 30-40x faster than the full approach.

%% Routines
% The code consists of a main routine, <gammatonegram.m gammatonegram>,
% which takes a waveform and other parameters and returns a
% spectrogram-like time-frequency matrix, and a helper function
% <fft2gammatonemx.m fft2gammatonemx>, which constructs the
% weighting matrix to convert FFT output spectra into gammatone
% approximations.

%% Example usage
% First, we calculate a Gammatone-based spectrogram-like image of
% a speech waveform using the fast approximation.  Then we do the
% same thing using the full filtering approach, for comparison.

% Load a waveform, calculate its gammatone spectrogram, then display:
[d,sr] = wavread('sa2.wav');
tic; D = gammatonegram(d,sr); toc
%Elapsed time is 0.140742 seconds.
subplot(211)
imagesc(20*log10(D)); axis xy
caxis([-90 -30])
colorbar
title('Gammatonegram - fast method')

% Now repeat with flag to use actual subband filters.
% Since it's the last argument, we have to include all the other
% arguments.  These are the default values for: summation window
% (0.025 sec), hop between successive windows (0.010 sec),
% number of gammatone channels (64), lowest frequency (50 Hz),
% and highest frequency (sr/2).  The last argument as zero
% means not to use the FFT approach.
tic; D2 = gammatonegram(d,sr,0.025,0.010,64,50,sr/2,0); toc
%Elapsed time is 3.165083 seconds.
subplot(212)
imagesc(20*log10(D2)); axis xy
caxis([-90 -30])
colorbar
title('Gammatonegram - accurate method')
% Actual gammatone filters appear somewhat narrower.  The fast
% version assumes coherence of addition of amplitude from
% different channels, whereas the actual subband energies will
% depend on how the energy in different frequencies combines.
% Also notice the visible time smearing in the low frequency
% channels that does not occur in the fast version.

%% Validation
% We can check the frequency responses of the filterbank
% simulated with the fast method against the actual filters
% from Malcolm's toolbox.  They match very closely, but of
% course this still doesn't mean the two approaches will give
% identical results - because the fast method ignores the phase
% of each frequency channel when summing up.

% Check the frequency responses to see that they match:
% Put an impulse through the Slaney ERB filters, then take the
% frequency response of each impulse response.
fcfs = flipud(MakeERBFilters(16000,64,50));
gtir = ERBFilterBank([1, zeros(1,1000)],fcfs);
H = zeros(64,512);
for i = 1:64; H(i,:) = abs(freqz(gtir(i,:),1,512)); end
% The weighting matrix for the FFT is the frequency response
% of each output filter
gtm = fft2gammatonemx(1024,16000,64,1,50,8000,512);
% Plot every 5th channel from both.  Offset by 3 dB just so we can
% see both
fs = [0:511]/512*8000;
figure
plot(fs,20*log10(H(5:5:64,:))','b',fs, -3 + 20*log10(gtm(5:5:64,:))','r')
axis([0 8000 -150 0])
grid
% Line up pretty well, apart from wiggles below -100 dB
% (from truncating the impulse response at 1000 samples?)

%% Download
% You can download all the code and data for these examples here:
% <gammatone.tgz gammatone.tgz>.

%% Referencing
% If you use this work in a publication, I would be grateful
% if you referenced this page as follows:
%
%  D. P. W. Ellis (2009).  "Gammatone-like spectrograms", web resource, http://www.ee.columbia.edu/~dpwe/resources/matlab/gammatonegram/ .

%% Acknowledgment
% This project was supported in part by the NSF under
% grant IIS-0535168. Any opinions, findings and conclusions
% or recommendations expressed in this material are those of the
% authors and do not necessarily reflect the views of the Sponsors.

% Last updated: $Date: 2009/02/22 01:46:42 $
% Dan Ellis <dpwe@ee.columbia.edu>

author: detly

repo: gammatone

save_path: github-repos/MATLAB/detly-gammatone


text:
function browse_topoplotVAR(cfg, data)

% BROWSE_TOPOPLOTVAR is a simple helper function for FT_DATABROWSER that
% computes the variance of band-pass filtered data and makes a topographic
% plot. It serves to make a quick-and-dirty power topography.
%
% See also BROWSE_MOVIEPLOTER, BROWSE_TOPOPLOTER, BROWSE_MULTIPLOTER, BROWSE_TOPOPLOTVAR, BROWSE_SIMPLEFFT

% Copyright (C) 2009, Robert Oostenveld

% compute the variance, i.e. the broad-band power
timelock        = [];
timelock.label  = data.label;
timelock.time   = mean(data.time{1});
timelock.avg    = sum(ft_preproc_baselinecorrect(data.trial{1}).^2, 2);
timelock.dimord = 'chan_time';

if isfield(data, 'grad')
  timelock.grad = data.grad;
end
if isfield(data, 'elec')
  timelock.elec = data.elec;
end

default             = [];
default.markers     = 'labels';
default.interactive = 'no';
cfg = mergestruct(cfg, default);

figure;
ft_topoplotER(cfg, timelock);

author: fieldtrip

repo: fieldtrip

save_path: github-repos/MATLAB/fieldtrip-fieldtrip


text:
function [varargout]=kolmogorov_max_flow(A,u,v,varargin)
% KOLMOGOROV_MAX_FLOW Kolmogorov's max flow algorithm
%
% Kolmogorov's algorithm implements a variation on the augmenting path idea
% to compute a maximum flow or minimum cut in a network.
%
% See the max_flow function for calling information and return parameters.
% This function just calls max_flow(...,struct('algname','kolmogorov'));
%
% Example:
%   load('graphs/max_flow_example.mat');
%   kolmogorov_max_flow(A,1,8)

% David Gleich
% Copyright, Stanford University, 2007-2008

%% History
%  2007-07-07: Initial version
%  2008-10-07: Changed options parsing
%%

algname = 'kolmogorov';
if ~isempty(varargin),
    options = merge_options(struct(),varargin{:});
    options.algname= algname;
else options = struct('algname',algname);
end

varargout = cell(1,max(nargout,1));

[varargout{:}] = max_flow(A,u,v,options);


text:
function[varargout]=polysmooth_sortfield(varargin)
%POLYSMOOTH_SORTFIELD  Returns sorted field values for a mapping problem.
%
%   ZS=POLYSMOOTH_SORTFIELD(INDEXS,Z), where INDEXS is a cell array of
%   indices as output by SPHERESORT or TWODSORT, applies INDEXS to the
%   input field Z, leading to the output ZS.
%
%   ZS contains the elements of Z, sorted in terms of increasing distance
%   from a set of grid points, as specified by SPHERESORT or TWODSORT.
%
%   [ZS1,ZS2,...,ZSK]=POLYSMOOTH_SORTFIELD(INDEXS,Z1,Z2,...,ZK) also works.
%
%   This is mainly an auxillary function used with SPHERESORT or TWODSORT.
%   However, it is useful to call directly for the "One grid, many fields"
%   case described in POLYSMOOTH.
%
%   Usage: zs=polysmooth_sortfield(indexs,z);
%          [zs1,zs2,zs3]=polysmooth_sortfield(indexs,z1,z2,z3);
%   __________________________________________________________________
%   This is part of JLAB --- type 'help jlab' for more information
%   (C) 2018--2020 J.M. Lilly --- type 'help jlab_license' for details

indexs=varargin{1};
varargin=varargin(2:end);

varargout=cell(length(varargin),1);
for k=1:length(varargin)
    temp=indexs;
    for j=1:length(indexs)
        nonnani=~isnan(indexs{j});
        temp{j}(nonnani)=varargin{k}(indexs{j}(nonnani));
        %adjustment for missing data in complex-valued fields
        if ~isreal(varargin{k})
            temp{j}(~nonnani)=temp{j}(~nonnani)+1i*nan;
        end
    end
    varargout{k}=temp;
end

author: jonathanlilly

repo: jLab

save_path: github-repos/MATLAB/jonathanlilly-jLab


text:
function calpak_test394 ( )

%*****************************************************************************80
%
%% CALPAK_TEST394 tests MONTH_TO_MONTH_NAME_EG_CIVIL.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    02 April 2013
%
%  Author:
%
%    John Burkardt
%
  fprintf ( 1, '\n' );
  fprintf ( 1, 'CALPAK_TEST394\n' );
  fprintf ( 1, '  For the Egyptian Civil calendar,\n' );
  fprintf ( 1, '  MONTH_TO_MONTH_NAME_EG_CIVIL names the months.\n' );
  fprintf ( 1, '\n' );

  y = 1;
  months = year_length_months_eg_civil ( y );

  for m = 1 : months
    month_name = month_to_month_name_eg_civil ( m );
    fprintf ( 1, '  %2d  %s\n', m, month_name );
  end

  return
end


text:
function visInferredShapes(class, jobID)
% Visualize fitted basis shapes cached in 'inferredShapes'

% Setup directories
globals;
inferredShapesDir = fullfile(cachedir, class, sprintf('inferredShapes%s',jobID));
meshDir = fullfile(cachedir, class, sprintf('meshes%s',jobID));
dmapDir = fullfile(cachedir, class, sprintf('depthMap%s',jobID));
fnames = getFileNamesFromDirectory(inferredShapesDir,'types',{'.mat'});
hFig = figure;
load(fullfile(datadir,'partNames',class));
for i = 1:length(fnames)
    fprintf('%s: VOC ID = %s\n',class, fnames{i}(1:end-4));
    load(fullfile(inferredShapesDir,fnames{i}));
    try
        mesh = load(fullfile(meshDir,fnames{i}));
    catch
        mesh = state2mesh(state);
    end

    subplot(231)
    imshow(state.im);
    title('Image');

    subplot(232)
    im = (color_seg(state.mask,state.im));
    %im = insertText(im,state.kps(1:2,:)',partNames,'FontSize',8);
    imshow(im);
    hold on;
    cmap = distinguishable_colors(size(state.kps,2),[0 1 1]);
    scatter(state.kps(1,:), state.kps(2,:),30,cmap,'Filled');
    hold off;
    title('Image Mask (and Keypoints)');

    subplot(233)
    imVerts = mesh.vertices;
    vertNormals = state.normals*state.cameraRot';

    %Flip z axis to be in the image frame
    imVerts(:,3) = -imVerts(:,3);
    imVerts(:,3) = imVerts(:,3)-min(0,min(imVerts(:,3)));
    zDepth = imVerts(:,3);
    vertNormals(:,3) = -vertNormals(:,3);

    imshow(state.im); hold on;

    %patch('vertices', imVerts, 'faces', mesh.faces, 'CData', zDepth,'VertexNormals',vertNormals,...
    %    'CDataMapping','scaled', 'FaceColor', 'interp', 'FaceAlpha', ...
    %     0.7, 'EdgeColor', 'None');

     patch('vertices', imVerts, 'faces', state.tri, 'FaceColor', 'red', 'VertexNormals',...
         vertNormals, 'FaceAlpha', 0.7, 'EdgeColor','None','FaceLighting','gouraud',...
         'BackFaceLighting','lit');

    camlight; lighting p;
    camproj('orthographic');
    hold off;
    title('Mesh overlayed on Image');

    subplot(234)
    showMeshTri(mesh);view(0,-90);
    title('Inferred Mesh');

    subplot(235)
    % Colored Point Cloud
    imsize = size(state.im); imsize = imsize(1:2);
    imVerts = mesh.vertices;

    % Mask project vertices only in mask
    maskInds = find(state.mask(:));
    imVerts(:,1) = max(1,min(imVerts(:,1),imsize(2)));
    imVerts(:,2) = max(1,min(imVerts(:,2),imsize(1)));
    linInds = sub2ind(imsize,round(imVerts(:,2)),round(imVerts(:,1)));
    [linInds, IA] = intersect(linInds,maskInds);
    imVerts = imVerts(IA,:);
    % Find colors for the points
    rCh = state.im(:,:,1); rCh = rCh(:);
    gCh = state.im(:,:,2); gCh = gCh(:);
    bCh = state.im(:,:,3); bCh = bCh(:);
    ptColors = [rCh(linInds) gCh(linInds) bCh(linInds)];
    % Plot
    scatter3(imVerts(:,1),imVerts(:,2),imVerts(:,3),15,single(ptColors)/255,'filled');
    axis off vis3d equal; view(0,-90);
    title('Colored Point Cloud');

    subplot(236)
    % Texture mapped depth map / Depth map
    try
        load(fullfile(dmapDir,fnames{i}));
    catch
        mesh = reducepatch(mesh,2000);
        dmap = meshToDepth(mesh.vertices,mesh.faces,size(state.mask));
        dmap(isinf(dmap)) = nan;
        dmap(~state.mask) = nan;
    end

    hIm = imagesc(imcrop(visualizeDEM(dmap),state.bbox));
    %set(hIm,'AlphaData',repmat(~isnan(imcrop(dmap,state.bbox)),[1 1 3]));
    %warp(dmap,state.im); view(0,90);
    axis equal off;
    title('Depth map');

    %subplotsqueeze(hFig, 1.2);
    pause; clf
end

repo: CategoryShapes


text:
function [VIn, MIn] = partition_distance(Cx, Cy)
%PARTITION_DISTANCE     Distance or similarity between community partitions
%
%   This function quantifies information-theoretic distance (normalized
%   variation of information) or similarity (normalized mutual information)
%   between community partitions.
%
%   VIn        = partition_distance(Cx);
%   VIn        = partition_distance(Cx, Cy);
%   [VIn, MIn] = partition_distance(Cx, Cy);
%
%   Inputs:
%       Cx,
%           Community partition vector or matrix of n rows and p columns,
%           n is the number of network nodes, and p is the number of input
%           community partitions (in the case of vector input p=1).
%
%       Cy (optional argument),
%           Community partition vector or matrix of n rows and q columns. n
%           is the number of nodes (must be equal to the number of nodes in
%           Cq) and q is the number of input community partitions (may be
%           different to the number of nodes in Cq). This argument may be
%           omitted, in which case, the partition distance is computed
%           between all pairwise partitions of Cx.
%
%   Outputs:
%       VIn,
%           Normalized variation of information ([p, q] matrix)
%
%       MIn,
%           Normalized mutual information ([p, q] matrix)
%
%   Notes:
%       Mathematical definitions.
%
%           VIn = [H(X) + H(Y) - 2MI(X, Y)]/log(n)
%           MIn = 2MI(X, Y) / [H(X) + H(Y)]
%
%           where H is the entropy and MI is the mutual information
%
%
%   Reference: Meila M (2007) J Multivar Anal 98, 873-895.
%
%
%   2011-2017, Mika Rubinov, UNSW, Janelia HHMI

%   Modification History:
%   Mar 2011: Original
%   Jan 2017: Added computation between input matrices.

s = (nargin==1);
if s
    Cy = Cx;
    d = 10.^ceil(log10(double(1 + max( Cx(:)) )));
else
    d = 10.^ceil(log10(double(1 + max([Cx(:);Cy(:)]) )));
end

if ~isequal([Cx(:);Cy(:)], int64([Cx(:);Cy(:)])) || min([Cx(:);Cy(:)])<=0
    error('Input partitions must contain only positive integers.')
end

[n, p] = size(Cx);
HX = zeros(p, 1);
for i = 1:p
    Px = nonzeros(accumarray(Cx(:, i), 1)) / n;                     % P(x)
    HX(i) = - sum(Px .* log(Px));                                   % H(x)
end

if s
    q = p;
    HY = HX;
else
    [n_, q] = size(Cy);
    assert(n == n_);
    HY = zeros(q, 1);
    for j = 1:q
        Py = nonzeros(accumarray(Cy(:, j), 1)) / n;                 % P(y)
        HY(j) = - sum(Py .* log(Py));                               % H(y)
    end
end

VIn = zeros(p, q);
MIn = zeros(p, q);
for i = 1:p
    j_idx = (s * (i - 1) + 1):q;
    for j = j_idx
        Pxy = nonzeros(accumarray(d*Cx(:, i) + Cy(:, j), 1)) / n; 	% P(x,y)
        Hxy = -sum(Pxy .* log(Pxy));                                % H(x,y)
        VIn(i, j) = (2 * Hxy - HX(i) - HY(j)) / log(n);             % VIn
        MIn(i, j) = 2 * (HX(i) + HY(j) - Hxy) / (HX(i) + HY(j));    % MIn
    end
    if s
        VIn(j_idx, i) = VIn(i, j_idx);
        MIn(j_idx, i) = MIn(i, j_idx);
    end
end


text:
function value = c8_sqrt ( c )

%*****************************************************************************80
%
%% C8_SQRT returns the principal square root of a C8.
%
%  Discussion:
%
%    A C8 is a complex double precision value.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    15 April 2007
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, complex double precision C, the number whose square root is desired.
%
%    Output, complex double precision VALUE, the square root of X.
%
  argument = c8_arg ( c );
  magnitude = c8_mag ( c );

  if ( magnitude == 0.0 )

    value = 0.0;

  else

    value = sqrt ( magnitude ) ...
      * ( cos ( argument / 2.0 ) + i * sin ( argument / 2.0 ) );

  end

  return
end


text:
function [ id, type, rep, field, symm, nrow, ncol, nnz, indx, jndx, value ] ...
  = mm_file_read ( input_unit )

%*****************************************************************************80
%
%% MM_FILE_READ reads data from a Matrix Market file.
%
%  Discussion:
%
%    The data may be either sparse coordinate format, or dense array format.
%
%    The unit input_unit must be open, and the file will be rewound on return.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    02 May 2004
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, integer INPUT_UNIT, the input unit identifier.
%
%    Output, character ( len = 14 ) ID, the Matrix Market identifier.
%    This value must be '%%MatrixMarket'.
%
%    Output, character ( len = 6 ) TYPE, the Matrix Market type.
%    This value must be 'matrix'.
%
%    Output, character ( len = 10 ) REP, the Matrix Market 'representation'
%    indicator.  Possible values include:
%    'coordinate'   (for sparse data)
%    'array'        (for dense data)
%    'elemental'    (to be added)
%
%    Output, character ( len = 7 ) FIELD, the Matrix Market 'field'.
%    Possible values include:
%    'real'
%    'double'
%    'complex'
%    'integer'
%    'pattern'
%
%    Output, character ( len = 19 ) SYMM, the Matrix Market symmetry.
%    Possible values include:
%    'symmetric'
%    'hermitian'
%    'skew-symmetric'
%    'general'
%
%    Output, integer NROW, the number of rows in the matrix.
%
%    Output, integer NCOL, the number of columns in the matrix.
%
%    Output, integer NNZ, the number of nonzero entries required to store
%    the matrix, if REP = 'coordinate'.
%
%    Output, integer INDX(NNZ), the row indices for coordinate format.
%    Not used if REP is 'array'.
%
%    Output, integer JNDX(NNZ), the column indices for coordinate format.
%    Not used if REP is 'array'.
%
%    Output, complex VALUE(NNZ), the matrix values.
%

%
%  Read and check the header line.
%
  [ id, type, rep, field, symm ] = mm_header_read ( input_unit );

  ok = mm_header_check ( id, type, rep, field, symm );

  if ( ~ok )
    fprintf ( 1, '\n' );
    fprintf ( 1, 'MM_FILE_READ - Fatal error!\n' );
    fprintf ( 1, '  The header information was not acceptable.\n' );
    error ( 'MM_FILE_READ - Fatal error!' );
  end
%
%  Read through the comment lines:
%
  while ( 1 )

    comment = mm_comment_read ( input_unit );

    if ( comment(1) ~= '%' )
      break;
    end

  end
%
%  The line we just read is not a comment.  It is presumably the SIZE line.
%
  [ nrow, ncol, nnz ] = mm_size_read_string ( comment, rep, symm );
%
%  Read the data values.
%
  [ value, indx, jndx ] = mm_values_read ( input_unit, rep, field, nnz );

  return
end


text:
function issym=issymmetric(X)
%ISSYMMETRIC Check if variable is symmetric

n = X.dim(1);
m = X.dim(2);
if (n==m)
    % What are the linar indicies to the transposed matrices
    if isa(X.basis,'lazybasis')
        issym = 1;
    else
        ind = reshape(reshape(1:n^2,n,n)',n^2,1);
        issym = norm(X.basis-X.basis(ind,:),1)<1e-10;
    end
else
    issym = 0;
end

author: zarathustr

repo: LibQPEP

save_path: github-repos/MATLAB/zarathustr-LibQPEP


text:
function M = threed_bilinear( kernel, phi, test, w_g )
%-----------------------------------------------------------------------
%  threed_bilinear.m - routine to compute \int{ kernel*phi*test }
%                      (same as twod_bilinear)
%
%  Copyright (c) 2002, Jeff Borggaard, Virginia Tech
%  Version: 1.0
%
%  Usage:    M = threed_bilinear(kernel, phi, test, w_g)
%
%  Variables:     kernel
%                        Kernel function in the integral evaluated
%                        at the Gauss points
%
%                 phi
%                        matrix of element test functions evaluated
%                        at the Gauss points (dim: n_gauss, n_dof)
%
%                 test
%                        matrix of test functions evaluated at the
%                        Gauss points (dim: n_gauss, n_dof)
%
%                 w_g
%                        Row vector of Gauss weights
%-----------------------------------------------------------------------
%   [n_gauss,n_row] = size(test);
%   [n_g1   ,n_col] = size(phi );
%
%   M = zeros(n_row,n_col);
%   for i=1:n_row
%     for j=1:n_col
%        M(i,j) = ( w_g'    .* test(:,i)' ) * ( kernel .* phi(:,j) );
%     end
%   end

  M = test'*diag( w_g.*kernel )*phi;


text:
function G = kronMatrix( varargin )

%
% G = kronMatrix( varargin )
%
% Constructor for the kronMatrix class.  This class is used to
% represent objects of the form:
%       \sum [A{i} (x) B{i}]
%
% Object's fields are:
%    G.a = cell array containing A{i} if i>1.
%    G.b = cell array containing B{i} if i>1.
%
% Calling Syntax:
%       G = kronMatrix
%       G = kronMatrix(G);
%       G = kronMatrix(A, B);
%
%   where
%       G = a pre-existing kronMatrix object
%       A, B = both matrices or both cell arrays (same length) such that
%              \sum [A{i} (x) B{i}] represents a bigger matrix.
%

%  L. Perrone, 5/13/03

%  Modifications:
%  6/4/03, J. Nagy
%          Cosmetic changes to incorporate into RestoreTools
%

%
% build the kronMatrix based on number and type of input arguments
%
switch nargin

  case 0
   G.a = cell(0);
   G.b = cell(0);
   G = class(G, 'kronMatrix');

  case 1
   if isa(varargin{1}, 'kronMatrix')
      G = varargin{1};
   else
      error('Wrong input argument')
   end

  case 2

   if isa(varargin{1}, 'double') & isa(varargin{2}, 'double')
     A = cell(1);, A{1} = varargin{1};
     B = cell(1);, B{1} = varargin{2};
     G.a = A;
     G.b = B;
     G = class(G, 'kronMatrix');

   elseif isa(varargin{1}, 'cell') & isa(varargin{2}, 'cell')
     if length(varargin{1}) == length(varargin{2})
        G.a = varargin{1};
        G.b = varargin{2};
        G = class(G, 'kronMatrix');
     else
        error('Input cell arrays must have the same length')
     end

   else
     error('Wrong input arguments')

   end % case 2

  otherwise
    error('Wrong input arguments')

end % switch nargin

repo: IRtools


text:
function k = indexardKernDiagCompute(kern, x)

% INDEXARDKERNDIAGCOMPUTE Compute diagonal of INDEXARD kernel.
% FORMAT
% DESC computes the diagonal of the kernel matrix for the index ard based covariance function kernel given a design matrix of inputs.
% ARG kern : the kernel structure for which the matrix is computed.
% ARG x : input data matrix in the form of a design matrix.
% RETURN k : a vector containing the diagonal of the kernel matrix
% computed at the given points.
%
% SEEALSO : indexKernParamInit, kernDiagCompute, kernCreate, indexKernCompute
%
% COPYRIGHT : Neil D. Lawrence, 2011

% KERN
  k = zeros(size(x, 1), 1);
  for i = 1:length(kern.indices)
    ind = find(round(x)==kern.indices(i));
    k(ind) = kern.indexScales(i);
  end
end

author: SheffieldML

repo: GPmat

save_path: github-repos/MATLAB/SheffieldML-GPmat


text:
function triangle_lyness_rule_test01 ( )

%*****************************************************************************80
%
%% TRIANGLE_LYNESS_RULE_TEST01 tests LYNESS_RULE_NUM, LYNESS_DEGREE, LYNESS_ORDER_NUM.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    30 September 2010
%
%  Author:
%
%    John Burkardt
%
  fprintf ( 1, '\n' );
  fprintf ( 1, 'TRIANGLE_LYNESS_RULE_TEST01\n' );
  fprintf ( 1, '  LYNESS_RULE_NUM returns the number of rules;\n' );
  fprintf ( 1, '  LYNESS_DEGREE returns the degree of a rule;\n' );
  fprintf ( 1, '  LYNESS_ORDER_NUM returns the order of a rule.\n' );

  rule_num = lyness_rule_num ( );

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

  for rule = 0 : rule_num
    order = lyness_order ( rule );
    precision = lyness_precision ( rule );
    fprintf ( 1, '  %8d  %8d  %8d\n', rule, order, precision );
  end

  return
end


text:
function writeMovie(path,vid)
% *************************************************************************
% Superresolution with Dictionary Technique
% writeMovie
%
% Saves an input video as an uncompressed grayscale video
%
% Version 1.0
%
% Created by:   Armin Kappeler
% Date:         03/27/2013
%
% Modifications:
%
% *************************************************************************
property = 'Uncompressed AVI';

writerObj = VideoWriter(path,property);
open(writerObj);

for k = 1:length(vid)
   writeVideo(writerObj,im2uint8(vid{k}));
end

% Create AVI file.
close(writerObj);

end

author: thomas-koehler

repo: SupER


text:
function[mu,sigma]=GMM_parameter(image,segmentation,class_number)
[n,d]=size(image);
mu=zeros(class_number,d);
sigma=zeros(d,d,class_number);
   for i=1:class_number
       Im_i=image(segmentation==i,:);
       [sigma(:,:,i),mu(i,:)]=covmatrix(Im_i);
    end
end

author: Sable

repo: mcbench-benchmarks


text:
%OPENFG Framegrabber device open function.
%
%SYNTAX devinfo = openfg('alias')
%
%DESCRIPTION This function will let the user open a DT
%framegrabber device. The user has to give the
%framegrabbers alias as input argument. The user
%already knows this alias. This function can also be used
%to get the device information structure of a device that
%was already opened.
%
%INPUT The alias (alias) of the framegrabber device.
%
%REMARKS It will be verified if the input is a valid string.
%If it's valid, it will try to open the framegrabber device
%with that specific alias. If more than one Framegrabber is
%installed in the computer, it is necessary to give unique
%aliases to these devices.
%
%OUTPUT If the input isn't valid or the framegrabber device
%can't be opened for some reason, some error message will be
%shown to the user. If everything's ok, this function will
%return a structure with the device information (devinfo) of
%the opened framegrabber device.
%
%EXAMPLE To open the framegrabber with alias 'DT3152' and
%store the device info strucutre in the variable m, or, if
%the device with alias 'DT3152' was already open, to store
%the device info structure in the variable m, use:
%
%	m = openfg('DT3152');

disp('Error: openfg not found.')


text:
function K = symbol_string(kern,dat1,dat2,ind1,ind2,kerParam)
% --- conversion to sequences------------------
D= []; Xf = get_x(dat1); Xs = get_x(dat2); X1 = {}; X2 = {};
for i = 1:size(Xf,1), tmp = Xf(i,:);  X1{i} = tmp(tmp > 0); end
for i = 1:size(Xs,1),  tmp = Xs(i,:);  X2{i} = tmp(tmp > 0);end

n = kerParam{1};
lambda = kerParam{2};

for i = 1:size(X1,1)
    for j = i:size(X2,1)
       K(i,j) =  k(X1{i},X2{j},n,lambda)
    end
end

%------------------------------------------------------
function kret = k(s,t,n,lambda)
        if min(length(s),length(t)) < n
            kret = 0; return
        else
            kret = k(s(1:end-1),t,n,lambda);
            ind = find(t == s(end));
            for i = 1:length(ind)
                kret = kret + k_prime();
            end
        end

%----------------------------------------------------------
function ret = k_prime(s,t,lambda,i)
    if n == 0
        ret = 1; return;
    elseif min(length(s),length(t)) < i
        ret = 0; return;
    end

author: canlab

repo: CanlabCore


text:
function d=ddiff(d1,d2), d=max(d1,-d2);

%   Copyright (C) 2004-2006 Per-Olof Persson. See COPYRIGHT.TXT for details.

repo: ifem