text:
function varargout = spm_update(update)
% Check (and install) SPM updates from the FIL server
% FORMAT spm_update
% This function will contact the FIL server, compare the version number of
% the updates with the one of the SPM installation currently in the MATLAB
% path and display the outcome.
%
% FORMAT spm_update(update)
% Invoking this function with any input parameter will do the same as
% above but will also attempt to download and install the updates.
% Note that it will close any open window and clear the workspace.
%
% FORMAT [sts, msg] = spm_update(update)
% sts  - status code:
%        NaN - SPM server not accessible
%        Inf - no updates available
%        0   - SPM installation up to date
%        n   - new revision <n> is available for download
% msg  - string describing outcome, that would otherwise be displayed.
%__________________________________________________________________________
% Copyright (C) 2010-2020 Wellcome Trust Centre for Neuroimaging

% Guillaume Flandin
% $Id: spm_update.m 7770 2020-01-13 11:26:14Z guillaume $

vspm = spm('Ver');
url  = ['http://www.fil.ion.ucl.ac.uk/spm/download/' lower(vspm) '_updates/'];

if ~nargin
    update = false;
else
    update = true;
end

%-Get list of updates from SPM server
%--------------------------------------------------------------------------
[s,status] = urlread(url);
if ~status
    sts = NaN;
    msg = 'Cannot access SPM server.';
    if ~nargout, error(msg); else varargout = {sts, msg}; end
    return
end

%-Get revision number of latest update
%--------------------------------------------------------------------------
n = regexp(s,[lower(vspm) '_updates_r(\d.*?)\.zip'],'tokens','once');
if isempty(n)
    sts = Inf;
    msg = 'There are no updates available yet.';
    if ~nargout, fprintf([blanks(9) msg '\n']);
    else varargout = {sts, msg}; end
    return
end
n = str2double(n{1});

%-Get revision number of SPM installation
%--------------------------------------------------------------------------
try
    [v,r] = spm('Ver','',1); r = str2double(r);
catch
    error('SPM cannot be found in MATLAB path.');
end
if ~strcmp(v,vspm), error('Your SPM version is %s and not %s',v,vspm); end
rs = [6225 6470 6685 6906 7219 7487 7771];
if isnan(r), r = rs(1); end
if floor(r) == str2double(vspm(4:end))
    try
        r = rs(round((r-floor(r))*10)+1);
    catch
        r = rs(end);
    end
end

%-Compare versions
%--------------------------------------------------------------------------
if n > r
    sts = n;
    msg = sprintf('         A new version of %s is available on:\n',vspm);
    msg = [msg sprintf('   %s\n',url)];
    msg = [msg sprintf('        (Your version: %d - New version: %d)\n',r,n)];
    if ~nargout, fprintf(msg); else varargout = {sts, msg}; end
else
    sts = 0;
    msg1 = sprintf('Your version of %s is up to date.',vspm);
    msg2 = sprintf('(Your version: %d - Online version: %d)',r,n);
    if ~nargout, fprintf([blanks(9) msg1 '\n' blanks(5) msg2 '\n']);
    else varargout = {sts, sprintf('%s\n%s',msg1,msg2)}; end
    return
end

%-and install...
%--------------------------------------------------------------------------
if update
    d = spm('Dir');
    delete(get(0,'Children')); spm('clean'); evalc('spm_rmpath'); drawnow
    try
        lastwarn('');
        m = '          Download and install in progress...\n';
        if ~nargout, fprintf(m); else varargout = {sts, [msg m]}; end
        s = unzip([url sprintf('%s_updates_r%d.zip',lower(vspm),n)], d);
        m = sprintf('         Success: %d files have been updated.\n',numel(s));
        if ~nargout, fprintf(m); else varargout = {sts, [msg m]}; end
    catch
        le = lasterror;
        switch le.identifier
            case 'MATLAB:checkfilename:urlwriteError'
                fprintf('          Update failed: cannot download update file.\n');
            otherwise
                fprintf('\n%s\n',le.message);
        end
    end
    [warnmsg, msgid] = lastwarn;
    switch msgid
        case ''
        case 'MATLAB:extractArchive:unableToCreate'
            fprintf('          Update failed: check folder permission.\n');
        case 'MATLAB:extractArchive:unableToOverwrite'
            fprintf('          Update failed: check file permissions.\n');
        otherwise
            fprintf('          Update failed: %s.\n',warnmsg);
    end
    addpath(d);
    rehash toolboxcache;
end

author: spm


text:
%--- help for splanar.rediff ---
%
%  differentiate - differentiates vectors of splanar objects
%
%  ::
%
%
%    derivs=splanar.differentiate(eqtns,nwrt)
%    derivs=splanar.differentiate(eqtns,nwrt,order)
%    derivs=splanar.differentiate(eqtns,nwrt,order,verbose)
%
%  Args:
%
%     - **eqtns** [splanar|cell array]: vector or cell array of splanar objects
%
%     - **nwrt** [integer]: number of variables for which differentiation is
%       taken
%
%     - **order** [integer|{1}]: order of differentiation
%
%     - **alien_list** [empty|cellstr|char]: list of alien functions
%       (functions that RISE does not recognize and that are to be
%       differentiated by the user himself).
%
%     - **verbose** [true|{false}]: displays information about the process e.g.
%       the amount of time it takes to differentiate each order
%
%  Returns:
%     :
%
%     - **derivs** [structure]: each element in the structure contains:
%       - **size** [2-colum vector]: number of rows and number of columns of
%         the compacted derivatives i.e. unique columns
%       - **derivatives** [vector of splanar]: derivatives
%       - **nwrt** [integer]: number of variables for which differentiation is
%         taken
%       - **order** [integer]: order of differentiation of the current
%         structure
%       - **nnz_derivs** [integer]: number of non-zero derivatives
%       - **partitions** [vector]: vector that help reconstruct the
%         expanded/grand derivative matrix.
%       - **map** []: empty field that will be used in a later stage
%
%  Note:
%
%  Example:
%
%     See also:
%

author: jmaih


text:
function l = legend( varargin )
%legend  Create a legend and make it match the background

%   Copyright 2011 The MathWorks, Inc.

l = legend( varargin{:} );
figh = ancestor( l, 'figure' );
set( l, 'Color', get( figh, 'Color' ) );

repo: Aero-matlab


text:
function X = matrandnorm(varargin)
%MATRANDNORM Normalizes columns of X so that each is unit 2-norm.
%
%   X = MATRANDNORM(M,N) creates a random M x N matrix with randomly using
%   normally distributed enries and then rescales the columsn so that each
%   has a unit 2-norm.
%
%   X = MATRANDNORM(X) rescales the columns of X so that each
%   column has a unit 2-norm.
%
%   Examples
%      X = MATRANDNORM(rand(5,5));
%      X = MATRANDNORM(3,2);
%      X = MATRANDNORM(ones(4));
%
%   See also MATRANDORTH, MATRANDNORM, CREATE_PROBLEM, CREATE_GUESS.
%
%MATLAB Tensor Toolbox.
%Copyright 2015, Sandia Corporation.

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

if nargin == 2
    X = randn(varargin{1}, varargin{2});
else
    X = varargin{1};
end

norms = sqrt(sum(X.^2,1));
X = bsxfun(@rdivide,X,norms);

author: zhangqianqianQQ

repo: MachineVisionAlgorithm


text:
%JC 7/16/08-sim_qpsk_video
%I thought some might be interested in this m-file, since we live in a
%world of data,text,audio and video.
%User has the option of using Data or an image(Lena1.bmp) to determine BER/SER for
%M=2,4 or 8 PSK.
%User also has the option of using Gray coding or no Gray coding.
%Viewing a still image file thru a channel with AWGN added gives an intuitive insight into the
%degradation of video at different levels of SNR. Running the image takes
%approx 8 minutes. The program showes, for (M=2,4), the image is reproduced
%(without errors) at approx. 10dB and M=8 requires approx 14 dB which
%agrees with theory. The bandwidth efficiency for M=2 is 1(bit/sec/Hz), M=4
%is 2(bits/sec/Hz), and M=8 is 3(bits/sec/Hz). Therefore M=4 and M=8 are good
%choices for band limited channels. I restructered this program from one
%written by
%Chen Zhifeng
%4/28/2007
%Search WEB for "Performance Analysis of Channel Estimation and
%Equalization in Slow Fading Channels"
%--------------------------------------------------------------------------
close all;
clear all;
clc;
%--------------------------------------------------------------------------
%Set parameters
%--------------------------------------------------------------------------
%randn('state',0);%keeps noise characteristics the same on reruns
%rand('state',0)%keeps intergers same on rerun
%Holding all values constant, shows a slight improvement on BER @minus5dB when using Gray
%coding and shows a slight degredation @plus 9dB. Is this valid according to theory?
%An interesting test would be to see the results with FEC(with and without Gray coding) if
%you have the Communications Toolbox(I don't) and also to test these
%results against it. Also,
%I challenge someone to write an easy to understand m-file Viterbi decoder
%that can be used in conjunction with the convolutional encoder m-file
%under JC files in author index. This would be helpful to many folks,
%including myself.
Max_dB =10  % input SNR in dB
modtype = 'psk'
M=8                    %M-ary
gray_encode = 1;      %1=yes 0=no
Ndata = 1000000;      %limit to 1000000 samples for matlab processing
Test_image = 1;       %1=yes 0=no
Image_name = 'Lena1.bmp' %color photo Lena1  im 128x128x3 uint8 array
                          %M=4   128x128x3x4=196608
k=log2(M);
ebn0=10.^(Max_dB/10);
%--------------------------------------------------------------------------
%produce Data or load image
%--------------------------------------------------------------------------
%use this method because the program may transmit data from file, such as a image
Data=My_randint(1,Ndata,0,3);%produce intergers from 0 to 3(QPSK-M=4)-Change for M=2,8
%Data=My_randint(1,Ndata,0,1);%M=8
%Data=My_randint(1,Ndata,0,1);%M=2
%use Data=randint(1,Ndata,M); if you have this function(Communications Toolbox)
if Test_image
    im = imread('Lena1.bmp');
    disp('Fig 1 showing original image before AWGN added');
    image(im);
    drawnow
    [Data, row_im, col_im, third_im] = image2data(im, M);
end
Ndata = length(Data);
Transmit = Data;
datalen = length(Transmit);
%--------------------------------------------------------------------------
%encode-Gray or no Gray
%--------------------------------------------------------------------------
if gray_encode
    % Create Gray encoding and decoding arrays-flip(2=3)and(3=2 for M=4)
    grayencod = bitxor(0:M-1, floor((0:M-1)/2));
    [dummy graydecod] = sort(grayencod); graydecod = graydecod - 1;

    % Gray encode symbols
    Transmit_gray = grayencod(Transmit+1);
end
%--------------------------------------------------------------------------
%Modulation
%--------------------------------------------------------------------------
if gray_encode
    step=2*pi/M;
    S=exp(j*Transmit_gray.*step);
else %no Gray encode
    step=2*pi/M;
    S=exp(j*Transmit.*step);
end
%-----------------------------------------------------------------------
%uniform the noise power by 1/var=SNR=Es/(N0/2)=2*k*ebn0
std=(1/2/k./ebn0).^0.5;
%-----------------------------------------------------------------------
%-----------------------------------------------------------------------
%add AWGN
%-----------------------------------------------------------------------
Z = std*( randn(1,datalen)+j*randn(1,datalen) );
R=S+Z;
%-------------------------------------------------------------------------
%demodulation
%-------------------------------------------------------------------------
step=2*pi/M;
Rphase = atan2(imag(R), real(R));%1 radian=57.296 degrees
%Sphase = atan2(imag(S), real(S));%for testing only
Receive = Rphase/step;%error is made if noise causes phase to fall outside +-45 degrees for M=4

Receive = round(Receive);
Receive = mod(Receive,M);
if gray_encode
% Gray decode message
Receive = graydecod(Receive+1);%Received image or data
end

ne=sum(Data~=Receive)%number of errors
BER=ne/Ndata%bit error rate
SER=ne/Ndata/k

if Test_image
   ima = data2image(Receive, row_im, col_im, third_im, M);
   ima = uint8(ima);
   imwrite(ima,'Received_image.bmp');
   disp(' ');
   disp('Fig 2 received image(after AWGN added) is saved as Received_image.bmp in the same directory');
   disp(' ');
   figure
   image(ima);
   drawnow
end

%INFO-The analytical expressions and references used in berawgn, bercoding, berfading, and
%BERTool are found in the Communications Toolbox documentation.

author: Sable

repo: mcbench-benchmarks


text:
function y = logWishart(Sigma, W, v)
% Compute log pdf of a Wishart distribution.
% Input:
%   Sigma: d x d covariance matrix
%   W: d x d covariance parameter
%   v: degree of freedom
% Output:
%   y: probability density in logrithm scale y=log p(Sigma)
% Written by Mo Chen (sth4nth@gmail.com).
d = length(Sigma);
B = -0.5*v*logdet(W)-0.5*v*d*log(2)-logmvgamma(0.5*v,d);
y = B+0.5*(v-d-1)*logdet(Sigma)-0.5*trace(W\Sigma);

author: PRML

repo: PRMLT

save_path: github-repos/MATLAB/PRML-PRMLT


text:
% Data file PILZ1
% Forced damped vibrations of a mechanical
% system with one degree of freedom, contained
% a reverse pendulum.
  Ek   = '3.6666*qt^2'; % Kinetic energy
  N    = '(10*sin(p*t)-c*q-k*qt+19.62*sin(10*q))*qt'; % power
  Tend = 10;    % upper bound of integration
  q0   = '0.1'; % initial coordinate
  qt0  = '0';   % initial velocity
  eps  = 1e-10; % desirable accuracy
  np   = 3;     % number of parameters
  P{1} = 'c';   % spring stiffness
  P{2} = 'k';   % coefficient of damping
  P{3} = 'p';   % disturbance frequency


text:
function s = ClenshawL(c,x)
% Clenshaw algorithmus for evaluating Legendre series
%
% Description
%
% $$ s = \sum_n=0^N c_n P_n(x)$$
%
% Input
%  c - coefficients
%  x - evaluation nodes
%
% Output
%  s - Legendre series
%

% init
c= [c(:);0;0];
dn = zeros(size(x));
d1 = zeros(size(x));
d2 = zeros(size(x));

% use backward tree term recurence
for l = length(c)-2:-1:1

  d1 = d2 + (2*l+1)/(l+1) * x .* dn;
  d2 = c(l) - l/(l+1) * dn;
  dn = d1;

end

s = d2 + x .* d1;

end

author: mtex-toolbox

repo: mtex

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


text:
function ft_realtime_downsample(cfg)

% FT_REALTIME_DOWNSAMPLE reads realtime data from one buffer and writes it after downsampling
% to another buffer.
%
% Use as
%   ft_realtime_downsample(cfg)
% with the following configuration options
%   cfg.channel              = cell-array, see FT_CHANNELSELECTION (default = 'all')
%   cfg.decimation           = integer, downsampling factor (default = 1, no downsampling)
%   cfg.order                = interger, order of butterworth lowpass filter (default = 4)
%   cfg.cutoff               = double, cutoff frequency of lowpass filter (default = 0.8*Nyquist-freq.)
%
% The source of the data is configured as
%   cfg.source.dataset       = string
% or alternatively to obtain more low-level control as
%   cfg.source.datafile      = string
%   cfg.source.headerfile    = string
%   cfg.source.eventfile     = string
%   cfg.source.dataformat    = string, default is determined automatic
%   cfg.source.headerformat  = string, default is determined automatic
%   cfg.source.eventformat   = string, default is determined automatic
%
% The target to write the data to is configured as
%   cfg.target.datafile      = string, target destination for the data (default = 'buffer://localhost:1972')
%   cfg.target.dataformat    = string, default is determined automatic
%
% To stop this realtime function, you have to press Ctrl-C

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

% set the defaults
if ~isfield(cfg, 'source'),               cfg.source = [];                                  end
if ~isfield(cfg, 'target'),               cfg.target = [];                                  end
if ~isfield(cfg.source, 'headerformat'),  cfg.source.headerformat = [];                     end % default is detected automatically
if ~isfield(cfg.source, 'dataformat'),    cfg.source.dataformat = [];                       end % default is detected automatically
if ~isfield(cfg.target, 'headerformat'),  cfg.target.headerformat = [];                     end % default is detected automatically
if ~isfield(cfg.target, 'dataformat'),    cfg.target.dataformat = [];                       end % default is detected automatically
if ~isfield(cfg.target, 'datafile'),      cfg.target.datafile = 'buffer://localhost:1972';  end
if ~isfield(cfg, 'minblocksize'),         cfg.minblocksize = 0;                             end % in seconds
if ~isfield(cfg, 'maxblocksize'),         cfg.maxblocksize = 1;                             end % in seconds
if ~isfield(cfg, 'channel'),              cfg.channel = 'all';                              end

if ~isfield(cfg, 'decimation'),           cfg.decimation = 1;                               end
if ~isfield(cfg, 'order'),                cfg.order = 4;                                    end

if cfg.decimation < 1 || cfg.decimation~=round(cfg.decimation)
  error 'Decimation factor must be integer quantity >= 1';
end

% translate dataset into datafile+headerfile
cfg.source = ft_checkconfig(cfg.source, 'dataset2files', 'yes');
cfg.target = ft_checkconfig(cfg.target, 'dataset2files', 'yes');
ft_checkconfig(cfg.source, 'required', {'datafile' 'headerfile'});
ft_checkconfig(cfg.target, 'required', {'datafile' 'headerfile'});

% ensure that the persistent variables related to caching are cleared
clear ft_read_header
% read the header for the first time
hdr = ft_read_header(cfg.source.headerfile);
fprintf('updating the header information, %d samples available\n', hdr.nSamples*hdr.nTrials);

targethdr = hdr;
targethdr.Fs = targethdr.Fs/cfg.decimation;

% Calculate cutoff frequency expressed in Nyquist freq of original signal
if ~isfield(cfg, 'cutoff')
  cutoff = 0.8/cfg.decimation;
else
  cutoff = cfg.cutoff / (0.5*hdr.Fs);
  if cutoff >= 1
    ft_error('Cutoff frequency too high');
  end
end
[B,A] = butter(cfg.order, cutoff, 'low');

% define a subset of channels for reading
cfg.channel = ft_channelselection(cfg.channel, hdr.label);
chanindx    = match_str(hdr.label, cfg.channel);
nchan       = length(chanindx);
if nchan==0
  ft_error('no channels were selected');
end

minblocksmp = round(cfg.minblocksize*hdr.Fs);
minblocksmp = max(minblocksmp, 1);
maxblocksmp = round(cfg.maxblocksize*hdr.Fs);
prevSample  = 0;
count       = 0;

% create filter model and downsampling model
FM = online_filter_init(B, A, zeros(nchan, 1));
DM = online_downsample_init(cfg.decimation);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% this is the general BCI loop where realtime incoming data is handled
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
while true

  % determine number of samples available in buffer
  hdr = ft_read_header(cfg.source.headerfile, 'cache', true);

  % see whether new samples are available
  newsamples = (hdr.nSamples*hdr.nTrials-prevSample);

  if newsamples>=minblocksmp

    % determine the samples to copy from source to target
    begsample = prevSample+1;
    endsample = prevSample + min(newsamples, maxblocksmp);

    % remember up to where the data was read
    prevSample  = endsample;
    count       = count + 1;
    fprintf('processing segment %d from sample %d to %d\n', count, begsample, endsample);

    % read data segment
    dat = ft_read_data(cfg.source.datafile, 'dataformat', cfg.source.dataformat, 'begsample', begsample, 'endsample', endsample, 'chanindx', chanindx, 'checkboundary', false);

    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    % from here onward it is specific to the downsampling
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

    [FM, dat] = online_filter_apply(FM, dat);
    [DM, xd] = online_downsample_apply(DM, dat);

    if count==1
      % flush the file, write the header and subsequently write the data segment
      ft_write_data(cfg.target.datafile, dat, 'header', targethdr, 'dataformat', cfg.target.dataformat, 'chanindx', chanindx, 'append', false);
    else
      % write the data segment
      ft_write_data(cfg.target.datafile, dat, 'header', targethdr, 'dataformat', cfg.target.dataformat, 'chanindx', chanindx, 'append', true);
    end

  end % if enough new samples
end % while true

author: fieldtrip

repo: fieldtrip

save_path: github-repos/MATLAB/fieldtrip-fieldtrip


text:
function  [PCA_idx, s_idx, seg]  =  Set_PCA_idx( im, par, codewords )

[h  w]     =  size(im);
cls_num    =  size( codewords, 2 );
LP_filter  =  fspecial('gaussian', 9, par.sigma);
lp_im      =  conv2( LP_filter, im );
lp_im      =  lp_im(4:h+3, 4:w+3);
hp_im      =  im - lp_im;

b       =  par.win;
s       =  par.step;
N       =  h-b+1;
M       =  w-b+1;

r       =  [1:s:N];
r       =  [r r(end)+1:N];
c       =  [1:s:M];
c       =  [c c(end)+1:M];
L       =  length(r)*length(c);
X       =  zeros(b*b, L, 'single');

% For the Y component
k    =  0;
for i  = 1:b
    for j  = 1:b
        k    =  k+1;
        blk  =  hp_im(r-1+i,c-1+j);
        X(k,:) =  blk(:)';
    end
end
PCA_idx   =  zeros(L, 1);

m     =  mean(X);
d     =  (X-repmat(m, size(X,1), 1)).^2;
v     =  sqrt( mean( d ) );

delta       =  sqrt(par.nSig^2+25);
[a, ind]    =  find( v<delta );

set         =  [1:L];
set(ind)    =  [];
L2          =  size(set,2);

for i = 1:L2

    wx            =   repmat( X(:,set(i)), 1, cls_num );
    dis           =   sum( (wx - codewords).^2 );
    [md, idx]     =   min(dis);
    PCA_idx( set(i) )  =   idx;

end

[s_idx, seg] =   Proc_cls_idx( PCA_idx );

author: lbasek

repo: image-denoising-benchmark


text:
function buzsakilab_database_submit(bz_datasets,database_name)
% Submit dataset to the Buzsakilab metadata-database
%
% v1.1
%
% INPUT
% bz_datasets :
% A structure with size n x 1, where n is number of sessions to submit.
%
% database_name :
% Name of database to submit to. Default 'buzsakid_submitdataset'
%
%
% By Peter Petersen
% petersen.peter@gmail.com

if nargin == 1
    database_name = 'buzsakid_submitdataset';
end

fprintf(['\nSubmitting datasets to database: ''', database_name, '''\n'])

% Make connection to database
bz_mySQL = bz_database_credentials; % Credentials loaded

conn = database(database_name,bz_mySQL.userid,bz_mySQL.password,'Vendor','MYSQL','Server',bz_mySQL.address,'PortNumber',3306);
curs = exec(conn,['SELECT * FROM ' database_name '.' bz_mySQL.table]);
curs = fetch(curs);

if ~strcmp(curs.Message,'Invalid connection.')
    bz_database = table2struct(curs.Data);
    % attributes = attr(curs);
    close(curs);

    % Execute query and fetch results
    if isfield(bz_datasets, 'Id')
        bz_datasets = rmfield(bz_datasets, 'Id');
    end

    for k = 1:size(bz_datasets,1)
        if sum(strcmp(bz_datasets(k).Session,{bz_database.Session}) & strcmp(bz_datasets(k).Animal,{bz_database.Animal}) & strcmp(bz_datasets(k).Investigator,{bz_database.Investigator}))
            bz_datasets(k).NewDataset = 0;
        else
            bz_datasets(k).NewDataset = 1;
        end
    end
    clear bz_database;
    n_datasets = size(bz_datasets,1);
    bz_datasets = bz_datasets(find([bz_datasets.NewDataset]))';

    if size(bz_datasets,2) > 0
        bz_datasets = rmfield(bz_datasets,{'NewDataset'});
        for j = 1:length(bz_datasets)
            % Checking for existing datasets in the database
            disp([num2str(j) '. dataset: ' bz_datasets(j).Investigator, ', ', bz_datasets(j).Animal,', ', bz_datasets(j).Session])
            if any(structfun(@isempty, bz_datasets(j)))
                idx = find(structfun(@isempty, bz_datasets(j)));
                names = fieldnames(bz_datasets);
                for field = 1:length(idx)
                    warning([names{idx(field)} ' was empty, filling with null.'])
                    bz_datasets(j) = setfield(bz_datasets(j),names{idx(field)},'null');
                end
            end
            % Data to insert
            colnames = fieldnames(bz_datasets(j))'; %
            data_table = struct2table(bz_datasets(j));
            % data_table = buzsakilab_database_verify_fieldtypes_and_dataformat(data_table,attributes)

            insert(conn,bz_mySQL.table,colnames,data_table);
        end
        disp([num2str(size(bz_datasets,2)), '/' num2str(n_datasets), ' dataset(s) submitted to database. ', num2str(n_datasets-size(bz_datasets,2)) ' dataset(s) already exist in database'])
    else
        disp('Datasets already exist in database')
    end
else
    disp('Connection failed. mySQL server not available')
    disp(curs.Message)
end
close(conn) % Close connection to database
clear conn curs % Clear variables

author: buzsakilab

repo: buzcode

save_path: github-repos/MATLAB/buzsakilab-buzcode


text:
function [DataMat, ChannelMat] = in_data_rhd( DataFile )
% IN_DATA_RHD: Read a RHD file.
%
% INPUT:
%    - DataFile : Full path to a recordings file (called 'data' files in Brainstorm)
% OUTPUT:
%    - DataMat : Brainstorm standard recordings ('data') structure

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

%% ===== GET OPTIONS =====
% Get file info
[sFile, ChannelMat] = in_fopen_intan(DataFile); % This call gets all of the file info (channels, events etc.)

%% Now Read entire recording and assign to DataMat
% Initialize returned structures
DataMat = db_template('DataMat');
DataMat.Comment  = 'EEG/Intan';
DataMat.DataType = 'recordings';
DataMat.Device   = 'Intan';

isProgress = bst_progress('isVisible');
if ~isProgress
    bst_progress('start', 'Intan importer', 'Reading entire recording');
end
% Read entire file
DataMat.F = in_fread_intan(sFile, [], [], []);

% Add the events
DataMat.Events = sFile.events;

% Build time vector
DataMat.Time = sFile.prop.times(1):1/sFile.prop.sfreq:sFile.prop.times(2);

% ChannelFlag
DataMat.ChannelFlag = ones(size(DataMat.F,1), 1);

% Replace NaN with 0
DataMat.F(isnan(DataMat.F)) = 0;

% Save number of trials averaged
DataMat.nAvg = 1;

% Build comment
[fPath, fBase, fExt] = bst_fileparts(DataFile);
DataMat.Comment = fBase;

isProgress = bst_progress('isVisible');
if isProgress
    bst_progress('stop');
end
end

author: brainstorm-tools


text:
%SETVERSION Set the version field of a dataset
%
%   A = SETVERSION(A,S)
%
% INPUT
%   A      Dataset
%   S      String to be stored in the version field
%
% OUTPUT
%   A      Updated dataset
%
% DESCRIPTION
% Set or reset the version field of A by S.

author: marianux

repo: ecg-kit

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


text:
function R = rotation_matrix(phi,dim,orientation)
%ROTATION_MATRIX 3D rotation matrix
%
%   Usage: R = rotation_matrix(phi,[dim,[orientation]])
%
%   Input parameters:
%       phi         - angle to rotate the given dim dimension vector / rad
%       dim         - dimension to turn around (default: 3, z-axis)
%       orientation - orientation of the rotation, 'clockwise' or
%                     'counterclockwise' (default: 'counterclockwise')
%
%   Output parameters:
%       R       - 3x3 rotation matrix to apply to your vector to
%                 rotate: R*y
%
%
%   ROTATION_MATRIX(phi,dimension,orientation) returns a rotation matrix R,
%   which is able to rotate a vector around the given dimension about phi.
%
%   See also: sin, cos

%*****************************************************************************
% The MIT License (MIT)                                                      *
%                                                                            *
% Copyright (c) 2010-2019 SFS Toolbox Developers                             *
%                                                                            *
% 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.                                                  *
%                                                                            *
% The SFS Toolbox  allows to simulate and  investigate sound field synthesis *
% methods like wave field synthesis or higher order ambisonics.              *
%                                                                            *
% https://sfs.readthedocs.io                            sfstoolbox@gmail.com *
%*****************************************************************************

%% ===== Checking of input  parameters ==================================
nargmin = 1;
nargmax = 3;
narginchk(nargmin,nargmax);
isargscalar(phi)
if nargin<3
    % Set defualt orientation of the rotation
    orientation = 'counterclockwise';
end
if nargin<2
    % set default rotation dimension to z-axis
    dim = 3;
end
isargchar(orientation);
isargpositivescalar(dim);

%% ===== Computation ====================================================
% Rotation matrix (see: https://en.wikipedia.org/wiki/Rotation_matrix)
% get single matrix entries
r1 = cos(phi);
r4 = cos(phi);
if strcmp('counterclockwise',orientation)
    r2 = -sin(phi);
    r3 =  sin(phi);
elseif strcmp('clockwise',orientation)
    r2 =  sin(phi);
    r3 = -sin(phi);
else
    error('%s: the given orientation "%s" is not known.', ...
        upper(mfilename),orientation);
end
% Fill up matrix to rotate around the given axis
if dim==1

    R = [1 0  0;  ...
         0 r1 r2; ...
         0 r3 r4];
elseif dim==2
    R = [r1 0 r2; ...
         0  1 0;  ...
         r3 0 r4];
elseif dim==3
    R = [r1 r2 0; ...
         r3 r4 0; ...
         0  0  1];
else
    error('%s: dim has to be 1,2, or 3 and not %i',upper(mfilename),dim);
end

author: sfstoolbox

repo: sfs-matlab


text:
function u=EOBurg(u,t,dx,dir,bndcond)
%%
% Solves Burgers' equation u_t + 0.5*u^2_x = 0 by the Engquist-Osher method.
%

%%
% Initial setup
if nargin<5,
  bndcond='Neumann';
end;
if nargin<4,
  dir=1;
end;
if (dir>2),
  error(strcat('You can easily extend this to',num2str(dir)','D yourself ...'));
end;
% Some things ensuring that you solve along the columns of u. Must
% transpose if u is a row vector or if dir = 2.
transpose=0;
S=size(u);
N=S(1);
if ((N==1)||(dir==2)),
  u=u';
  transpose=1;
end;
maxspeed=max(max(abs(u)));
CFL=0.95;
dt=CFL*dx/maxspeed;
nstep=ceil(t/dt);
dt=t/nstep;
lambda=dt/dx;
%%
% Solving by EO metod
for j=1:nstep,
  uext=set_extend('iden',u,bndcond);
  u=u-lambda*diff(eoflux(uext));
end;
if (transpose)
  u=u';
end;

%% eoflux
function f=eoflux(u)
  S=size(u);
  N=S(1);
  MM=max(u(1:N-1,:),0.0); mm=min(u(2:N,:),0.0);
  f=0.5*(MM.^2+mm.^2);


text:
function example4 ( )

%*****************************************************************************80
%
%% EXAMPLE4 demonstrates the use of PDEPE on a variable coefficient PDE system.
%
%  Discussion:
%
%    Solve the system of convection-diffusion equations.
%
%    PDE:
%      ut - 2 ux -   vx = uxx
%      vt -   ux - 2 vx - ( 3 (ubar(x))^2 u ) = vxx
%    BC:
%      u(t,-oo) = u(t,+oo) = 0
%      v(t,-oo) = v(t,+oo) = 0
%    IC:
%      u(0,x) = exp ( - (x-5)^2 )
%      v(0,x) = exp ( - (x+5)^2 )
%
%    Here, ubar(x) is the first component in the solution of the boundary value
%    problem:
%
%      ubarx = - 2 ( ubar + 2 ) - vbar
%      vbarx = - ( ubar + 2 ) - 2 vbar - ( ubar^3 + 8 )
%      ubar(-oo) = -2, ubar(+oo) = 1
%      vbar(-oo) = 0, vbar(+00) = -6
%
%    Although the mathematical problem has boundary conditions at infinity,
%    we approximate this by the interval [-25,+25].
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    06 September 2013
%
%  Author:
%
%    Original formulation by P Howard.
%    This version by John Burkardt.
%
%  Reference:
%
%    P Howard,
%    Partial Differential Equations in Matlab 7.0,
%    Spring 2005.
%
  timestamp ( );
  fprintf ( 1, '\n' );
  fprintf ( 1, 'EXAMPLE4:\n' );
  fprintf ( 1, '  A convection-diffusion system:\n' );
  fprintf ( 1, '  ut - 2 ux -   vx = uxx\n' );
  fprintf ( 1, '  vt -   ux - 2 vx - ( 3 (ubar(x))^2 u ) = vxx\n' );
  fprintf ( 1, '  u(t,-oo) = u(t,+oo) = 0\n' );
  fprintf ( 1, '  v(t,-oo) = v(t,+oo) = 0\n' );
  fprintf ( 1, '  u(0,x) = exp ( - (x-5)^2 )\n' );
  fprintf ( 1, '  v(0,x) = exp ( - (x+5)^2 )\n' );
%
%  M defines the coordinate system:
%  0: cartesian
%  1: cylindrical
%  2: spherical
%
  m = 0;
%
%  Define the spatial mesh.
%
  nx = 101;
  xmesh = linspace ( -25.0, +25.0, nx );
%
%  Define the time mesh.
%
  nt = 21;
  tspan = linspace ( 0.0, 2.0, nt );
%
%  Call PDEPE() for the solution.
%
  sol = pdepe ( m, @pdefun, @icfun, @bcfun, xmesh, tspan );
%
%  Copy out the two components of the solution.
%
  u = sol(:,:,1);
  v = sol(:,:,2);

  figure ( 1 )
  subplot(2,1,1)
  surf ( xmesh, tspan, u, 'EdgeColor', 'None' );
  title ( 'Example 4: Solution U Over Time', 'Fontsize', 16 );
  xlabel ( '<--- X --->' )
  ylabel ( '<--- T --->' );
  zlabel ( '<---U(X,T)--->' );
  subplot(2,1,2)
  surf ( xmesh, tspan, v, 'EdgeColor', 'None' );
  title ( 'Example 4: Solution V Over Time', 'Fontsize', 16 );
  xlabel ( '<--- X --->' )
  ylabel ( '<--- T --->' );
  zlabel ( '<---V(X,T)--->' );
  filename = 'example4.png';
  print ( '-dpng', filename );
  fprintf ( 1, '\n' );
  fprintf ( 1, '  Saved solution plot in file "%s"\n', filename );
%
%  Plot the initial condition, U at time 0.
%
  figure ( 2 )
  subplot(2,1,1)
  plot ( xmesh, u(1,:), 'LineWidth', 3 );
  grid on
  title ( 'Example 4: Initial Condition U(X,0)', 'Fontsize', 16 );
  xlabel ( '<--- X --->' )
  ylabel ( '<--- U(X,T0) --->' );
  subplot(2,1,2)
  plot ( xmesh, v(1,:), 'LineWidth', 3 );
  grid on
  title ( 'Example 4: Initial Condition V(X,0)', 'Fontsize', 16 );
  xlabel ( '<--- X --->' )
  ylabel ( '<--- V(X,T0) --->' );
  filename = 'example4_ic.png';
  print ( '-dpng', filename );
  fprintf ( 1, '  Saved initial condition plot in file "%s"\n', filename );
%
%  Plot the solution at a fixed point, with time varying.
%
  mid = ( nx + 1 ) / 2;

  figure ( 3 )
  subplot(2,1,1)
  plot ( tspan, u(:,mid), 'LineWidth', 3 );
  grid on
  title ( 'Example 4: Time evolution of U(0,T)', 'Fontsize', 16 );
  xlabel ( '<--- T --->' )
  ylabel ( '<--- U(0.0,T) --->' );
  subplot(2,1,2)
  plot ( tspan, v(:,mid), 'LineWidth', 3 );
  grid on
  title ( 'Example 4: Time evolution of V(0,T)', 'Fontsize', 16 );
  xlabel ( '<--- T --->' )
  ylabel ( '<--- V(0.0,T) --->' );
  filename = 'example4_profile.png';
  print ( '-dpng', filename );
  fprintf ( 1, '  Saved time evolution plot in file "%s"\n', filename );
%
%  Terminate.
%
  fprintf ( 1, '\n' );
  fprintf ( 1, 'EXAMPLE4:\n' );
  fprintf ( 1, '  Normal end of execution.\n' );
  fprintf ( 1, '\n' );
  timestamp ( );

  return
end
function xprime = degode ( t, x )

%*****************************************************************************80
%
%% DEGODE stores an ODE for a standing wave.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    06 September 2013
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, real T, the current time.
%
%    Input, real X, the spatial location.
%
%    Output, real XPRIME(2,1), the right hand side of the differential
%    equation satisfied by ( UBAR, VBAR ).
%
  xprime(1,1) =  - 2.0 * ( x(1) + 2.0 ) - x(2);
  xprime(2,1) = - ( x(1) + 2.0 ) - 2.0 * x(2) - ( x(1)^3 + 8.0 );

  return
end
function ubar = pdegwave ( x )

%*****************************************************************************80
%
%% PDEGWAVE determines the value of UBAR in the equation.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    06 September 2013
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, real X, the spatial location.
%
%    Output, real UBAR, the value of UBAR at X.
%
  small = 0.000001;

  if ( x < -20.0 )
    ubar = -2.0;
    vbar = 0.0;
  else
    tspan = [ -20, x ];
    x0 = [ -2.0 + small, - small ];
    [ t, x ] = ode45 ( @degode, tspan, x0 );
    ubar = x(end,1);
    vbar = x(end,2);
  end

  return
end
function [ c, f, s ] = pdefun ( x, t, u, dudx )

%*****************************************************************************80
%
%% PDEFUN defines the components of the PDE.
%
%  Discussion:
%
%    The PDE has the form:
%
%      c * du/dt = x^(-m) d/dx ( x^m f ) + s
%
%    where m is 0, 1 or 2,
%    c, f and s are functions of x, t, u, and dudx,
%    and most typically, f = dudx.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    06 September 2013
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, real X, the spatial location.
%
%    Input, real T, the current time.
%
%    Input, real U(:,1), the estimated solution at T and X.
%
%    Input, real DUDX(:,1), the estimated spatial derivative of U at T and X.
%
%    Output, real C(:,1), the coefficients of du/dt.
%
%    Output, real F(:,1), the flux terms.
%
%    Output, real S(:,1), the source terms.
%
  c = [ 1.0; 1.0 ];
  ubar = pdegwave ( x );
  f(1) = dudx(1) + 2.0 * u(1) + u(2);
  f(2) = dudx(2) + u(1) + 2.0 * u(2) + 3.0 * ubar^2 * u(1);
  s = [ 0.0; 0.0 ];

  return
end
function u0 = icfun ( x )

%*****************************************************************************80
%
%% ICFUN defines the initial conditions.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    06 September 2013
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, real X, the spatial location.
%
%    Output, real U0(:,1), the value of the solution at the initial time,
%    and location X.
%
  u0(1) = exp ( - ( x - 5.0 )^2 );
  u0(2) = exp ( - ( x + 5.0 )^2 );

  return
end
function [ pl, ql, pr, qr ] = bcfun ( xl, ul, xr, ur, t )

%*****************************************************************************80
%
%% BCFUN defines the boundary conditions.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    06 September 2013
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, real XL, the spatial coordinate of the left boundary.
%
%    Input, real UL(:,1), the solution estimate at the left boundary.
%
%    Input, real XR, the spatial coordinate of the right boundary.
%
%    Input, real UR(:,1), the solution estimate at the right boundary.
%
%    Output, real PL(:,1), the Dirichlet portion of the left boundary condition.
%
%    Output, real QL(:,1), the coefficient of the flux portion of the left
%    boundary condition.
%
%    Output, real PR(:,1), the Dirichlet portion of the right boundary condition.
%
%    Output, real QR(:,1), the coefficient of the flux portion of the right
%    boundary condition.
%
  pl = [ ul(1); ul(2) ];
  ql = [ 0.0; 0.0 ];
  pr = [ ur(1); ur(2) ];
  qr = [ 0.0; 0.0 ];

  return
end
function timestamp ( )

%*****************************************************************************80
%
%% TIMESTAMP prints the current YMDHMS date as a timestamp.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    14 February 2003
%
%  Author:
%
%    John Burkardt
%
  t = now;
  c = datevec ( t );
  s = datestr ( c, 0 );
  fprintf ( 1, '%s\n', s );

  return
end

author: johannesgerer

repo: jburkardt-m


text:
function simOutput = mpcSim(nx, ny, nu, Nsim, e0, umin, umax, H, G, Lineq, Mineq, Leq, Meq, phi, gamma, lambda, chkInputBoundFlag, sivdMode, Xref)
% mpcSim Simulate the system with MPC Controller
%
% Author : Ajinkya Khade, askhade@ncsu.edu

[J, sivd,exitflag] = deal(zeros(1, Nsim));

% Initializing matrices to log simulation data
empc = zeros(nx, Nsim);         % stores error vectors at all time steps
ympc = zeros(ny, Nsim);         % stores output vectors at all time steps
umpc = zeros(nu, Nsim);         % stores calculated input vectors
uact = umpc;                    % stores applied input vectors

% Assigning values at the beginning of simulation
empc(:,1) = e0;
umpc(:,1) = 0;
sivd(1)   = Xref(2);

% Options structure for the QP Solver
options = optimoptions('quadprog','Algorithm','active-set','Display','off','TolFun',0.001);

% Simulating over N time steps
for k = 1:Nsim

    F = (2*G(empc(:,k)))';

    % Solving the QP at current time step
    [U, J(k),exitflag(k)] = quadprog(H, F, Lineq, Mineq(empc(:,k),sivd(k),Xref(2,k)), ...
                                           Leq,   Meq(empc(:,k)), ...
                                           [], [], [], options);

	% Extracting and storing calculated input at current time step
    umpc(:,k)   = U(1:nu);
    % Checking for saturation and storing applied input at current time step
    uact(:,k)   = getActualU(umpc(:,k), umin, umax, chkInputBoundFlag);
	% Calculating and storing error vector at current time step
    empc(:,k+1) = phi*empc(:,k) + gamma*uact(:,k);
	% Calculating and storing ouput at current time step
    ympc(:,k)   = lambda*empc(:,k);

    % Calculating SIVD
    if sivdMode == 0
        % SIVD varies linearly with reference trajectory
        sivd(k+1) = Xref(2,k);
    elseif sivdMode == 1
        % SIVD varies dynamically with range rate
        sivd(k+1) = getsivd(empc(2,k), umin);
    end

    % Update range and range rate to make it relative to SIVD
    empc(1,k+1) = empc(1,k+1) + (sivd(k+1) - sivd(k));
    empc(2,k+1) = empc(2,k+1) - (Xref(2,k+1) - Xref(2,k));
end

empc = empc(:,1:end-1);
sivd = sivd(:,1:end-1);

% Storing all simulation information in a structure for returning
simOutput.empc = empc;
simOutput.ympc = ympc;
simOutput.umpc = umpc;
simOutput.uact = uact;
simOutput.sivd = sivd;
simOutput.exitflag = exitflag;
simOutput.J = J;

end

function uact = getActualU(umpc, umin, umax, uConstraintFlag)
%getActualU  Calculate actual input applied.
%   This function checks if input violates bounds and calculates the actual
%   input to be applied.

uact = umpc;

if uConstraintFlag
    if umpc > umax
        uact = umax;
    elseif umpc < umin
        uact = umin;
    end
end

end

author: ajinkya-khade