text:
function bi = evalBSpline(xi,deg)
% evalBSpline: evaluate BSpline basis function
% usage: bi = evalBSpline(xi,deg);
%
% arguments:
%   xi (N-D arrary) - positions at which to evaluate BSpline basis function
%   deg (scalar) - degree of desired basis function
%
% Note: no error checking is performed for speed
%

% author: Nathan D. Cahill
% email: ndcahill@gmail.com
% date: 18 April 2008

switch deg
    case 0
        bi = (xi>=(-1/2))&(xi<(1/2));
    case 1
        bi = zeros(size(xi));
        k = (xi>=0) & (xi<1);
        bi(k) = 1 - xi(k);
        k(:) = (xi>=-1) & (xi<0);
        bi(k) = 1 + xi(k);
    case 2
        bi = zeros(size(xi));
        x2 = xi.*xi;
        k = (xi>=1/2) & (xi<3/2);
        bi(k) = (9/8) - (3/2).*xi(k) + (1/2).*x2(k);
        k = (xi>=-1/2) & (xi<1/2);
        bi(k) = (3/4) - x2(k);
        k = (xi>=-3/2) & (xi<-1/2);
        bi(k) = (9/8) + (3/2).*xi(k) + (1/2).*x2(k);
    case 3
        bi = zeros(size(xi));
        x2 = xi.*xi; x3 = x2.*xi;
        k = (xi>=1) & (xi<2);
        bi(k) = 4/3 - 2.*xi(k) + x2(k) - (1/6).*x3(k);
        k(:) = (xi>=0) & (xi<1);
        bi(k) = 2/3 - x2(k) + (1/2).*x3(k);
        k(:) = (xi>=-1) & (xi<0);
        bi(k) = 2/3 - x2(k) - (1/2).*x3(k);
        k(:) = (xi>=-2) & (xi<-1);
        bi(k) = 4/3 + 2.*xi(k) + x2(k) + (1/6).*x3(k);
    case 4
        bi = zeros(size(xi));
        x2 = xi.*xi; x3 = x2.*xi; x4 = x3.*xi;
        k = (xi>=3/2) & (xi<5/2);
        bi(k) = 625/384 - (125/48).*xi(k) + (25/16).*x2(k) - (5/12).*x3(k) + (1/24).*x4(k);
        k = (xi>=1/2) & (xi<3/2);
        bi(k) = 55/96 + (5/24).*xi(k) - (5/4).*x2(k) + (5/6).*x3(k) - (1/6).*x4(k);
        k = (xi>=-1/2) & (xi<1/2);
        bi(k) = 115/192 - (5/8).*x2(k) + (1/4).*x4(k);
        k = (xi>=-3/2) & (xi<-1/2);
        bi(k) = 55/96 - (5/24).*xi(k) - (5/4).*x2(k) - (5/6).*x3(k) - (1/6).*x4(k);
        k = (xi>=-5/2) & (xi<-3/2);
        bi(k) = 625/384 + (125/48).*xi(k) + (25/16).*x2(k) + (5/12).*x3(k) + (1/24).*x4(k);
    case 5
        bi = zeros(size(xi));
        x2 = xi.*xi; x3 = x2.*xi; x4 = x3.*xi; x5 = x4.*xi;
        k = (xi>=2) & (xi<3);
        bi(k) = 81/40 - (27/8).*xi(k) + (9/4).*x2(k) - (3/4).*x3(k) + (1/8).*x4(k) - (1/120).*x5(k);
        k = (xi>=1) & (xi<2);
        bi(k) = 17/40 + (5/8).*xi(k) - (7/4).*x2(k) + (5/4).*x3(k) - (3/8).*x4(k) + (1/24).*x5(k);
        k = (xi>=0) & (xi<1);
        bi(k) = 11/20 - (1/2).*x2(k) + (1/4).*x4(k) - (1/12).*x5(k);
        k = (xi>=-1) & (xi<0);
        bi(k) = 11/20 - (1/2).*x2(k) + (1/4).*x4(k) + (1/12).*x5(k);
        k = (xi>=-2) & (xi<-1);
        bi(k) = 17/40 - (5/8).*xi(k) - (7/4).*x2(k) - (5/4).*x3(k) - (3/8).*x4(k) - (1/24).*x5(k);
        k = (xi>=-3) & (xi<-2);
        bi(k) = 81/40 + (27/8).*xi(k) + (9/4).*x2(k) + (3/4).*x3(k) + (1/8).*x4(k) + (1/120).*x5(k);
    case 6
        bi = zeros(size(xi));
        x2 = xi.*xi; x3 = x2.*xi; x4 = x3.*xi; x5 = x4.*xi; x6 = x5.*xi;
        k = (xi>=5/2) & (xi<7/2);
        bi(k) = 117649/46080 - (16807/3840).*xi(k) + (2401/768).*x2(k) - (343/288).*x3(k) + (49/192).*x4(k) - (7/240).*x5(k) + (1/720).*x6(k);
        k = (xi>=3/2) & (xi<5/2);
        bi(k) = 1379/7680 + (1267/960).*xi(k) - (329/128).*x2(k) + (133/72).*x3(k) - (21/32).*x4(k) + (7/60).*x5(k) - (1/120).*x6(k);
        k = (xi>=1/2) & (xi<3/2);
        bi(k) = 7861/15360 - (7/768).*xi(k) - (91/256).*x2(k) - (35/288).*x3(k) + (21/64).*x4(k) - (7/48).*x5(k) + (1/48).*x6(k);
        k = (xi>=-1/2) & (xi<1/2);
        bi(k) = 5887/11520 - (77/192).*x2(k) + (7/48).*x4(k) - (1/36).*x6(k);
        k = (xi>=-3/2) & (xi<-1/2);
        bi(k) = 7861/15360 + (7/768).*xi(k) - (91/256).*x2(k) + (35/288).*x3(k) + (21/64).*x4(k) + (7/48).*x5(k) + (1/48).*x6(k);
        k = (xi>=-5/2) & (xi<-3/2);
        bi(k) = 1379/7680 - (1267/960).*xi(k) - (329/128).*x2(k) - (133/72).*x3(k) - (21/32).*x4(k) - (7/60).*x5(k) - (1/120).*x6(k);
        k = (xi>=-7/2) & (xi<-5/2);
        bi(k) = 117649/46080 + (16807/3840).*xi(k) + (2401/768).*x2(k) + (343/288).*x3(k) + (49/192).*x4(k) + (7/240).*x5(k) + (1/720).*x6(k);
    case 7
        bi = zeros(size(xi));
        x2 = xi.*xi; x3 = x2.*xi; x4 = x3.*xi; x5 = x4.*xi; x6 = x5.*xi; x7 = x6.*xi;
        k = (xi>=3) & (xi<4);
        bi(k) = 6405119470038039/1970324836974592 - (672537544353994073/118219490218475520).*xi(k) + (64/15).*x2(k) - (16/9).*x3(k) + (4/9).*x4(k) - (1/15).*x5(k) + (1/180).*x6(k) - (1/5040).*x7(k);
        k = (xi>=2) & (xi<3);
        bi(k) = -2173612320154509/9851624184872960 + (855120979246972939/354658470655426560).*xi(k) - (23/6).*x2(k) + (49/18).*x3(k) - (19/18).*x4(k) + (7/30).*x5(k) - (1/36).*x6(k) + (1/720).*x7(k);
        k = (xi>=1) & (xi<2);
        bi(k) = 103/210 - (7/90).*xi(k) - (1/10).*x2(k) - (7/18).*x3(k) + (1/2).*x4(k) - (7/30).*x5(k) + (1/20).*x6(k) - (1/240).*x7(k);
        k = (xi>=0) & (xi<1);
        bi(k) = 151/315 - (1/3).*x2(k) + (1/9).*x4(k) - (1/36).*x6(k) + (1/144).*x7(k);
        k = (xi>=-1) & (xi<0);
        bi(k) = 151/315 - (1/3).*x2(k) + (1/9).*x4(k) - (1/36).*x6(k) - (1/144).*x7(k);
        k = (xi>=-2) & (xi<-1);
        bi(k) = 103/210 + (7/90).*xi(k) - (1/10).*x2(k) + (7/18).*x3(k) + (1/2).*x4(k) + (7/30).*x5(k) + (1/20).*x6(k) + (1/240).*x7(k);
        k = (xi>=-3) & (xi<-2);
        bi(k) = -2173612320154509/9851624184872960 - (855120979246972939/354658470655426560).*xi(k) - (23/6).*x2(k) - (49/18).*x3(k) - (19/18).*x4(k) - (7/30).*x5(k) - (1/36).*x6(k) - (1/720).*x7(k);
        k = (xi>=-4) & (xi<-3);
        bi(k) = 6405119470038039/1970324836974592 + (672537544353994073/118219490218475520).*xi(k) + (64/15).*x2(k) + (16/9).*x3(k) + (4/9).*x4(k) + (1/15).*x5(k) + (1/180).*x6(k) + (1/5040).*x7(k);
    otherwise
        bi = ((xi + ((deg+1)/2)).*evalBSpline(xi+(1/2),deg-1) + ...
            (((deg+1)/2) - xi).*evalBSpline(xi-(1/2),deg-1))./deg;
end

author: Sable

repo: mcbench-benchmarks


text:
function varargout = ft_error(varargin)

% FT_ERROR prints an error message on screen, just like the standard ERROR function.
%
% Use as
%   ft_error(...)
% with arguments similar to fprintf, or
%   ft_error(msgId, ...)
% with arguments similar to error.
%
% See also FT_ERROR, FT_WARNING, FT_NOTICE, FT_INFO, FT_DEBUG, ERROR, WARNING

% Copyright (C) 2017, Robert Oostenveld
%
% 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$

if nargout
  varargout{:} = ft_notification(varargin{:});
elseif isequal(varargin, {'last'})
  % return an answer anyway
  varargout{1} = ft_notification(varargin{:});
else
  ft_notification(varargin{:});
end

author: spm


text:
function [PQ, PV, REF, NONE, BUS_I, BUS_TYPE, PD, QD, GS, BS, BUS_AREA, VM, ...
    VA, BASE_KV, ZONE, VMAX, VMIN, LAM_P, LAM_Q, MU_VMAX, MU_VMIN] = idx_bus
%IDX_BUS   Defines constants for named column indices to bus matrix.
%   Example:
%
%   [PQ, PV, REF, NONE, BUS_I, BUS_TYPE, PD, QD, GS, BS, BUS_AREA, VM, ...
%   VA, BASE_KV, ZONE, VMAX, VMIN, LAM_P, LAM_Q, MU_VMAX, MU_VMIN] = idx_bus;
%
%   Some examples of usage, after defining the constants using the line above,
%   are:
%
%    Pd = bus(4, PD);       % get the real power demand at bus 4
%    bus(:, VMIN) = 0.95;   % set the min voltage magnitude to 0.95 at all buses
%
%   The index, name and meaning of each column of the bus matrix is given
%   below:
%
%   columns 1-13 must be included in input matrix (in case file)
%    1  BUS_I       bus number (positive integer)
%    2  BUS_TYPE    bus type (1 = PQ, 2 = PV, 3 = ref, 4 = isolated)
%    3  PD          Pd, real power demand (MW)
%    4  QD          Qd, reactive power demand (MVAr)
%    5  GS          Gs, shunt conductance (MW demanded at V = 1.0 p.u.)
%    6  BS          Bs, shunt susceptance (MVAr injected at V = 1.0 p.u.)
%    7  BUS_AREA    area number, (positive integer)
%    8  VM          Vm, voltage magnitude (p.u.)
%    9  VA          Va, voltage angle (degrees)
%    10 BASE_KV     baseKV, base voltage (kV)
%    11 ZONE        zone, loss zone (positive integer)
%    12 VMAX        maxVm, maximum voltage magnitude (p.u.)
%    13 VMIN        minVm, minimum voltage magnitude (p.u.)
%
%   columns 14-17 are added to matrix after OPF solution
%   they are typically not present in the input matrix
%                   (assume OPF objective function has units, u)
%    14 LAM_P       Lagrange multiplier on real power mismatch (u/MW)
%    15 LAM_Q       Lagrange multiplier on reactive power mismatch (u/MVAr)
%    16 MU_VMAX     Kuhn-Tucker multiplier on upper voltage limit (u/p.u.)
%    17 MU_VMIN     Kuhn-Tucker multiplier on lower voltage limit (u/p.u.)
%
%   additional constants, used to assign/compare values in the BUS_TYPE column
%    1  PQ    PQ bus
%    2  PV    PV bus
%    3  REF   reference bus
%    4  NONE  isolated bus
%
%   See also DEFINE_CONSTANTS.

%   MATPOWER
%   Copyright (c) 1996-2016, Power Systems Engineering Research Center (PSERC)
%   by Ray Zimmerman, PSERC Cornell
%
%   This file is part of MATPOWER.
%   Covered by the 3-clause BSD License (see LICENSE file for details).
%   See https://matpower.org for more info.

%% define bus types
PQ      = 1;
PV      = 2;
REF     = 3;
NONE    = 4;

%% define the indices
BUS_I       = 1;    %% bus number (1 to 29997)
BUS_TYPE    = 2;    %% bus type (1 - PQ bus, 2 - PV bus, 3 - reference bus, 4 - isolated bus)
PD          = 3;    %% Pd, real power demand (MW)
QD          = 4;    %% Qd, reactive power demand (MVAr)
GS          = 5;    %% Gs, shunt conductance (MW at V = 1.0 p.u.)
BS          = 6;    %% Bs, shunt susceptance (MVAr at V = 1.0 p.u.)
BUS_AREA    = 7;    %% area number, 1-100
VM          = 8;    %% Vm, voltage magnitude (p.u.)
VA          = 9;    %% Va, voltage angle (degrees)
BASE_KV     = 10;   %% baseKV, base voltage (kV)
ZONE        = 11;   %% zone, loss zone (1-999)
VMAX        = 12;   %% maxVm, maximum voltage magnitude (p.u.)      (not in PTI format)
VMIN        = 13;   %% minVm, minimum voltage magnitude (p.u.)      (not in PTI format)

%% included in opf solution, not necessarily in input
%% assume objective function has units, u
LAM_P       = 14;   %% Lagrange multiplier on real power mismatch (u/MW)
LAM_Q       = 15;   %% Lagrange multiplier on reactive power mismatch (u/MVAr)
MU_VMAX     = 16;   %% Kuhn-Tucker multiplier on upper voltage limit (u/p.u.)
MU_VMIN     = 17;   %% Kuhn-Tucker multiplier on lower voltage limit (u/p.u.)

author: MATPOWER

repo: matpower

save_path: github-repos/MATLAB/MATPOWER-matpower


text:
function calc_prob

%%

nGames = 50; % number of games
bet = 100; % units: escudos

pwin = 0.8; % probability of winning
ploss = 1 - pwin; % loosing probability

rate = 0.25;  % paying bet (what the house pays). This will be the average of
% bet payment with multiple bets

earnings =  ( (rate * win) -  ploss) * bet * nGames

repo: BeatTheBookie


text:
function [latb,lonb] = bufferm2(varargin) %lat,lon,dist,direction,npts,outputformat)
%BUFFERM2 Computes buffer zone around a polygon
%
% [latb,lonb] = bufferm2(lat,lon,dist,direction)
% [latb,lonb] = bufferm2(lat,lon,dist,direction,npts)
% [latb,lonb] = bufferm2(lat,lon,dist,direction,npts,outputformat)
% [xb,  yb]   = bufferm2('xy',x,y,dist,direction,npts,outputformat)
%
% This function was originally designed as a replacement for the Mapping
% Toolbox function bufferm, which calculates a buffer zone around a
% polygon. The original bufferm function had some serious bugs that could
% result in incorrect buffer results and/or errors, and was also very slow.
% As of R2006b, those bugs have been fixed.  However, this version still
% maintains a few advantages over the original:
%
%   - Can be applied to polygons in either geographical space (as in
%   bufferm) or in cartesian coordinates.
%
%   - Better treatment of polygon holes.  The original function simply
%   filled in all holes; this version trims or pads holes according to the
%   buffer width given.
%
% Input and output format is identical to bufferm unless the 'xy' option is
% specified, so it can be used interchangeably.
%
% Input variables:
%
%   lat:            Latitude values defining the polygon to be buffered.
%                   This can be either a NaN-delimited vector, or a cell
%                   array containing individual polygonal contours (each of
%                   which is a vector). External contours should be listed
%                   in a clockwise direction, and internal contours (holes)
%                   in a counterclockwise direction.
%
%   lon:            Longitude values defining the polygon to be buffered.
%                   Same format as lat.
%
%   dist:           Width of buffer, in degrees of arc along the surface
%                   (unless 'xy' is used, in which case units correspond to
%                   x-y coordinates)
%
%   direction:      'in' or 'out'
%
%   npts:           Number of points used to contruct the circles around
%                   each polygon vertex.  If omitted, default is 13.
%
%   outputformat:   'vector' (NaN-delimited vectors), 'cutvector'
%                   (NaN-clipped vectors with cuts connecting holes to the
%                   exterior of the polygon), or 'cell' (cell arrays in
%                   which each element of the cell array is a separate
%                   polygon), defining format of output.  If omitted,
%                   default is 'vector'.
%
%   'xy':           If first input is 'xy', then data will be assumed to
%                   lie on a cartesian plane rather than on a sphere.  Use
%                   x and y coordinates as first two inputs rather than lat
%                   and lon.  Units of x, y, and distance should be the
%                   same.
%
% Output variables:
%
%   latb:           Latitude values for buffer polygon
%
%   lonb:           Longitude values for buffer polygon
%
% Example:
%
%   load conus
%   tol = 0.1; % Tolerance for simplifying polygon outlines
%   [reducedlat, reducedlon] = reducem(gtlakelat, gtlakelon, tol);
%   dist = 1;  % Buffer distance in degrees
%   [latb, lonb] = bufferm2(reducedlat, reducedlon, dist, 'out');
%   figure('Renderer','painters')
%   usamap({'MN','NY'})
%   geoshow(latb, lonb, 'DisplayType', 'polygon', 'FaceColor', 'yellow')
%   geoshow(gtlakelat, gtlakelon,...
%                       'DisplayType', 'polygon', 'FaceColor', 'blue')
%   geoshow(uslat, uslon)
%   geoshow(statelat, statelon)
%
% See also:
%
%   bufferm, polybool

% Copyright 2010 Kelly Kearney

%---------------------------
% Check input
%---------------------------

error(nargchk(3,7,nargin));

% Determine if geographic or cartesian

if ischar(varargin{1}) && strcmp(varargin{1}, 'xy')
    geo = false;
    param = varargin(2:end);
else
    geo = true;
    param = varargin;
end

% Set defaults if not provided as input

nparam = length(param);

if geo
    [lat, lon, dist] = deal(param{1:3});
else
    [lon, lat, dist] = deal(param{1:3}); % lon = x, lat = y for mental clarity, will switch back at end
end

if nparam < 4
    direction = 'out';
else
    direction = param{4};
end

if nparam < 5
    npts = 13;
else
    npts = param{5};
end

if nparam < 6
    outputformat = 'vector';
else
    outputformat = param{6};
end

% Check format and dimensions of input

if ~ismember(direction, {'in', 'out'})
    error('Direction must be either ''in'' or ''out''.');
end

if ~ismember(outputformat, {'vector', 'cutvector', 'cell'})
    error('Unrecognized output format flag.');
end

if ~isnumeric(dist) || numel(dist) > 1
    error('Distance must be a scalar.')
end

if ~isnumeric(npts) || numel(npts) > 1
    error('Number of points must be a scalar.')
end

if iscell(lat)
    for il = 1:numel(lat)
        if ~isvector(lat{il}) | ~isvector(lon{il}) | ~isequal(length(lat{il}), length(lon{il}))
            error('Lat (or x) and lon (or y) must be vectors or cells of vectors with identical dimensions');
        end
        lat{il} = lat{il}(:);
        lon{il} = lon{il}(:);
    end
else
    if ~isvector(lat) || ~isvector(lon) || ~isequal(length(lat), length(lon))
        error('Lat (or x) and lon (or y) must be vectors or cells of vectors with identical dimensions');
    end
    lat = lat(:);
    lon = lon(:);
end

%---------------------------
% Split polygon(s) into
% separate faces
%---------------------------

if iscell(lat)
    [lat, lon] = polyjoin(lat, lon);  % In case multiple faces in one cell.
end

[latcells, loncells] = polysplit(lat, lon);

%---------------------------
% Create buffer shapes
%---------------------------

plotflag = 0;

if plotflag

    Plt.x = lon;
    Plt.y = lat;

end

latcrall = cell(0);
loncrall = cell(0);

for ipoly = 1:length(latcells)

    % Circles around each vertex

    if geo
        [latc, lonc] = calccircgeo(latcells{ipoly}, loncells{ipoly}, dist, npts);
    else
        [lonc, latc] = calccirccart(loncells{ipoly}, latcells{ipoly}, dist, npts);
    end

    % Rectangles around each edge

    if geo
        [latr, lonr] = calcrecgeo(latcells{ipoly}, loncells{ipoly}, dist);
    else
        [lonr, latr] = calcreccart(loncells{ipoly}, latcells{ipoly}, dist);
    end

    % Union of circles and rectangles

    if plotflag
        Plt.rectx = lonr;
        Plt.recty = latr;
        Plt.circx = lonc;
        Plt.circy = latc;
    end

    [latc, lonc] = multipolyunion(latc, lonc);
    [latr, lonr] = multipolyunion(latr, lonr);

    if plotflag
        Plt.rectcombox = lonr;
        Plt.rectcomboy = latr;
        Plt.circcombox = lonc;
        Plt.circcomboy = latc;
    end

    [loncr, latcr] = polybool('+', lonr, latr, lonc, latc);

    % Union of new circle/rectangle combo with that from other faces

    [loncrall, latcrall] = polybool('+', loncrall, latcrall, loncr, latcr);

    % Plotting (for debugging only)

    if plotflag

        Plt.allx = loncrall;
        Plt.ally = latcrall;

        if ipoly == 1
            figure;
            plot(Plt.x, Plt.y, 'k', 'linewidth', 2);
            hold on
        end

        plot(cat(2, Plt.rectx{:}), cat(2, Plt.recty{:}), 'b');
        plot(cat(2, Plt.circx{:}), cat(2, Plt.circy{:}), 'r');
        plot(Plt.allx{1}, Plt.ally{1}, 'g', 'linewidth', 2);

    end

end

%---------------------------
% Calculate union/difference
%---------------------------

switch direction
    case 'out'
        [lonb, latb] = polybool('+', loncells, latcells, loncrall, latcrall);
    case 'in'
        [lonb, latb] = polybool('-', loncells, latcells, loncrall, latcrall);
end

if plotflag
    [Plt.yfinal, Plt.xfinal] = polyjoin(latb, lonb);
    plot(Plt.xfinal, Plt.yfinal, 'linestyle', '--', 'color', [0 .5 0], 'linewidth', 2);
end

%---------------------------
% Reformat output
%---------------------------

if ~geo
    y = latb; % Switch, since cartesion uses opposite order
    x = lonb;
    latb = x;
    lonb = y;
end

switch outputformat
    case 'vector'
        [latb, lonb] = polyjoin(latb, lonb);
    case 'cutvector'
        [latb, lonb] = polycut(latb, lonb);
    case 'cell'
end

%**************************************************************************

function [latc, lonc] = calccircgeo(lat, lon, radius, npts)
% lat and lon: n x 1 vectors
% radius: scalar

radius = ones(length(lat),1) * radius;
[latc, lonc] = scircle1(lat, lon, radius, [], [], [], npts);
latc = num2cell(latc, 1);
lonc = num2cell(lonc, 1);

function [latr, lonr] = calcrecgeo(lat, lon, halfwidth)
% lat and lon: n x 1 vectors
% halfwidth: scalar

range = halfwidth * ones(length(lat)-1, 1);

az = azimuth(lat(1:end-1), lon(1:end-1), lat(2:end), lon(2:end));

[latbl1,lonbl1] = reckon(lat(1:end-1), lon(1:end-1), range, az-90);
[latbr1,lonbr1] = reckon(lat(1:end-1), lon(1:end-1), range, az+90);
[latbl2,lonbl2] = reckon(lat(2:end),   lon(2:end),   range, az-90);
[latbr2,lonbr2] = reckon(lat(2:end),   lon(2:end),   range, az+90);

latr = [latbl1 latbl2 latbr2 latbr1 latbl1]';
lonr = [lonbl1 lonbl2 lonbr2 lonbr1 lonbl1]';
latr = num2cell(latr, 1);
lonr = num2cell(lonr, 1);

function [latu, lonu] = multipolyunion(lat, lon)
% lat and lon are n x 1 cell arrays of vectors

latu = lat{1};
lonu = lon{1};

for ip = 2:length(lat)
    [lonu, latu] = polybool('+', lonu, latu, lon{ip}, lat{ip});
end
[latu, lonu] = polysplit(latu, lonu);

function [xc, yc] = calccirccart(x, y, radius, npts)

ang = linspace(0, 2*pi, npts+1);
ang = ang(end-1:-1:1);
xc = bsxfun(@plus, x, radius * cos(ang));
yc = bsxfun(@plus, y, radius * sin(ang));
xc = num2cell(xc', 1);
yc = num2cell(yc', 1);

% if ~ispolycw(x,y)
%     [xc,yc] = poly2ccw(xc,yc);
% end

function [xrec, yrec] = calcreccart(x, y, halfwidth)

dx = diff(x);
dy = diff(y);

is1 = dx >= 0 & dy >= 0;
is2 = dx < 0 & dy >= 0;
is3 = dx < 0 & dy < 0;
is4 = dx >= 0 & dy < 0;

ish1 = dy == 0 & dx > 0;
ish2 = dy == 0 & dx < 0;

theta = zeros(5,1);
theta(is1 | is3) = atan(dy(is1 | is3)./dx(is1 | is3));
theta(is2 | is4) = -atan(dy(is2 | is4)./dx(is2 | is4));

[xl,xr,yl,yr] = deal(zeros(size(dx)));

xl(is1) = -halfwidth * sin(theta(is1));
xr(is1) =  halfwidth * sin(theta(is1));
yl(is1) =  halfwidth * cos(theta(is1));
yr(is1) = -halfwidth * cos(theta(is1));

xl(is2) = -halfwidth * sin(theta(is2));
xr(is2) =  halfwidth * sin(theta(is2));
yl(is2) = -halfwidth * cos(theta(is2));
yr(is2) =  halfwidth * cos(theta(is2));

xl(is3) =  halfwidth * sin(theta(is3));
xr(is3) = -halfwidth * sin(theta(is3));
yl(is3) = -halfwidth * cos(theta(is3));
yr(is3) =  halfwidth * cos(theta(is3));

xl(is4) =  halfwidth * sin(theta(is4));
xr(is4) = -halfwidth * sin(theta(is4));
yl(is4) =  halfwidth * cos(theta(is4));
yr(is4) = -halfwidth * cos(theta(is4));

xrec = [xl+x(1:end-1) xl+x(2:end) xr+x(2:end) xr+x(1:end-1) xl+x(1:end-1)];
yrec = [yl+y(1:end-1) yl+y(2:end) yr+y(2:end) yr+y(1:end-1) yl+y(1:end-1)];

xrec = num2cell(xrec, 2);
yrec = num2cell(yrec, 2);


text:
function f = eval(SO3F,g,varargin)
% evaluate sum of unimodal components at orientation g
%
% Syntax
%   f = SO3F.eval(g)
%
% Input
%  SO3F - @SO3FunUnimodal
%  rot  - @rotation
%
% Options
%  exact   -
%  epsilon -
%
% Description
% general formula:
%
% $$ f(r) = sum_j w_j \psi(r,c_j) $$

% if isa(g,'orientation')
%   ensureCompatibleSymmetries(SO3F,g)
% end

% decide along which dimension to split the summation matrix
if isa(g,'SO3Grid')
  lg1 = length(g);
else
  lg1 = -length(g);
end

if isa(SO3F.center,'SO3Grid')
  lg2 = length(SO3F.center);
else
  lg2 = -length(SO3F.center);
end

along = (lg1 > lg2 && lg1 > 0) || (abs(lg1) > abs(lg2) && lg2 < 0);
if along
  num = abs(lg2);
else
  num = abs(lg1);
end

% init variables
f = SO3F.c0 * ones(size(g));
iter = 0; numiter = 1; ind = 1; %for first run

% now iterate along the splitting
while iter <= numiter
  if iter > 0% split
    ind = 1 + (1+(iter-1)*diter:min(num-1,iter*diter));
    if isempty(ind), return; end
  end

  %eval the kernel
  if along
    M = SO3F.psi.K_symmetrised(g,SO3F.center(ind),SO3F.CS,SO3F.SS,'nocubictrifoldaxis',varargin{:});
    if SO3F.antipodal
      M = 0.5*(M + SO3F.psi.K_symmetrised(g,inv(SO3F.center(ind)),SO3F.CS,SO3F.SS,'nocubictrifoldaxis',varargin{:}));
    end
    f = f + reshape(full(M * reshape(SO3F.weights(ind),[],1)),size(f));
  else
    M = SO3F.psi.K_symmetrised(g(ind),SO3F.center,SO3F.CS,SO3F.SS,'nocubictrifoldaxis',varargin{:});
    if SO3F.antipodal
      M = 0.5*(M + SO3F.psi.K_symmetrised(inv(g(ind)),SO3F.center,SO3F.CS,SO3F.SS,'nocubictrifoldaxis',varargin{:}));
    end
    f(ind) = f(ind) + reshape(full(M * SO3F.weights(:)),size(f(ind)));
  end

  if num == 1
    return
  elseif iter == 0 % iterate due to memory restrictions?
    numiter = ceil( max(1,nnz(M))*num / getMTEXpref('memory',512 * 1024) / 256 );
    diter = ceil(num / numiter);
  end

  if numiter > 1 && ~check_option(varargin,'silent'), progress(iter,numiter); end

  iter = iter + 1;
end

if isalmostreal(f)
  f = real(f);
end

end

author: mtex-toolbox

repo: mtex

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


text:
function K = rbf_of_dist(kern,dat1,dat2,ind1,ind2,kerParam),
  %
  % introduce nonlinearity by performing POLYNOMIAL on input kernel matrix
  %
  K=get_x(dat1,ind1);
  K = exp(-K/(2*kerParam^2));

author: canlab

repo: CanlabCore

save_path: github-repos/MATLAB/canlab-CanlabCore


text:
function [ A, D ] = SIC( Fopt, H, N , SNR)
[Nt,Ns] = size(Fopt);
M = Nt/N;
T = eye(N);
A = [];
D = [];
P = [];
R = [eye(M),zeros(M,M*(N-1))];
G = R*H'*H*R';
for n = 1:N
    [U,S,V] = svd(G);
    v1 = V(:,1);

    A = blkdiag( A, exp(1i*angle(v1))/sqrt(M) );
    D = blkdiag( D, norm(v1,1)/sqrt(M) );

%     P = [P, [zeros(M*(n-1),1);norm(v1,1)/M*exp(1i*angle(v1));zeros(M*(N-n),1)] ];
%     T = eye(size(H,1)) + SNR/Ns * H*P*P'*H';
%     R = [zeros(M,M*(n-1)),eye(M),zeros(M,M*(N-n))];
%     G = R*H'*inv(T)*H*R';

    G = G - SNR/Ns * S(1)^2 * v1*v1' / (1+SNR/Ns*S(1));
end

author: yuxianghao


text:
function varargout = ellipsoid( a, b, c )
%ELLIPSOID  Generate an ellipsoid-like surface. (Not necessarily an ellipsoid!)
%   ELLIPSOID(A,B,C), where A, B, and C are CHEBFUN2 objects on the domain [0
%   pi]x[0 2*pi] plots the "ELLIPSOID" of semi axis lengths A(th,phi),
%   B(th,phi), and C(th,phi).
%
%   [X Y Z]=ELLIPSOID(A,B,C) returns X, Y, and Z as CHEBFUN2 objects such that
%   SURF(X,Y,Z) plots an ELLIPSOID of semi axis lengths A(th,phi), B(th,phi),
%   and C(th,phi).
%
%   F = ELLIPSOID(A,B,C) returns the CHEBFUN2V representing the ELLIPSOID
%   SURF(F) plots the ELLIPSOID.
%
%   Omitting output arguments causes the ELLIPSOID command to be displayed with
%   a SURF command and no outputs are returned.
%
% For the ellipsoid:
%   a = chebfun2(@(th,phi) 1+0*th,[0 pi 0 2*pi]);
%   F = ellipsoid(a,2*a,3*a); surf(F)
%
% For a badly shaped ship:
%   a = chebfun2(@(th,phi) th,[0 pi 0 2*pi])
%   F = ellipsoid(a,2*a,cos(2*a)+.25); surf(F)
%
% See also SPHERE, CYLINDER.

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

if ( nargin < 3 )
    error('CHEBFUN:CHEBFUN2:ellipsoid:inputs', 'Incorrect input arguments.');
end

% ELLIPSOID with axis lengths a(th,phi), b(th,phi), c(th,phi).
dom = [0 pi 0 2*pi];
th = chebfun2(@(th, phi) th, dom);
phi = chebfun2(@(th, phi) phi, dom);

x = a.*sin(th).*cos(phi);
y = b.*sin(th).*sin(phi);
z = c.*cos(th);

if ( nargout == 0 )
    surf(x, y, z)
    axis equal
elseif ( nargout == 1 )
    varargout = { [x ; y ; z] };
else
    varargout = { x, y, z };
end

end

author: chebfun

repo: chebfun

save_path: github-repos/MATLAB/chebfun-chebfun


text:
function vl_test_binsum
% VL_TEST_BINSUM  Test VL_BINSUM function

testh({[0 0],   1, 2},                  [0 1]   ) ;
testh({[1 7],  -1, 1},                  [0 7]   ) ;
testh({[1 7],  -1, [1 2 2 2 2 2 2 2]},  [0 0]   ) ;
testh({eye(3), [1 1 1],  [1 2 3],  1 }, 2*eye(3)) ;
testh({eye(3), [1 1 1]', [1 2 3]', 2 }, 2*eye(3)) ;
testh({eye(3), 1, [1 2 3],  1 },        2*eye(3)) ;
testh({eye(3), 1, [1 2 3]', 2 },        2*eye(3)) ;

Z = zeros(3,3,3) ;
B = 3*ones(3,1,3) ;
R = Z ; R(:,3,:) = 17 ;

testh({Z, 17, B, 2}, R) ;

Z = zeros(3,3,3) ;
B = 3*ones(3,3,1) ;
X = zeros(3,3,1) ; X(:,:,1) = 17 ;
R = Z ; R(:,:,3) = 17 ;

testh({Z, X, B, 3}, R) ;

function testh(args, H_)
H__ = vl_binsum(args{:}) ;
if any(any(any(H_ ~= H__)))
  fprintf('H:\n') ; disp(args{1});
  fprintf('X:\n') ; disp(args{2});
  fprintf('B:\n') ; disp(args{3});
  if length(args) > 3,
    fprintf('d:\n') ; disp(args{4}) ;
  end
  fprintf('R computed:\n') ; disp(H__) ;
  fprintf('R correct:\n') ; disp(H_) ;
  error('vl_binsum regression test failed') ;
end

author: jianxiongxiao

repo: ProfXkit


text:
function stats = applyStatsfun(problem, x, storedb, options, stats)
% Apply the statsfun function to a stats structure (for solvers).
%
% function stats = applyStatsfun(problem, x, storedb, options, stats)
%
% Applies the options.statsfun user supplied function to the stats
% structure, if it was provided, with the appropriate inputs, and returns
% the (possibly) modified stats structure.
%
% See also:

% This file is part of Manopt: www.manopt.org.
% Original author: Nicolas Boumal, April 3, 2013.
% Contributors:
% Change log:

	if isfield(options, 'statsfun')

		is_octave = exist('OCTAVE_VERSION', 'builtin');
		if ~is_octave
			narg = nargin(options.statsfun);
		else
			narg = 4;
		end

        switch narg
            case 3
                stats = options.statsfun(problem, x, stats);
            case 4
                store = getStore(problem, x, storedb);
                stats = options.statsfun(problem, x, stats, store);
            otherwise
                warning('manopt:statsfun', ...
                        'statsfun unused: wrong number of inputs');
        end
	end

end

author: andrewssobral

repo: lrslibrary


text:
function cvt_movie5 ( )

%*****************************************************************************80
%
%% CVT_MOVIE5 creates a movie from data files.
%
%  Discussion:
%
%    This simple example creates a movie by reading in a set of point
%    values at each step in an iteration.
%
%    The point values are 389 points inside the P08 "holey pie" region.
%    They were computed in a CVT calculation external to this procedure.
%
%    This routine simply reads in the data, adds a picture of the boundary
%    of the region, and computes the Voronoi diagram, before "taking a picture"
%    which creates a JPEG image.  The JPEG images can easily be made into a movie.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    03 July 2005
%
%  Author:
%
%    John Burkardt
%
%  Local parameters:
%
%    Local parameter FILE_NAME, the name of the data file being read.
%
  data_file = 'p08_hbf_000.txt';

  timestamp ( );

  fprintf ( 1, '\n' );
  fprintf ( 1, 'CVT_MOVIE5:\n' );
  fprintf ( 1, '  Create a movie of a CVT computation.\n' );
  fprintf ( 1, '  This example uses the "holey pie" geometry.\n' );
  fprintf ( 1, '\n' );
  fprintf ( 1, '  The data has been precomputed.\n' );
  fprintf ( 1, '  The first data file to read is %s\n', data_file );
%
%  This switch controls whether the Delaunay triangulation is displayed.
%
  delaunay_display = 0;

  if ( delaunay_display )
    fprintf ( 1, '\n' );
    fprintf ( 1, 'The Delaunay triangulation will be displayed.\n' );
  else
    fprintf ( 1, '\n' );
    fprintf ( 1, 'The Delaunay triangulation will NOT be displayed.\n' );
  end

  frame_num = 0;

  frame_file = 'cvt_movie5_000.jpg';
  frame_title = 'cvt\_movie5\_000.jpg';

  while ( 1 )

    if ( ~file_exist ( data_file ) )
      break
    end

    p = load ( data_file );
    frame_num = frame_num + 1;
%
%  Compute the Delaunay triangle information T for the current nodes.
%
    t = delaunay ( p(:,1), p(:,2) );
%
%  Display a graph of the generators.
%
%  This was a little tricky.  The TRIMESH command must come before
%  the VORONOI command in order for both items to appear.
%
%  There are obvious discrepancies between the Delaunay and Voronoi
%  graphs which suggest that there are mistakes in the Voronoi diagram
%  as computed by MATLAB's VORONOI command!
%
    if ( delaunay_display )
      trimesh ( t, p(:,1), p(:,2), zeros(n,1) )
      hold on
    end
%
%  Here, we use a new, unreleased version of the VORONOI command, which
%  does a better job of handling the semi-infinite sides.
%
    voronoi ( p(:,1), p(:,2), t );

    title ( frame_title );
    axis equal
    axis ( [ 0.0, 1.0, -0.29, +0.29 ] )
    segment_num = p08_boundary_segment_num ( 'DUMMY' );
    h = 0.05;
    for segment_index = 1 : segment_num
      segment_length = p08_boundary_segment_length ( segment_index, h );
      segment = p08_boundary_segment ( segment_index, 2, segment_length );
      line ( segment(1,:), segment(2,:), 'Color', 'r' );
    end
    view ( 2 )
    drawnow

    if ( delaunay_display )
      hold off
    end

    F = getframe;
    [ X, map ] = frame2im ( F );
    imwrite ( X, frame_file, 'JPEG' );
%
%  Prepare for next step.
%
    data_file = file_name_inc ( data_file );
    frame_file = file_name_inc ( frame_file );
    frame_title = file_name_inc ( frame_title );

  end

  fprintf ( 1, '\n' );
  fprintf ( 1, '  Number of frames created is %d\n', frame_num );

  fprintf ( 1, '\n' );
  fprintf ( 1, 'CVT_MOVIE5:\n' );
  fprintf ( 1, '  Normal end of execution.\n' );

  fprintf ( 1, '\n' );
  timestamp ( );

  return
end

author: johannesgerer

repo: jburkardt-m


text:
function A = create_2d_image_upsample_matrix(insz, outsz, varargin)
    % (insz, outsz, {interp_type})
    interp_type = 1;
    if nargin >= 3
        interp_type = varargin{1};
    end
    [m1, m2] = ndgrid(1:outsz(1), 1:outsz(2));
    kk = (insz - 1) ./ (outsz - 1);
    dd = 1 - kk;
    t1 = kk(1) * m1 + dd(1);
    t2 = kk(2) * m2 + dd(2);
    T = cat(3, t1 - m1, t2 - m2);
    if interp_type == 0
        A = create_2d_bilinear_interp_matrix(T, insz);
    elseif interp_type == 1
        A = create_2d_bicubic_interp_matrix(T, insz);
    end
end

repo: pTVreg


text:
function mrVistaData = dtiXformGetMrVistaDataForFibers(handles, view, scanNum, fgNum)
%
% mrVistaData = dtiXformGetMrVistaDataForFibers(handles, [view], [scanNum], [fgNum])
%
% Uses the xformVAnatToAcpc (see dtiXformVanatCompute) to convert the specified
% fiber group to mrVista vAnatomy coords and pulls data from the selected
% gray view. Note that only fiber end-points are used, and they are
% projected to the nearest gray-matter coord.
%
% HISTORY:
% 2005.07.28 RFD (bob@white.stanford.edu) wrote it.

if(~exist('view','var') | isempty(view))
    view = getSelectedGray;
end
if(~exist('scanNum','var') | isempty(scanNum))
    scanNum = getCurScan(view);
end
if(~exist('fgNum','var') | isempty(fgNum))
    fgNum = handles.curFiberGroup;
end

% This is in vAnatomy voxels. We should probably scale by mmPerVox.
distThresh = 3;

fg = handles.fiberGroups(fgNum);

% Select only layer 1 for the grayCoords.
grayCoords = view.coords(:,view.nodes(6,:)==1);

distSqThresh = distThresh.^2;
fibNum = 0;
allCoords = [];
requireBothEnds = 1;

% We are only grabbing data for fiber endpoints, so we need num fibers * 2.
% We will apply a 4x4 xform, so we need homogeneous coords (nx4).
coords = zeros(length(fg.fibers)*2, 3);
for(ii=1:length(fg.fibers))
    % Get coords for fiber endpoints (first and last point)
    coords((ii-1)*2+1,:) = fg.fibers{ii}(:,1)';
    coords((ii-1)*2+2,:) = fg.fibers{ii}(:,end)';
end
coords = mrAnatXformCoords(inv(handles.xformVAnatToAcpc), coords);

% Now project endpoint coords to nearest gray matter coord
[indices, bestSqDist] = nearpoints(coords', grayCoords);
tooFarInds = bestSqDist>distSqThresh;
if(requireBothEnds)
    % Sometimes we want to delete fibers that have only one endpoint close
    % enough to the gray matter (eg. if we want to test connectivity).
    notBoth = tooFarInds(1:2:end)|tooFarInds(2:2:end);
    tooFarInds = repmat(notBoth,2,1);
    tooFarInds = tooFarInds(:)';
end
indices(tooFarInds) = [];

%mrVistaData.coords = grayCoords(:,indices);

mrVistaData.co = getCurDataROI(view,'co',scanNum,grayCoords(:,indices));
noData = isnan(mrVistaData.co);
mrVistaData.ph = getCurDataROI(view,'ph',scanNum,grayCoords(:,indices));
mrVistaData.amp = getCurDataROI(view,'amp',scanNum,grayCoords(:,indices));
mrVistaData.co(noData) = [];
mrVistaData.ph(noData) = [];
mrVistaData.amp(noData) = [];
if(~isempty(view.map) && ~isempty(view.map{scanNum}))
    mrVistaData.map = view.map{curScan}(indices);
    mrVistaData.map(noData) = [];
else
    mrVistaData.map = [];
end
mrVistaData.grayInds = indices(~noData);

fibNum = length(indices);

disp([num2str(fibNum) ' fibers had endpoints within ' num2str(sqrt(distSqThresh)) ' units of gray coords.']);
return;

author: vistalab

repo: vistasoft

save_path: github-repos/MATLAB/vistalab-vistasoft


text:
function [essentialRxn4Models, dataStruct] = essentialRxn4MultipleModels(allModels, objFun)
% Performs single reactions deletions across multiple models and identifies which reactions
% have variable essentiality across all models for the chosen objective function
% (eg. it identifies reactions whos deletion would result in a zero flux through ATP
% consumption reaction - ATPM).
%
% USAGE:
%
%    [essentialRxn4Models, dataStruct] = essentialRxn4MultipleModels(allModels, objFun)
%
% INPUTS:
%    allModels:    directory of the structure with multiple COBRA model structures
%    objFun:       string with objective function reaction (e.g. 'ATPM')
%
% OUTPUT:
%    essentialRxn4Models:    Table with reaction fluxes (within the objective function reaction)
%                            after single deletion of model reaction (rows) across models (columns)
%    dataStruct:             Structure with all models
%
% EXAMPLE:
%
%    [essentialRxn4Models, dataStruct] = essentialRxn4MultipleModels(allModels, 'ATPM')
%
% .. Authors:
%	    - Dr. Miguel A.P. Oliveira, 13/11/2017, Luxembourg Centre for Systems Biomedicine, University of Luxembourg
%   	- Diana C. El Assal-Jordan, 13/11/2017, Luxembourg Centre for Systems Biomedicine, University of Luxembourg

if isstr(allModels)  % Locate COBRA models in a directory and load them into a structure:
    modelsDir = allModels;
    clear allModels
    allModelFilenames = dir(strcat(modelsDir, '/', '*.mat'));
    for i = 1:size(allModelFilenames, 1)
        % Extract model name from filenames
        match = {'.mat', '.'};
        str = allModelFilenames(i).name;
        for j = 1:size(match, 2)
            str = strrep(str, match{j}, {''});
        end
        newFilename{1, i} = horzcat(str{1}, '_model');

        loadedFile = load(strcat(modelsDir, '/', allModelFilenames(i).name));
        fields = fieldnames(loadedFile);
        model = loadedFile.(fields{1, 1});
        allModels.(newFilename{1, i}) = model;
    end
end

% Load structure with models and perform singleRxnDeletion in all COBRA models
modelNames = fieldnames(allModels);
numModels = size(modelNames, 1);
sumRxnSubsystems = {};

for j = 1:numModels
    model = changeObjective(allModels.(modelNames{j}), objFun);
    fprintf(strcat(' \nAnalysing model: \n', modelNames{j}, '\n'))

    [~, grRateKO, ~, ~, delRxn, fluxSolution] = singleRxnDeletion(model);

    delRxnSubsystems(:, 1) = model.rxns;
    % delRxnSubsystems(:,2) = model.subSystems;
    test = setdiff(model.rxns, delRxn);
    if ~isempty(test)
        fprintf('Warning: different list of reactions found for:')
        modelNames{j};
    end
    sumRxnSubsystems = vertcat(sumRxnSubsystems, delRxnSubsystems);
    dataStruct.(modelNames{j}).rxnSubsystems = delRxnSubsystems;
    dataStruct.(modelNames{j}).grRateKO = grRateKO;
    dataStruct.(modelNames{j}).fluxSolution = fluxSolution;
    clear model
    clear delRxnSubsystems
end

% Identify unique reactions across all models
uniqueRxns = unique(sumRxnSubsystems(:, 1));
allRxns = {};
for i = 1:size(uniqueRxns, 1)
    allRxns(i, 1) = sumRxnSubsystems(find(strcmp(uniqueRxns{i, 1}, sumRxnSubsystems(:, 1)), 1), :);
end

% Build essential reaction table for all models
essentialRxn4Models = cell2table(allRxns, 'VariableNames', {'rxn'});  % ,'subsystem'

for j = 1:size(modelNames, 1)
    for i = 1:size(allRxns, 1)
        idx = find(strcmp(allRxns{i, 1}, dataStruct.(modelNames{j}).rxnSubsystems(:, 1)));
        if idx ~= 0
            essentialRxn4Models.(modelNames{j}){i} = dataStruct.(modelNames{j}).grRateKO(idx, 1);
        else
            essentialRxn4Models.(modelNames{j}){i} = 'NotIncluded';
        end
    end
end

end

author: opencobra

repo: cobratoolbox


text:
% COMPUTATION TIME FIGURES
function [currentFigure,figureHandle] = GetFigure_computationTimes(SIM,objectIndex,DATA,currentFigure)
% This function gets the computation time figures for all agents using the
% objectIndex.DATA fields.
% INPUTS:
% SIM     - Local copy of the META structure
% OBJECTS - The object class objects.
% currentFigure - The current figure number
% OUTPUTS:
% currentFigure - Updated figure number
% figureHandle - Handle to the created figure

figureHandle = [];

% SANITY CHECK 1
if DATA.totalAgents == 0
    warning('No agent data available.');
    return
else
    % LOGICALLY SELECT THE AGENTS FROM THE OBJECT INDEX
    agentObjectIndex = objectIndex([SIM.OBJECTS.type] == OMAS_objectType.agent); % The agents themselves
end

iter = 0;
for ID1 = 1:SIM.totalAgents
    % TRY TO GET THE REQUIRED PROPERTIES
    try
       agentObjectIndex{ID1}.DATA.dt;
    catch
       iter = iter + 1;
    end
end

% SANITY CHECK 2
if iter == SIM.totalAgents
    warning('No agents have computational timeseries data stored in (agent)."DATA.dt".');
    return
end

% CONFIGURE THE PLOT ATTRIBUTES
figurePath = strcat(SIM.outputPath,'computation-times');
figureHandle = figure('Name','OpenMAS computation timeseries');
ax = axes(figureHandle);
set(figureHandle,'Position',DATA.figureProperties.windowSettings);         % [x y width height]
set(figureHandle,'Color',DATA.figureProperties.figureColor);               % Background colour
setappdata(figureHandle, 'SubplotDefaultAxesLocation', [0.08, 0.08, 0.90, 0.88]);
hold on;
legendEntries = cell(SIM.totalAgents,1);
for ind = 1:SIM.totalAgents
    % Get object data
    entity = objectIndex{ind};
    % Get the equivalent SIM object
    SIMobject = SIM.OBJECTS([SIM.OBJECTS.objectID] == entity.objectID);
    % Legend entry
    legendEntries{ind} = sprintf('[ID-%s] %s',num2str(SIMobject.objectID),SIMobject.name);
    % Generate the figures
    stairs(ax,DATA.timeVector(1:SIM.TIME.endStep),entity.DATA.dt*1E3,...
         'Color',SIMobject.colour,...
         'LineWidth',DATA.figureProperties.lineWidth,...
         'DisplayName',legendEntries{ind});
end

% Title
title(ax,...
    sprintf('Agent computation times over a period of %ss',num2str(SIM.TIME.endTime)),...
    'Interpreter',DATA.figureProperties.interpreter,...
    'fontname',DATA.figureProperties.fontName,...
    'Fontweight',DATA.figureProperties.fontWeight,...
    'FontSize',DATA.figureProperties.titleFontSize,...
    'FontSmoothing','on');
% X-Label
xlabel(ax,'t (s)',...
    'Interpreter',DATA.figureProperties.interpreter,...
    'fontname',DATA.figureProperties.fontName,...
    'Fontweight',DATA.figureProperties.fontWeight,...
    'FontSize',DATA.figureProperties.axisFontSize,...
    'FontSmoothing','on');
% Y-Label
ylabel(ax,'Computation time (ms)',...
    'Interpreter',DATA.figureProperties.interpreter,...
    'fontname',DATA.figureProperties.fontName,...
    'Fontweight',DATA.figureProperties.fontWeight,...
    'FontSize',DATA.figureProperties.axisFontSize,...
    'FontSmoothing','on');
% Legend
legend(...
    'Location','northeast',...
    'Interpreter',DATA.figureProperties.interpreter,...
    'fontname',DATA.figureProperties.fontName);
% Axes
set(ax,...
    'TickLabelInterpreter',DATA.figureProperties.interpreter,...
    'fontname',DATA.figureProperties.fontName,...
    'Fontweight',DATA.figureProperties.fontWeight,...
    'FontSize',DATA.figureProperties.axisFontSize,...
    'FontSmoothing','on',...
    'Color',DATA.figureProperties.axesColor,...
    'GridLineStyle','--',...
    'GridAlpha',0.25,...
    'GridColor','k');
grid on; box on; grid minor;

% Show the timestep difference in the figure
ax.XAxis.MinorTickValues = ax.XAxis.Limits(1):SIM.TIME.dt:ax.XAxis.Limits(2);

% SAVE THE OUTPUT FIGURE
savefig(figureHandle,figurePath);

% Publish as .pdf if requested
if DATA.figureProperties.publish
	GetFigurePDF(figureHandle,figurePath);
end

% INCREMENT THE FIGURE INDEX
DATA.figureProperties.alignment = DATA.figureProperties.alignment + DATA.figureProperties.spacing;
currentFigure = currentFigure + 1;
end

repo: OpenMAS


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

% figure(1);
% path = points_to_idx(map, path);

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

hold on;
for i = 1:size(map.blocks,1)
    block = map.blocks(i, :);

    x = [ones(4,1) * block(1); ones(4,1) * block(4)];
    y = [ones(2,1) * block(5); ones(2,1) * block(2); ones(2,1) * block(5); ones(2,1) * block(2)];
    z = [block(3);block(6);block(3);block(6);block(3);block(6);block(3);block(6)];

    vert = [x, y, z];
    fac = [1 2 6 5;2 3 7 6;3 4 8 7;4 1 5 8;1 2 3 4;5 6 7 8];
    c = block(7:9)/255;
    patch('Vertices',vert,'Faces',fac,...
          'FaceVertexCData',hsv(6),'FaceColor',c);

    x = [ones(4,1) * block(1); ones(4,1) * block(4)];
    y = [block(2);block(5);block(2);block(5);block(2);block(5);block(2);block(5)];
    z = [ones(2,1) * block(3); ones(2,1) * block(6); ones(2,1) * block(3); ones(2,1) * block(6)];

    vert = [x, y, z];
    fac = [1 2 6 5;2 3 7 6;3 4 8 7;4 1 5 8;1 2 3 4;5 6 7 8];
    c = block(7:9)/255;
    patch('Vertices',vert,'Faces',fac,...
          'FaceVertexCData',hsv(6),'FaceColor',c);
end

if size(path,1) > 0
% pcshow(path, [0,0,0],'MarkerSize', 0.1);
end
% axis([map.boundary(1)-1, map.boundary(4)-1, map.boundary(2)-1,map.boundary(5)+1,map.boundary(3)+1,map.boundary(6)+1])
hold off;
% view(3);
end

author: yorgoon

repo: minimum-snap-geometric-control


text:
%+========================================================================+
%|                                                                        |
%|            This script uses the GYPSILAB toolbox for Matlab            |
%|                                                                        |
%| COPYRIGHT : Matthieu Aussal (c) 2017-2018.                             |
%| PROPERTY  : Centre de Mathematiques Appliquees, Ecole polytechnique,   |
%| route de Saclay, 91128 Palaiseau, France. All rights reserved.         |
%| LICENCE   : This program is free software, distributed in the hope that|
%| it will be useful, but WITHOUT ANY WARRANTY. Natively, you can use,    |
%| redistribute and/or modify it under the terms of the GNU General Public|
%| License, as published by the Free Software Foundation (version 3 or    |
%| later,  http://www.gnu.org/licenses). For private use, dual licencing  |
%| is available, please contact us to activate a "pay for remove" option. |
%| CONTACT   : matthieu.aussal@polytechnique.edu                          |
%| WEBSITE   : www.cmap.polytechnique.fr/~aussal/gypsilab                 |
%|                                                                        |
%| Please acknowledge the gypsilab toolbox in programs or publications in |
%| which you use it.                                                      |
%|________________________________________________________________________|
%|   '&`   |                                                              |
%|    #    |   FILE       : nrtRayTheatre.m                               |
%|    #    |   VERSION    : 0.41                                          |
%|   _#_   |   AUTHOR(S)  : Matthieu Aussal                               |
%|  ( # )  |   CREATION   : 14.03.2017                                    |
%|  / 0 \  |   LAST MODIF : 01.04.2018                                    |
%| ( === ) |   SYNOPSIS   : Ray tracing with theatre                      |
%|  `---'  |                                                              |
%+========================================================================+

% Cleaning
clear all
close all
clc

% Gypsilab path
run('../../addpathGypsilab.m')

% Parameters
Xsrc  = [0 8 1.7];
Xmes  = [0 -20 8];
Nray  = 1e5
rad   = 1

% Read mesh
mesh        = msh('theatre.mesh');
mesh.col(:) = 100; % Rough concrete (Bobran, 1973)

% Extract toit
ctr = mesh.ctr;
mesh = mesh.sub(ctr(:,3)<14);

% Initialize ray
ray1 = ray(mesh,Xsrc,Nray);

% Graphical sphere
[X,Y,Z] = sphere(50);
X = Xmes(1) + rad*X; Y = Xmes(2) + rad*Y; Z = Xmes(3) + rad*Z;

% Graphical representation
figure
plot(mesh,'w')
hold on
plot(ray1)
surf(X,Y,Z,ones(size(X)))
axis equal
xlabel('X');   ylabel('Y');   zlabel('Z');
view(0,45)
alpha(0.5)

% Maximum distances
r1000 = rad/2 * sqrt(Nray/1000);
r100  = rad/2 * sqrt(Nray/100);
r10   = rad/2 * sqrt(Nray/10);
rMax  = rad/2 * sqrt(Nray/2);

% Ray-tracing
tic
ray1 = ray1.tracer(30,rMax);
toc
% plot(ray1)

% Images sources
tic
[img,nrg] = ray1.image(Xmes,rad,rMax);
toc

% Data
r   = sqrt(sum(img.^2,2));
sol = mean(nrg,2);

% Impacts
ray2 = ray(mesh,Xmes,img);
ray2 = ray2.tracer;
plot(ray2)
plot(msh(ray2.pos{2},(1:length(ray2))',10*log10(sol)));
axis equal
colorbar

% Energy in dB
figure
plot(r,10*log10(sol),'+b')
hold on
plot([r1000 r1000],[-60,0],'k--')
text(r1000,-55,' n = 1000')
plot([r100 r100],[-60,0],'k--')
text(r100,-57,' n = 100')
plot([r10 r10],[-60,0],'k--')
text(r10,-55,' n = 10')
plot([rMax rMax],[-60,0],'k--')
text(rMax,-55,' n = 1')
xlabel('Distance source mesure')
ylabel('Energie mesuree (dB)')
title('Energie mesuree selon la distance de mesure')
grid on

% Audio file
[audio,fs] = audioread('anechoicVoice.wav');

% Fir from images
T = floor(r/340*fs) + 1;
for i = 1:size(nrg,2)
    fir8(:,i) = accumarray(T,nrg(:,i),[max(T) 1]);
end

% Bank filtering
dir8 = ray1.bank(256,fs);
rir  = 0;
for i = 1:size(dir8,2)
    rir = rir + fftfilt(dir8(:,i),fir8(:,i));
end

% Graphical representation
figure
plot(rir)

% Audio rendering (5 seconds)
out = fftfilt(rir,audio(1:5*fs));
% sound(out,fs)

disp('~~> Michto gypsilab !')

author: SwanLab

repo: Swan

save_path: github-repos/MATLAB/SwanLab-Swan


text:
%  Figure 7.53      Feedback Control of Dynamic Systems, 5e
%                        Franklin, Powell, Emami
%
% script for fig. 7.53
%
clf;
dp=[1 1 0];
np=[1];
nc=conv([1 1.001],[8.32 0.8]);
dc=conv([1 4.08],[1 0.0196]);
num=conv(np,nc);
den=conv(dp,dc);
dcl=[0 0 num]+den; % closed-loop denominator
t=0:.1:5;
y=step(num,dcl,t);
plot(t,y);
xlabel('Time (sec)');
ylabel('y(t)');
grid;
title('Fig.7.53: Step response for lag compensation design')


text:
function Y = sym_sne(X, d, perplexity)
%SNE Implementation of symmetric Stochastic Neighbor Embedding
%
%   mappedX = sym_sne(X, no_dims, perplexity)
%
% Runs the symmetric Stochastic Neighbor Embedding algorithm. The
% high-dimensional datapoints are specified by X. The target dimensionality
% is specified in no_dims (default = 2), and the perplexity of the Gaussian
% kernel is specified in perplexity (default = 30).
% The function returns the embedded points in mappedX.
%

% This file is part of the Matlab Toolbox for Dimensionality Reduction.
% The toolbox can be obtained from http://homepage.tudelft.nl/19j49
% You are free to use, change, or redistribute this code in any way you
% want for non-commercial purposes. However, it is appreciated if you
% maintain the name of the original author.
%
% (C) Laurens van der Maaten, Delft University of Technology

    if ~exist('d', 'var') || isempty(d)
        d = 2;
    end
    if ~exist('perplexity', 'var') || isempty(perplexity)
        perplexity = 30;
    end

    % Initialize some variables
    n = size(X, 1);                 % number of instances
    eta = .05;                      % learning rate
    max_iter = 2000;                % maximum number of iterations
    jitter = 0.3;                   % initial jitter
    jitter_decay = 0.99;            % jitter decay
    momentum = 0.5;                 % initial momentum
    final_momentum = 0.8;           % final momentum
    mom_switch_iter = 750;          % iteration where momentum changes

    % Initialize embedding coordinates randomly (close to origin)
    Y = 0.0001 * rand(n, d);
    dC = zeros(n, d);
    y_incs = zeros(n, d);

    % Compute Gaussian kernel for high-dimensional data representation
    P = x2p(X, perplexity, 1e-5);                                        % use fixed perplexity
    P = P + P';
    P = P ./ sum(P(:));
    P = max(P, eps);

    % Iterating loop
    for iter=1:max_iter

        % Compute Gaussian kernel for low-dimensional data representation
        sum_Y = sum(Y .^ 2, 2);                                                     % precomputation for pairwise distances
        Q = exp(-bsxfun(@plus, sum_Y, bsxfun(@plus, sum_Y', -2 * Y * Y')) ./ 2);    % Gaussian probabilities
        Q = Q ./ repmat(sum(Q, 2), [1 n]);
        Q = max(Q, eps);

        % Compute cost function between P and Q
        if ~rem(iter, 20)
            costs = sum(P .* log((P + eps) ./ (Q + eps)), 2) ./ n;              % division by n corrects for # of datapoints
            cost = sum(costs);
            disp(['Iteration ' num2str(iter) ': error is ' num2str(cost)]);
        end

        % Compute gradient
        PQ = P - Q;
        for i=1:n
            dC(i,:) = sum(bsxfun(@times, bsxfun(@minus, Y(i,:), Y), PQ(i,:)'), 1);
        end

        % Perform the gradient search
        y_incs = momentum * y_incs - eta * dC;
        Y = Y + y_incs;
		Y = Y + jitter * randn(size(Y));
        Y = bsxfun(@minus, Y, mean(Y, 1));

        % Reduce jitter over time and change momentum
        jitter = jitter * jitter_decay;
        if iter == mom_switch_iter
            momentum = final_momentum;
        end
    end

repo: cluster


text:
function d = fd03 ( p )

%*****************************************************************************80
%
%% FD03 is a signed distance function for problem 3.
%
%  Discussion:
%
%    The formula used here is not quite correct.  In particular, it is wrong
%    for points exterior to the cube whose nearest point on the cube is at a corner.
%
%    For DISTMESH_3D's purposes, though, this computation is accurate enough.
%
%  Modified:
%
%    12 September 2005
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, real P(N,3), one or more points.
%
%    Output, real D(N), the signed distance of each point to the boundary of the region.
%
  d = - min ( min ( min ( min ( min ( -0.0+p(:,3), 1.0-p(:,3) ), ...
                                     -0.0+p(:,2) ), ...
                                      1.0-p(:,2) ), ...
                                     -0.0+p(:,1) ), ...
                                      1.0-p(:,1) );

  return
end


text:
function [ mO ] = CompressImageSvd( mI, energyThr, blockRadius )
% ----------------------------------------------------------------------------------------------- %
% [ mO ] = CompressImageSvd( mI, energyThr, blockRadius )
%   Compresses the image using SVD.
% Input:
%   - mI            -   Input Image.
%                       Structure: Image Matrix (Signle Channel or RGB).
%                       Type: 'Single' / 'Double'.
%                       Range: [0, 1].
%   - energyThr     -   Energy Threshold.
%                       Sets the threshold for Singular Value kept energy.
%                       Structure: Scalar.
%                       Type: 'Single' / 'Double'.
%                       Range: [0, 1].
%   - blockRadius   -   Block Radius.
%                       Sets the block radius for the compression.
%                       Structure: Scalar.
%                       Type: 'Single' / 'Double'.
%                       Range: {1, 2, ...}.
% Output:
%   - mO            -   Output Image.
%                       Structure: Image Matrix (Signle Channel or RGB).
%                       Type: 'Single' / 'Double'.
%                       Range: [0, 1].
% References
%   1.  SVD Wikipedia - https://en.wikipedia.org/wiki/Singular_value_decomposition.
% Remarks:
%   1.  a
% TODO:
%   1.  U.
% Release Notes:
%   -   1.0.000     01/09/2017  Royi
%       *   First release version.
% ----------------------------------------------------------------------------------------------- %

FALSE   = 0;
TRUE    = 1;

OFF     = 0;
ON      = 1;

numRows = size(mI, 1);
numCols = size(mI, 2);
numChan = size(mI, 3); %<! Num Channels

vImageDim = [numRows, numCols];

blockLength = (2 * blockRadius) + 1;
vBlockDim   = [blockLength, blockLength];

mO = zeros([numRows, numCols, numChan]);

for ii = 1:numChan

    mII     = mI(:, :, ii);
    dcLevel = mean(mII(:)); %<! Extracting DC Level
    mII     = mII - dcLevel;

    % Decomposing the image into blocks. Each block becomes a vector in the
    % Columns Images.
    mColImage   = im2col(mII, vBlockDim, 'distinct');

    % The SVD Step
    [mU, mS, mV] = svd(mColImage);

    vSingularValues = diag(mS);

    vSingularValueEnergy = cumsum(vSingularValues) / sum(vSingularValues);
    lastIdx = find(vSingularValueEnergy >= energyThr, 1, 'first');

    vSingularValues(lastIdx + 1:end) = 0;
    % mS isn't necessarily square matrix. Hence only work on its main
    % diagonal.
    mS(1:length(vSingularValues), 1:length(vSingularValues)) = diag(vSingularValues);

    % Reconstruction of the image using "Less Energy".
    mColImage = mU * mS * mV.';

    % Restorig the original structure and the DC Level
    mO(:, :, ii) = col2im(mColImage, vBlockDim, vImageDim, 'distinct') + dcLevel;

end

end

author: RoyiAvital

repo: StackExchangeCodes


text:
%
% D = gpdist(x1, x2, funcDist)
%
% Returns matrix of pair-wise distances between x1 and x2 as a matrix.
% x1 is K x M matrix
% x2 is K x N matrix
% where K is the dimensionality, M and N are the number of samples.
% The resulting matrix D is M x N matrix.
% If no function is given for distance calculation, 2-norm is used.
function D = gpdist(x1, x2, funcDist)

if nargin < 2
  x2 = x1;
end
if isvector(x1)
  x1 = x1(:)';
  x2 = x2(:)';
end
if nargin < 3
  funcDist = @(z1,z2) norm(z1-z2);
%  funcDist = @(z1,z2) abs(z1-z2);
end

n1 = size(x1,2);
n2 = size(x2,2);
D = zeros(n1,n2);
for i=1:n1
  for j=1:n2
    D(i,j) = funcDist(x1(:,i),x2(:,j));
  end
end

author: jluttine

repo: matlab

save_path: github-repos/MATLAB/jluttine-matlab


text:
function mldate = x2mdate(xlsdate, type)
% X2MDATE provides a simple method to convert between excel dates and MATLAB dates.
%
% USAGE:
%   [MLDATE] = x2mdate(XLSDATE)
%   [MLDATE] = x2mdate(XLSDATE, TYPE)
%
% INPUTS:
%   XLSDATE   - A scalar or vector of Excel dates.
%   TYPE      - [OPTIONAL] A scalar or vector of the same size as XLSDATE that describes the Excel
%                 basedate.  Can be either 0 or 1.  If 0 (default), the base date of Dec-31-1899 is
%                 used.  If 1, the base date is Jan 1, 1904.
%
% OUTPUTS:
%   MLDATE    - A vector with the same size as XLSDATE consisting of MATLAB dates.
%
% EXAMPLE:
%   XLSDATE = [35000 40000 41000];
%   MLDATE  = x2mdate(XLSDATE);
%   datestr(MLDATE)
%       28-Oct-1995
%       06-Jul-2009
%       01-Apr-2012
%
% COMMENTS:
%   This is a reverse engineered clone of the MATLAB function x2mdate and should behave the same.
%   You only need it if you do not have the financial toolbox installed.
%
% See also C2MDATE

% Copyright: Kevin Sheppard
% kevin.sheppard@economics.ox.ac.uk
% Revision: 1    Date: 10/27/2006

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Input Validation
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if nargin==1
    type=0;
end

if any(ischar(xlsdate))
    error('XLSDATE must be numeric')
end

if ndims('xlsdate')~=2 &&  min(size(type))~=1
    error('XLSDATE must be a T by 1 or 1 by T vector');
end

if nargin>1
    if nargin>2
        error('1 or 2 inputs only');
    end
    if min(size(type))~=1 || ndims(type)~=2
        error('TYPE must be either a scalar or a vector conformable to xlsdate');
    end
    if max(size(type))~=1
        if max(size(type))~=max(size(xlsdate))
            error('TYPE must be either a scalar or a vector conformable to xlsdate');
        end
    end
    if any(~ismember(type,[0 1]))
        error('TYPE must be either 0 or 1.')
    end
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Input Validation
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if isscalar(type)
    type=type*ones(size(xlsdate));
end
type = logical(type);

mldate=zeros(size(xlsdate));
mldate(~type) = xlsdate(~type) + datenum('30-Dec-1899') ;
mldate(type) = xlsdate(type) + datenum('1-Jan-1904') ;

author: bashtage

repo: mfe-toolbox


text:
function [ x ] = hyperdistance( geo1, geo2 )
%HYPERDISTANCE Calculate hyperdistance between two geometry structures.
%   x = HYPERDISTANCE(geometry1,geometry2) evaluates the hyperdistance
%   between two geometries. The hyperdistance is defined as the root mean
%   square of the displacements of each atom. The returned value is only
%   meaningful if the two geometries have the same lattice vectors. The
%   minimum image convention is used to calculate atomic displacements. See
%   IMPORT_POSCAR for a detailed description of the geometry structure.
%
%   See also IMPORT_POSCAR.

       x = sqrt(sum(sum(((mod(abs(geo2.coords-geo1.coords)+0.5,1)-0.5)*geo1.lattice).^2)));

end

author: Sable

repo: mcbench-benchmarks


text:
function rf=lpcis2rf(is)
%LPCRF2IS Convert inverse sines to reflection coefficients RF=(IS)

%      Copyright (C) Mike Brookes 1997
%      Version: $Id: lpcis2rf.m,v 1.4 2007/05/04 07:01:38 dmb Exp $
%
%   VOICEBOX is a MATLAB toolbox for speech processing.
%   Home page: http://www.ee.ic.ac.uk/hp/staff/dmb/voicebox/voicebox.html
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%   This program is free software; you can redistribute it and/or modify
%   it under the terms of the GNU General Public License as published by
%   the Free Software Foundation; either version 2 of the License, or
%   (at your option) any later version.
%
%   This program is distributed in the hope that it will be useful,
%   but WITHOUT ANY WARRANTY; without even the implied warranty of
%   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
%   GNU General Public License for more details.
%
%   You can obtain a copy of the GNU General Public License from
%   http://www.gnu.org/copyleft/gpl.html or by writing to
%   Free Software Foundation, Inc.,675 Mass Ave, Cambridge, MA 02139, USA.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

rf=sin(is*pi/2);

author: decouples


text:
function Rt = cameraAngleAxis2Rt(vec)

Rt = AngleAxisToRotationMatrix(vec(1:3));
Rt(1:3,4) = vec(4:6);


text:
function DebugBatch(mfile,bfile,item,varargin)

%DebugBatch - Help debug a batch job.
%
%  This function can be used to help debug a batch. It retrieves the parameters
%  corresponding to the required item from the batch file, and either starts
%  the batch function with these parameters in debug mode, or simply assigns
%  these parameters to new variables in matlab's workspace so they can be
%  directly manipulated.
%
%  Set breakpoints in the batch function before running DebugBatch. If no
%  breakpoints have been defined, DebugBatch will stop at the first line.
%
%  USAGE
%
%    DebugBatch(mfile,bfile,item,<options>)
%
%    mfile          batch function (M-file name or function handle)
%    bfile          batch file listing the parameters for each iteration
%    item           item number in batch function
%    <options>      optional list of property-value pairs (see table below)
%
%    =========================================================================
%     Properties    Values
%    -------------------------------------------------------------------------
%     'mode'        run the batch function ('run', default) or simply set
%                   variables in matlab's workspace ('set')
%    =========================================================================
%
%  OUTPUT
%
%    Instantiates the variables of the batch function using the values in the
%    batch file ('set' mode).
%
%  SEE
%
%    See also StartBatch.
%

% Copyright (C) 2007-2013 by Michaël Zugaro
%
% This program is free software; you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation; either version 3 of the License, or
% (at your option) any later version.

% Default values
mode = 'run';

% Check number of parameters
if nargin < 3,
	error(['Incorrect number of parameters (type ''help DebugBatch'' for details).']);
end

% Batch function name
if isa(mfile,'function_handle'),
	mfile = func2str(mfile);
end

% Check batch file and function are valid
if ~isstring(bfile) || ~exist(bfile,'file'),
	error(['Batch file not found (type ''help DebugBatch'' for details).']);
end
if isempty(which(mfile)),
	error(['Batch function not found (type ''help DebugBatch'' for details).']);
end

% Item
if ~isiscalar(item),
	error(['Incorrect item (type ''help DebugBatch'' for details).']);
end

% Parse parameter list
for i = 1:2:length(varargin),
  if ~ischar(varargin{i}),
    error(['Parameter ' num2str(i+1) ' is not a property (type ''help DebugBatch'' for details).']);
  end
  switch(lower(varargin{i})),
    case 'mode',
      mode = lower(varargin{i+1});
      if ~isstring(mode,'run','set'),
        error('Incorrect value for property ''mode'' (type ''help DebugBatch'' for details).');
      end
    otherwise,
      error(['Unknown property ''' num2str(varargin{i}) ''' (type ''help DebugBatch'' for details).']);
  end
end

% Open batch file
f = fopen(bfile,'r');
if f == -1, error(['Could not open file ''' bfile '''.']); end

% Parse batch file
b = ParseBatch(bfile);

if strcmp(mode,'run'),

	% Check breakpoints
	if length(dbstatus(mfile)) == 0,
		eval(['dbstop in ' mfile ' at 1;']);
	end

	% Run batch function with appropriate parameters
	for i = 1:size(b.field,2),
		args{i} = b.field{item,i};
	end
	outputs = cell(1,nargout(mfile));
	[outputs{:}] = feval(mfile,args{:})

else

	% Open batch function
	if isa(mfile,'function_handle'),
		mfile = func2str(mfile);
	end
	f = fopen(which(mfile));

	% Find first line containing the (uncommented) 'function' keyword
	found = [];
	while isempty(found),
		line = fgets(f);
		if line == -1, error('Could not find function definition in batch function.'); end
		line = regexprep(line,'%.*function.*','');
		found = regexp(line,'.*function[^(]*\(([^)]*)\).*');
	end
	fclose(f);

	% Extract parameter names, and assign them in 'base' workspace
	parameters = regexprep(line,'.*function[^(]*\(([^)]*)\).*','$1');
	parameters = regexp(parameters,'[^,]*','match');
	for i = 1:length(parameters),
		assignin('base',parameters{i},b.field{item,i});
	end

end

author: buzsakilab

repo: buzcode

save_path: github-repos/MATLAB/buzsakilab-buzcode


text:
function test_ft_getopt

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

success = true;

success = success && isequal(ft_getopt({'key1', 'val1'}, 'key1'), 'val1');
success = success && isequal(ft_getopt({'key1', 'val1'}, 'key2'), []);
success = success && isequal(ft_getopt({'key1', 'val1'}, 'key2', 'default'), 'default');
success = success && isequal(ft_getopt({'key1', 'val1'}, 'key2', 'default'), 'default');

cfg      = [];
cfg.key1 = 'val1';
success = success && isequal(ft_getopt(cfg, 'key1'), 'val1');
success = success && isequal(ft_getopt(cfg, 'key2'), []);
success = success && isequal(ft_getopt(cfg, 'key2', 'default'), 'default');
success = success && isequal(ft_getopt(cfg, 'key2', 'default'), 'default');
success = success && isequal(ft_getopt({'key1', []}, 'key1'), []);
success = success && isequal(ft_getopt({'key1', []}, 'key1', 'default'), 'default');

cfg      = [];
cfg.key1 = [];
success = success && isequal(ft_getopt(cfg, 'key1'), []);
success = success && isequal(ft_getopt(cfg, 'key1', 'default'), 'default');
success = success && isequal(ft_getopt({'key1', 'val1', 'key2', 'val2', 'key3', 'val3'}, 'key3'), 'val3');
success = success && isequal(ft_getopt({'key1', 'val1', 'key2', 'val2', 'key3', 'val3'}, 'key4'), []);
success = success && isequal(ft_getopt({'key1', 'val1', 'key2', 'val2', 'key3', 'val3'}, 'key4', 'default'), 'default');

cfg      = [];
cfg.key1 = 'val1';
cfg.key2 = 'val2';
cfg.key3 = 'val3';
success = success && isequal(ft_getopt(cfg, 'key3'), 'val3');
success = success && isequal(ft_getopt(cfg, 'key4'), []);
success = success && isequal(ft_getopt(cfg, 'key4', 'default'), 'default');

% test the emptymeaningful option
success = success &&  isempty(ft_getopt({'key',[]},'key','default',1));
success = success && ~isempty(ft_getopt({'key',[]},'key','default',0));

if ~success
  error('there is a problem with ft_getopt');
end

author: fieldtrip

repo: fieldtrip

save_path: github-repos/MATLAB/fieldtrip-fieldtrip


text:
function backDiffUy = backdiffy(u)
% backward diff in y direction in equation (2.9a)
% \Delta^y_{+}u_{ij} = u_{i,j+1} - u_{ij}

cols = size(u,2);
% u([2:rows rows], :) is u_{i+1,j} in (2.9a)
backDiffUy = u(:, [2:cols cols]) - u;

author: YimianDai

path: github-repos/MATLAB/YimianDai-Image-Processing-Codes-for-Easier-Understanding/Image-Processing-Codes-for-Easier-Understanding-874302799e48852624bc3760b58b46bd9360f238/src/(Physica D 1992) Nonlinear Total Variation based noise removal algorithms/version_1/support/backdiffy.m


text:
%% housekeeping
close all
clear
clc
%% load the data

load tut01_data

load tut02_estimation

load tut09_bootstrap

load tut04_identification

ci=[30,50,68,90];

%% pick a model

model='ve_lr'; % ve

endog=models.(model).endogenous;

%% compute decompositions

% choose identification scheme
Rfunc=[];

hd=variance_decomposition(models.(model),params.(model),Rfunc,shock_names);

%% plot decompositions
shock_tex=shock_names;

myrange='1:50';

for iv=1:numel(endog)

    vname=tex.(endog{iv});

    titel=['Variance Decomposition (in %) of ',vname];

    figure('name',titel);

    d=hd.conditional.(endog{iv});

    d.varnames=shock_names;

    d=pages2struct(d);

    contributors=fieldnames(d); % = shock_names

    for ii=1:numel(contributors)

        subplot(3,2,ii)

        % note we are multiplying by 100, this just by pure convenience
        %--------------------------------------------------------------
        out=fanchart(100*d.(contributors{ii})(myrange),ci);

        plot_fanchart(out)

        title(contributors{ii})

        axis tight

    end

    [~,h]=sup_label(titel,'t');

    set(h,'fontsize',12)

end

author: jmaih


text:
function obj = Jacob_pnp_func_new(in1,in2,in3,in4,in5)
coef_f0_q_sym1 = in1(:,1);
coef_f0_q_sym2 = in1(:,2);
coef_f1_q_sym1 = in2(:,1);
coef_f0_q_sym3 = in1(:,3);
coef_f1_q_sym2 = in2(:,2);
coef_f2_q_sym1 = in3(:,1);
coef_f0_q_sym4 = in1(:,4);
coef_f1_q_sym3 = in2(:,3);
coef_f2_q_sym2 = in3(:,2);
coef_f3_q_sym1 = in4(:,1);
coef_f0_q_sym5 = in1(:,5);
coef_f1_q_sym4 = in2(:,4);
coef_f2_q_sym3 = in3(:,3);
coef_f3_q_sym2 = in4(:,2);
coef_f0_q_sym6 = in1(:,6);
coef_f1_q_sym5 = in2(:,5);
coef_f2_q_sym4 = in3(:,4);
coef_f3_q_sym3 = in4(:,3);
coef_f0_q_sym7 = in1(:,7);
coef_f1_q_sym6 = in2(:,6);
coef_f2_q_sym5 = in3(:,5);
coef_f3_q_sym4 = in4(:,4);
coef_f0_q_sym8 = in1(:,8);
coef_f1_q_sym7 = in2(:,7);
coef_f2_q_sym6 = in3(:,6);
coef_f3_q_sym5 = in4(:,5);
coef_f0_q_sym9 = in1(:,9);
coef_f1_q_sym8 = in2(:,8);
coef_f2_q_sym7 = in3(:,7);
coef_f3_q_sym6 = in4(:,6);
coef_f1_q_sym9 = in2(:,9);
coef_f2_q_sym8 = in3(:,8);
coef_f3_q_sym7 = in4(:,7);
coef_f2_q_sym9 = in3(:,9);
coef_f3_q_sym8 = in4(:,8);
coef_f3_q_sym9 = in4(:,9);
coef_f0_q_sym10 = in1(:,10);
coef_f0_q_sym11 = in1(:,11);
coef_f0_q_sym20 = in1(:,20);
coef_f1_q_sym10 = in2(:,10);
coef_f0_q_sym12 = in1(:,12);
coef_f0_q_sym21 = in1(:,21);
coef_f1_q_sym11 = in2(:,11);
coef_f1_q_sym20 = in2(:,20);
coef_f2_q_sym10 = in3(:,10);
coef_f0_q_sym13 = in1(:,13);
coef_f0_q_sym22 = in1(:,22);
coef_f1_q_sym12 = in2(:,12);
coef_f1_q_sym21 = in2(:,21);
coef_f2_q_sym11 = in3(:,11);
coef_f2_q_sym20 = in3(:,20);
coef_f3_q_sym10 = in4(:,10);
coef_f0_q_sym14 = in1(:,14);
coef_f0_q_sym23 = in1(:,23);
coef_f1_q_sym13 = in2(:,13);
coef_f1_q_sym22 = in2(:,22);
coef_f2_q_sym12 = in3(:,12);
coef_f2_q_sym21 = in3(:,21);
coef_f3_q_sym11 = in4(:,11);
coef_f3_q_sym20 = in4(:,20);
coef_f0_q_sym15 = in1(:,15);
coef_f0_q_sym24 = in1(:,24);
coef_f1_q_sym14 = in2(:,14);
coef_f1_q_sym23 = in2(:,23);
coef_f2_q_sym13 = in3(:,13);
coef_f2_q_sym22 = in3(:,22);
coef_f3_q_sym12 = in4(:,12);
coef_f3_q_sym21 = in4(:,21);
coef_f0_q_sym16 = in1(:,16);
coef_f1_q_sym15 = in2(:,15);
coef_f1_q_sym24 = in2(:,24);
coef_f2_q_sym14 = in3(:,14);
coef_f2_q_sym23 = in3(:,23);
coef_f3_q_sym13 = in4(:,13);
coef_f3_q_sym22 = in4(:,22);
coef_f0_q_sym17 = in1(:,17);
coef_f1_q_sym16 = in2(:,16);
coef_f2_q_sym15 = in3(:,15);
coef_f2_q_sym24 = in3(:,24);
coef_f3_q_sym14 = in4(:,14);
coef_f3_q_sym23 = in4(:,23);
coef_f0_q_sym18 = in1(:,18);
coef_f1_q_sym17 = in2(:,17);
coef_f2_q_sym16 = in3(:,16);
coef_f3_q_sym15 = in4(:,15);
coef_f3_q_sym24 = in4(:,24);
coef_f0_q_sym19 = in1(:,19);
coef_f1_q_sym18 = in2(:,18);
coef_f2_q_sym17 = in3(:,17);
coef_f3_q_sym16 = in4(:,16);
coef_f1_q_sym19 = in2(:,19);
coef_f2_q_sym18 = in3(:,18);
coef_f3_q_sym17 = in4(:,17);
coef_f2_q_sym19 = in3(:,19);
coef_f3_q_sym18 = in4(:,18);
coef_f3_q_sym19 = in4(:,19);
q0 = in5(1,:);
q1 = in5(2,:);
q2 = in5(3,:);
q3 = in5(4,:);
t114 = q0.^2;
t115 = q1.^2;
t116 = q2.^2;
t117 = q3.^2;
t118 = coef_f0_q_sym11.*q0;
t119 = coef_f0_q_sym24.*q3;
t120 = coef_f0_q_sym1.*q0.*t114;
t121 = coef_f0_q_sym23.*q3.*t117;
t122 = coef_f0_q_sym4.*q3.*t114;
t123 = coef_f0_q_sym10.*q0.*t117;
t124 = coef_f0_q_sym6.*q0.*q1.*q2.*2.0;
t125 = coef_f0_q_sym16.*q1.*q2.*q3.*2.0;
t126 = coef_f0_q_sym18.*q1;
t127 = coef_f0_q_sym22.*q2;
t128 = coef_f0_q_sym12.*q1.*t115;
t129 = coef_f0_q_sym19.*q2.*t116;
t130 = coef_f0_q_sym2.*q1.*t114;
t131 = coef_f0_q_sym5.*q0.*t115;
t132 = coef_f0_q_sym3.*q2.*t114;
t133 = coef_f0_q_sym8.*q0.*t116;
t134 = coef_f0_q_sym7.*q0.*q1.*q3.*2.0;
t135 = coef_f0_q_sym9.*q0.*q2.*q3.*2.0;
obj = reshape([coef_f1_q_sym11.*q0.*-2.0+coef_f0_q_sym11.*q1-coef_f1_q_sym18.*q1-coef_f1_q_sym22.*q2-coef_f1_q_sym24.*q3-coef_f1_q_sym1.*q0.*t114.*4.0+coef_f0_q_sym2.*q0.*t115.*2.0-coef_f1_q_sym5.*q0.*t115.*2.0-coef_f1_q_sym8.*q0.*t116.*2.0-coef_f1_q_sym10.*q0.*t117.*2.0+coef_f0_q_sym1.*q1.*t114.*3.0-coef_f1_q_sym2.*q1.*t114.*3.0+coef_f0_q_sym5.*q1.*t115+coef_f0_q_sym8.*q1.*t116+coef_f0_q_sym10.*q1.*t117-coef_f1_q_sym12.*q1.*t115-coef_f1_q_sym15.*q1.*t116-coef_f1_q_sym17.*q1.*t117-coef_f1_q_sym3.*q2.*t114.*3.0+coef_f0_q_sym6.*q2.*t115-coef_f1_q_sym13.*q2.*t115-coef_f1_q_sym21.*q2.*t117-coef_f1_q_sym19.*q2.*t116-coef_f1_q_sym4.*q3.*t114.*3.0+coef_f0_q_sym7.*q3.*t115-coef_f1_q_sym20.*q3.*t116-coef_f1_q_sym14.*q3.*t115-coef_f1_q_sym23.*q3.*t117+coef_f0_q_sym3.*q0.*q1.*q2.*2.0-coef_f1_q_sym6.*q0.*q1.*q2.*2.0+coef_f0_q_sym4.*q0.*q1.*q3.*2.0-coef_f1_q_sym7.*q0.*q1.*q3.*2.0-coef_f1_q_sym9.*q0.*q2.*q3.*2.0+coef_f0_q_sym9.*q1.*q2.*q3-coef_f1_q_sym16.*q1.*q2.*q3,coef_f2_q_sym11.*q0.*-2.0-coef_f2_q_sym18.*q1+coef_f0_q_sym11.*q2-coef_f2_q_sym22.*q2-coef_f2_q_sym24.*q3-coef_f2_q_sym1.*q0.*t114.*4.0+coef_f0_q_sym3.*q0.*t116.*2.0-coef_f2_q_sym5.*q0.*t115.*2.0-coef_f2_q_sym8.*q0.*t116.*2.0-coef_f2_q_sym10.*q0.*t117.*2.0-coef_f2_q_sym2.*q1.*t114.*3.0+coef_f0_q_sym6.*q1.*t116-coef_f2_q_sym12.*q1.*t115-coef_f2_q_sym15.*q1.*t116-coef_f2_q_sym17.*q1.*t117+coef_f0_q_sym1.*q2.*t114.*3.0-coef_f2_q_sym3.*q2.*t114.*3.0+coef_f0_q_sym5.*q2.*t115+coef_f0_q_sym8.*q2.*t116+coef_f0_q_sym10.*q2.*t117-coef_f2_q_sym13.*q2.*t115-coef_f2_q_sym21.*q2.*t117-coef_f2_q_sym19.*q2.*t116-coef_f2_q_sym4.*q3.*t114.*3.0+coef_f0_q_sym9.*q3.*t116-coef_f2_q_sym20.*q3.*t116-coef_f2_q_sym14.*q3.*t115-coef_f2_q_sym23.*q3.*t117+coef_f0_q_sym2.*q0.*q1.*q2.*2.0-coef_f2_q_sym6.*q0.*q1.*q2.*2.0-coef_f2_q_sym7.*q0.*q1.*q3.*2.0+coef_f0_q_sym4.*q0.*q2.*q3.*2.0-coef_f2_q_sym9.*q0.*q2.*q3.*2.0+coef_f0_q_sym7.*q1.*q2.*q3-coef_f2_q_sym16.*q1.*q2.*q3,coef_f3_q_sym11.*q0.*-2.0-coef_f3_q_sym18.*q1-coef_f3_q_sym22.*q2+coef_f0_q_sym11.*q3-coef_f3_q_sym24.*q3-coef_f3_q_sym1.*q0.*t114.*4.0+coef_f0_q_sym4.*q0.*t117.*2.0-coef_f3_q_sym5.*q0.*t115.*2.0-coef_f3_q_sym8.*q0.*t116.*2.0-coef_f3_q_sym10.*q0.*t117.*2.0-coef_f3_q_sym2.*q1.*t114.*3.0+coef_f0_q_sym7.*q1.*t117-coef_f3_q_sym12.*q1.*t115-coef_f3_q_sym15.*q1.*t116-coef_f3_q_sym17.*q1.*t117-coef_f3_q_sym3.*q2.*t114.*3.0+coef_f0_q_sym9.*q2.*t117-coef_f3_q_sym13.*q2.*t115-coef_f3_q_sym21.*q2.*t117-coef_f3_q_sym19.*q2.*t116+coef_f0_q_sym1.*q3.*t114.*3.0+coef_f0_q_sym5.*q3.*t115-coef_f3_q_sym4.*q3.*t114.*3.0+coef_f0_q_sym8.*q3.*t116+coef_f0_q_sym10.*q3.*t117-coef_f3_q_sym20.*q3.*t116-coef_f3_q_sym14.*q3.*t115-coef_f3_q_sym23.*q3.*t117-coef_f3_q_sym6.*q0.*q1.*q2.*2.0+coef_f0_q_sym2.*q0.*q1.*q3.*2.0-coef_f3_q_sym7.*q0.*q1.*q3.*2.0+coef_f0_q_sym3.*q0.*q2.*q3.*2.0-coef_f3_q_sym9.*q0.*q2.*q3.*2.0+coef_f0_q_sym6.*q1.*q2.*q3-coef_f3_q_sym16.*q1.*q2.*q3,q0.*2.0,t118+t119+t120+t121+t122+t123+t124+t125+t127+t129+t132+t133+t134-coef_f1_q_sym18.*q0+coef_f0_q_sym18.*q1.*2.0-coef_f1_q_sym2.*q0.*t114+coef_f0_q_sym5.*q0.*t115.*3.0-coef_f1_q_sym12.*q0.*t115.*3.0-coef_f1_q_sym15.*q0.*t116-coef_f1_q_sym17.*q0.*t117+coef_f0_q_sym2.*q1.*t114.*2.0-coef_f1_q_sym5.*q1.*t114.*2.0+coef_f0_q_sym12.*q1.*t115.*4.0+coef_f0_q_sym15.*q1.*t116.*2.0+coef_f0_q_sym17.*q1.*t117.*2.0-coef_f1_q_sym6.*q2.*t114+coef_f0_q_sym13.*q2.*t115.*3.0+coef_f0_q_sym21.*q2.*t117-coef_f1_q_sym7.*q3.*t114+coef_f0_q_sym20.*q3.*t116+coef_f0_q_sym14.*q3.*t115.*3.0-coef_f1_q_sym13.*q0.*q1.*q2.*2.0-coef_f1_q_sym14.*q0.*q1.*q3.*2.0+coef_f0_q_sym9.*q0.*q2.*q3-coef_f1_q_sym16.*q0.*q2.*q3,-coef_f2_q_sym18.*q0+coef_f0_q_sym18.*q2-coef_f2_q_sym2.*q0.*t114+coef_f0_q_sym6.*q0.*t116-coef_f2_q_sym12.*q0.*t115.*3.0-coef_f2_q_sym15.*q0.*t116-coef_f2_q_sym17.*q0.*t117-coef_f2_q_sym5.*q1.*t114.*2.0+coef_f0_q_sym13.*q1.*t116.*2.0+coef_f0_q_sym2.*q2.*t114-coef_f2_q_sym6.*q2.*t114+coef_f0_q_sym12.*q2.*t115.*3.0+coef_f0_q_sym15.*q2.*t116+coef_f0_q_sym17.*q2.*t117-coef_f2_q_sym7.*q3.*t114+coef_f0_q_sym16.*q3.*t116+coef_f0_q_sym5.*q0.*q1.*q2.*2.0-coef_f2_q_sym13.*q0.*q1.*q2.*2.0-coef_f2_q_sym14.*q0.*q1.*q3.*2.0+coef_f0_q_sym7.*q0.*q2.*q3-coef_f2_q_sym16.*q0.*q2.*q3+coef_f0_q_sym14.*q1.*q2.*q3.*2.0,-coef_f3_q_sym18.*q0+coef_f0_q_sym18.*q3-coef_f3_q_sym2.*q0.*t114+coef_f0_q_sym7.*q0.*t117-coef_f3_q_sym12.*q0.*t115.*3.0-coef_f3_q_sym15.*q0.*t116-coef_f3_q_sym17.*q0.*t117-coef_f3_q_sym5.*q1.*t114.*2.0+coef_f0_q_sym14.*q1.*t117.*2.0-coef_f3_q_sym6.*q2.*t114+coef_f0_q_sym16.*q2.*t117+coef_f0_q_sym2.*q3.*t114-coef_f3_q_sym7.*q3.*t114+coef_f0_q_sym12.*q3.*t115.*3.0+coef_f0_q_sym15.*q3.*t116+coef_f0_q_sym17.*q3.*t117-coef_f3_q_sym13.*q0.*q1.*q2.*2.0+coef_f0_q_sym5.*q0.*q1.*q3.*2.0-coef_f3_q_sym14.*q0.*q1.*q3.*2.0+coef_f0_q_sym6.*q0.*q2.*q3-coef_f3_q_sym16.*q0.*q2.*q3+coef_f0_q_sym13.*q1.*q2.*q3.*2.0,q1.*2.0,-coef_f1_q_sym22.*q0+coef_f0_q_sym22.*q1-coef_f1_q_sym3.*q0.*t114+coef_f0_q_sym6.*q0.*t115-coef_f1_q_sym13.*q0.*t115-coef_f1_q_sym21.*q0.*t117-coef_f1_q_sym19.*q0.*t116.*3.0+coef_f0_q_sym3.*q1.*t114-coef_f1_q_sym6.*q1.*t114+coef_f0_q_sym13.*q1.*t115+coef_f0_q_sym21.*q1.*t117+coef_f0_q_sym19.*q1.*t116.*3.0-coef_f1_q_sym8.*q2.*t114.*2.0+coef_f0_q_sym15.*q2.*t115.*2.0-coef_f1_q_sym9.*q3.*t114+coef_f0_q_sym16.*q3.*t115+coef_f0_q_sym8.*q0.*q1.*q2.*2.0-coef_f1_q_sym15.*q0.*q1.*q2.*2.0+coef_f0_q_sym9.*q0.*q1.*q3-coef_f1_q_sym16.*q0.*q1.*q3-coef_f1_q_sym20.*q0.*q2.*q3.*2.0+coef_f0_q_sym20.*q1.*q2.*q3.*2.0,t118+t119+t120+t121+t122+t123+t124+t125+t126+t128+t130+t131+t135-coef_f2_q_sym22.*q0+coef_f0_q_sym22.*q2.*2.0-coef_f2_q_sym3.*q0.*t114+coef_f0_q_sym8.*q0.*t116.*3.0-coef_f2_q_sym13.*q0.*t115-coef_f2_q_sym21.*q0.*t117-coef_f2_q_sym19.*q0.*t116.*3.0-coef_f2_q_sym6.*q1.*t114+coef_f0_q_sym15.*q1.*t116.*3.0+coef_f0_q_sym17.*q1.*t117+coef_f0_q_sym3.*q2.*t114.*2.0-coef_f2_q_sym8.*q2.*t114.*2.0+coef_f0_q_sym13.*q2.*t115.*2.0+coef_f0_q_sym21.*q2.*t117.*2.0+coef_f0_q_sym19.*q2.*t116.*4.0-coef_f2_q_sym9.*q3.*t114+coef_f0_q_sym20.*q3.*t116.*3.0+coef_f0_q_sym14.*q3.*t115-coef_f2_q_sym15.*q0.*q1.*q2.*2.0+coef_f0_q_sym7.*q0.*q1.*q3-coef_f2_q_sym16.*q0.*q1.*q3-coef_f2_q_sym20.*q0.*q2.*q3.*2.0,-coef_f3_q_sym22.*q0+coef_f0_q_sym22.*q3-coef_f3_q_sym3.*q0.*t114+coef_f0_q_sym9.*q0.*t117-coef_f3_q_sym13.*q0.*t115-coef_f3_q_sym21.*q0.*t117-coef_f3_q_sym19.*q0.*t116.*3.0-coef_f3_q_sym6.*q1.*t114+coef_f0_q_sym16.*q1.*t117-coef_f3_q_sym8.*q2.*t114.*2.0+coef_f0_q_sym20.*q2.*t117.*2.0+coef_f0_q_sym3.*q3.*t114-coef_f3_q_sym9.*q3.*t114+coef_f0_q_sym13.*q3.*t115+coef_f0_q_sym21.*q3.*t117+coef_f0_q_sym19.*q3.*t116.*3.0-coef_f3_q_sym15.*q0.*q1.*q2.*2.0+coef_f0_q_sym6.*q0.*q1.*q3-coef_f3_q_sym16.*q0.*q1.*q3+coef_f0_q_sym8.*q0.*q2.*q3.*2.0-coef_f3_q_sym20.*q0.*q2.*q3.*2.0+coef_f0_q_sym15.*q1.*q2.*q3.*2.0,q2.*2.0,-coef_f1_q_sym24.*q0+coef_f0_q_sym24.*q1-coef_f1_q_sym4.*q0.*t114+coef_f0_q_sym7.*q0.*t115-coef_f1_q_sym20.*q0.*t116-coef_f1_q_sym14.*q0.*t115-coef_f1_q_sym23.*q0.*t117.*3.0+coef_f0_q_sym4.*q1.*t114-coef_f1_q_sym7.*q1.*t114+coef_f0_q_sym20.*q1.*t116+coef_f0_q_sym14.*q1.*t115+coef_f0_q_sym23.*q1.*t117.*3.0-coef_f1_q_sym9.*q2.*t114+coef_f0_q_sym16.*q2.*t115-coef_f1_q_sym10.*q3.*t114.*2.0+coef_f0_q_sym17.*q3.*t115.*2.0+coef_f0_q_sym9.*q0.*q1.*q2-coef_f1_q_sym16.*q0.*q1.*q2+coef_f0_q_sym10.*q0.*q1.*q3.*2.0-coef_f1_q_sym17.*q0.*q1.*q3.*2.0-coef_f1_q_sym21.*q0.*q2.*q3.*2.0+coef_f0_q_sym21.*q1.*q2.*q3.*2.0,-coef_f2_q_sym24.*q0+coef_f0_q_sym24.*q2-coef_f2_q_sym4.*q0.*t114+coef_f0_q_sym9.*q0.*t116-coef_f2_q_sym20.*q0.*t116-coef_f2_q_sym14.*q0.*t115-coef_f2_q_sym23.*q0.*t117.*3.0-coef_f2_q_sym7.*q1.*t114+coef_f0_q_sym16.*q1.*t116+coef_f0_q_sym4.*q2.*t114-coef_f2_q_sym9.*q2.*t114+coef_f0_q_sym20.*q2.*t116+coef_f0_q_sym14.*q2.*t115+coef_f0_q_sym23.*q2.*t117.*3.0-coef_f2_q_sym10.*q3.*t114.*2.0+coef_f0_q_sym21.*q3.*t116.*2.0+coef_f0_q_sym7.*q0.*q1.*q2-coef_f2_q_sym16.*q0.*q1.*q2-coef_f2_q_sym17.*q0.*q1.*q3.*2.0+coef_f0_q_sym10.*q0.*q2.*q3.*2.0-coef_f2_q_sym21.*q0.*q2.*q3.*2.0+coef_f0_q_sym17.*q1.*q2.*q3.*2.0,t118+t120+t125+t126+t127+t128+t129+t130+t131+t132+t133+t134+t135-coef_f3_q_sym24.*q0+coef_f0_q_sym24.*q3.*2.0-coef_f3_q_sym4.*q0.*t114+coef_f0_q_sym10.*q0.*t117.*3.0-coef_f3_q_sym20.*q0.*t116-coef_f3_q_sym14.*q0.*t115-coef_f3_q_sym23.*q0.*t117.*3.0-coef_f3_q_sym7.*q1.*t114+coef_f0_q_sym15.*q1.*t116+coef_f0_q_sym17.*q1.*t117.*3.0-coef_f3_q_sym9.*q2.*t114+coef_f0_q_sym13.*q2.*t115+coef_f0_q_sym21.*q2.*t117.*3.0+coef_f0_q_sym4.*q3.*t114.*2.0+coef_f0_q_sym20.*q3.*t116.*2.0-coef_f3_q_sym10.*q3.*t114.*2.0+coef_f0_q_sym14.*q3.*t115.*2.0+coef_f0_q_sym23.*q3.*t117.*4.0+coef_f0_q_sym6.*q0.*q1.*q2-coef_f3_q_sym16.*q0.*q1.*q2-coef_f3_q_sym17.*q0.*q1.*q3.*2.0-coef_f3_q_sym21.*q0.*q2.*q3.*2.0,q3.*2.0],[4, 4]);

author: zarathustr

repo: LibQPEP

save_path: github-repos/MATLAB/zarathustr-LibQPEP


text:
function ft_plot_layout(layout, varargin)

% FT_PLOT_LAYOUT plots a two-dimensional channel layout
%
% Use as
%   ft_plot_layout(layout, ...)
% where the layout is a FieldTrip structure obtained from FT_PREPARE_LAYOUT.
%
% Additional options should be specified in key-value pairs and can be
%   'chanindx'    = list of channels to plot (default is all)
%   'point'       = yes/no
%   'box'         = yes/no
%   'label'       = yes/no
%   'labeloffset' = offset of label from point (default = 0)
%   'labelrotate' = scalar, vector with rotation angle (in degrees) per label (default = 0)
%   'labelalignh' = string, or cell-array specifying the horizontal alignment of the text (default = 'center')
%   'labelalignv' = string, or cell-array specifying the vertical alignment of the text (default = 'middle')
%   'mask'        = yes/no
%   'outline'     = yes/no
%   'verbose'     = yes/no
%   'pointsymbol' = string with symbol (e.g. 'o') - all three point options need to be used together
%   'pointcolor'  = string with color (e.g. 'k')
%   'pointsize'   = number indicating size (e.g. 8)
%   'fontcolor'   = string, color specification (default = 'k')
%   'fontsize'    = number, sets the size of the text (default = 10)
%   'fontunits'   =
%   'fontname'    =
%   'fontweight'  =
%
% It is possible to plot the object in a local pseudo-axis (c.f. subplot), which is specfied as follows
%   'hpos'        = horizontal position of the lower left corner of the local axes
%   'vpos'        = vertical position of the lower left corner of the local axes
%   'width'       = width of the local axes
%   'height'      = height of the local axes
%
% See also FT_PREPARE_LAYOUT, FT_PLOT_TOPO

% Copyright (C) 2009, Robert Oostenveld
%
% 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$

ws = warning('on', 'MATLAB:divideByZero');

% get the optional input arguments
chanindx     = ft_getopt(varargin, 'chanindx',     []);
hpos         = ft_getopt(varargin, 'hpos',         0);
vpos         = ft_getopt(varargin, 'vpos',         0);
width        = ft_getopt(varargin, 'width',        []);
height       = ft_getopt(varargin, 'height',       []);
point        = ft_getopt(varargin, 'point',        true);
box          = ft_getopt(varargin, 'box',          true);
label        = ft_getopt(varargin, 'label',        true);
labeloffset  = ft_getopt(varargin, 'labeloffset',  0);
labelxoffset = ft_getopt(varargin, 'labelxoffset', labeloffset);
labelyoffset = ft_getopt(varargin, 'labelyoffset', labeloffset*1.5);
mask         = ft_getopt(varargin, 'mask',         true);
outline      = ft_getopt(varargin, 'outline',      true);
verbose      = ft_getopt(varargin, 'verbose',      false);
pointsymbol  = ft_getopt(varargin, 'pointsymbol');
pointcolor   = ft_getopt(varargin, 'pointcolor');
pointsize    = ft_getopt(varargin, 'pointsize');

% these have to do with the font
fontcolor   = ft_getopt(varargin, 'fontcolor', 'k'); % default is black
fontsize    = ft_getopt(varargin, 'fontsize',   get(0, 'defaulttextfontsize'));
fontname    = ft_getopt(varargin, 'fontname',   get(0, 'defaulttextfontname'));
fontweight  = ft_getopt(varargin, 'fontweight', get(0, 'defaulttextfontweight'));
fontunits   = ft_getopt(varargin, 'fontunits',  get(0, 'defaulttextfontunits'));
% these have to do with the font
interpreter  = ft_getopt(varargin, 'interpreter', 'tex');

% some stuff related to some refined label plotting
labelrotate   = ft_getopt(varargin, 'labelrotate',  0);
labelalignh   = ft_getopt(varargin, 'labelalignh',  'center');
labelalignv   = ft_getopt(varargin, 'labelalignv',  'middle');
labelcolor    = ft_getopt(varargin, 'labelcolor', 'k');

% convert between true/false/yes/no etc. statements
point   = istrue(point);
box     = istrue(box);
label   = istrue(label);
mask    = istrue(mask);
outline = istrue(outline);
verbose = istrue(verbose);

% color management
if ischar(pointcolor) && exist([pointcolor '.m'], 'file')
  pointcolor = eval(pointcolor);
end

if ~(point || box || label || mask || outline)
  % there is nothing to be plotted
  return;
end

% everything is added to the current figure
holdflag = ishold;
if ~holdflag
  hold on
end

% make a selection of the channels
if ~isempty(chanindx)
  layout.pos    = layout.pos(chanindx,:);
  layout.width  = layout.width(chanindx);
  layout.height = layout.height(chanindx);
  layout.label  = layout.label(chanindx);
end

% the units can be arbitrary (e.g. relative or pixels), so we need to compute the right scaling factor and offset
% create a matrix with all coordinates from positions, mask, and outline
allCoords = layout.pos;
if isfield(layout, 'mask') && ~isempty(layout.mask)
  for k = 1:numel(layout.mask)
    allCoords = [allCoords; layout.mask{k}];
  end
end
if isfield(layout, 'outline') &&~isempty(layout.outline)
  for k = 1:numel(layout.outline)
    allCoords = [allCoords; layout.outline{k}];
  end
end

naturalWidth = (max(allCoords(:,1))-min(allCoords(:,1)));
naturalHeight = (max(allCoords(:,2))-min(allCoords(:,2)));

if isempty(width) && isempty(height)
  xScaling = 1;
  yScaling = 1;
elseif isempty(width) && ~isempty(height)
  % height specified, auto-compute width while maintaining aspect ratio
  yScaling = height/naturalHeight;
  xScaling = yScaling;
elseif ~isempty(width) && isempty(height)
  % width specified, auto-compute height while maintaining aspect ratio
  xScaling = width/naturalWidth;
  yScaling = xScaling;
else
  % both width and height specified
  xScaling = width/naturalWidth;
  yScaling = height/naturalHeight;
end

X      = layout.pos(:,1)*xScaling + hpos;
Y      = layout.pos(:,2)*yScaling + vpos;
Width  = layout.width*xScaling;
Height = layout.height*yScaling;
Lbl    = layout.label;

if point
  if ~isempty(pointsymbol) && ~isempty(pointcolor) && ~isempty(pointsize) % if they're all non-empty, don't use the default
    plot(X, Y, 'marker', pointsymbol, 'color', pointcolor, 'markersize', pointsize, 'linestyle', 'none');
  else
    plot(X, Y, 'marker', '.', 'color', 'b', 'linestyle', 'none');
    plot(X, Y, 'marker', 'o', 'color', 'y', 'linestyle', 'none');
  end
end

if label
  % the MATLAB text function fails if the position for the string is specified in single precision
  X = double(X);
  Y = double(Y);

  % check whether fancy label plotting is needed, this requires a for loop,
  % otherwise print text in a single shot
  if numel(labelrotate)==1
    text(X+labelxoffset, Y+labelyoffset, Lbl , 'interpreter', interpreter, 'horizontalalignment', labelalignh, 'verticalalignment', labelalignv, 'color', fontcolor, 'fontunits', fontunits, 'fontsize', fontsize, 'fontname', fontname, 'fontweight', fontweight);
  else
    n = numel(Lbl);
    if ~iscell(labelalignh)
      labelalignh = repmat({labelalignh},[n 1]);
    end
    if ~iscell(labelalignv)
      labelalignv = repmat({labelalignv},[n 1]);
    end
    if numel(Lbl)~=numel(labelrotate)||numel(Lbl)~=numel(labelalignh)||numel(Lbl)~=numel(labelalignv)
      error('there is something wrong with the input arguments');
    end
    for k = 1:numel(Lbl)
      h = text(X(k)+labelxoffset, Y(k)+labelyoffset, Lbl{k}, 'interpreter', interpreter, 'horizontalalignment', labelalignh{k}, 'verticalalignment', labelalignv{k}, 'rotation', labelrotate(k), 'color', fontcolor, 'fontunits', fontunits, 'fontsize', fontsize, 'fontname', fontname, 'fontweight', fontweight);
    end
  end
end

if box
  line([X-Width/2 X+Width/2 X+Width/2 X-Width/2 X-Width/2]',[Y-Height/2 Y-Height/2 Y+Height/2 Y+Height/2 Y-Height/2]', 'color', [0 0 0]);
end

if outline && isfield(layout, 'outline')
  if verbose
    fprintf('solid lines indicate the outline, e.g. head shape or sulci\n');
  end
  for i=1:length(layout.outline)
    if ~isempty(layout.outline{i})
      X = layout.outline{i}(:,1)*xScaling + hpos;
      Y = layout.outline{i}(:,2)*yScaling + vpos;
      h = line(X, Y);
      set(h, 'color', 'k');
      set(h, 'linewidth', 2);
    end
  end
end

if mask && isfield(layout, 'mask')
  if verbose
    fprintf('dashed lines indicate the mask for topograpic interpolation\n');
  end
  for i=1:length(layout.mask)
    if ~isempty(layout.mask{i})
      X = layout.mask{i}(:,1)*xScaling + hpos;
      Y = layout.mask{i}(:,2)*yScaling + vpos;
      % the polygon representing the mask should be closed
      X(end+1) = X(1);
      Y(end+1) = Y(1);
      h = line(X, Y);
      set(h, 'color', 'k');
      set(h, 'linewidth', 1.5);
      set(h, 'linestyle', ':');
    end
  end
end

axis auto
axis equal
axis off

if ~holdflag
  hold off
end

warning(ws); %revert to original state

author: spm


text:
function [ A,E,Z,iter1 ] = core_RMAMR( D,lambda, I, J, tol, maxiter )
%MC_RPCA_MIXED_NOISE implements the inexact augmented lagrange multiplier
% method for matrix recovery with erase and sparse noise
%   D    - m*n matrix of observations
%
% lambda - weight on sparse error term in the cost function
%
% I, J   - two dimension index of the regions containing values
%
% tol    - tolerance for stopping criterion
%    -DEFAULT 1e-3 if omitted or -1
%
% maxiter - maximum number of iterations
%    -DEFAULT 500 if omitted or -1
%
% Model:
%   min |A|_* +lambda*|ProjectionOnOmega(E)|_1 + gamma*|ProjectionOnOmega(Z)|_F^2
%   subj A+E+Z=D;
% Copyright:Xinchen YE, Tianjin University, 2014

[m,n] = size(D);
if nargin < 5
    tol = 1e-3;
elseif tol == -1
    tol = 1e-3;
end
if nargin < 6
    maxiter = 500;
elseif maxiter == -1
    maxiter = 500;
end

rho = 1.2;%1.1+2.5*rho_s;
% lambda = 10;%1/sqrt(m);
gamma =1; % weight on noise term in the cost function
norm_two = lansvd(D, 1, 'L');   %computes the 1 largest singular value
muk = 10/norm_two; %can be tuned
d_norm=norm(D,'fro');

Ek=zeros(m,n);Yk=zeros(m,n);
Zk=zeros(m,n);
iter1=0;
converged1=false;
sv = 5;
while ~converged1
    iter1 = iter1+1;
    [U, S, V]=lansvd(D-Ek-Zk+(1/muk)*Yk,sv,'L');
    diagS = diag(S);
    diagS = diagS(1:sv);
    svn = length(find(diagS > 1/muk));
    svp = svn;

    ratio = diagS(1:end-1)./diagS(2:end);
    [max_ratio, max_idx] = max(ratio);
    if max_ratio > 2
        svp = min(svn, max_idx);
    end
    if svp < sv %|| iter < 10
        sv = min(svp + 1, n);
    else
        sv = min(svp + 10, n);
    end
   Ak=U(:,1:svp)*diag(diagS(1:svp)-1/muk)*V(:,1:svp)';
%     [U,S,V] = svd(D-Ek-Zk+(1/muk)*Yk);
%     Ak = U*(shrink(S,1/muk))*V';

    Ek = MtOmega(shrink(D-Ak-Zk+(1/muk)*Yk,lambda/muk),I,J,m,n)+...
        D-Ak-Zk+(1/muk)*Yk-MtOmega(D-Ak-Zk+(1/muk)*Yk,I,J,m,n);

    Zk = (muk/(muk+2*gamma))*MtOmega(D-Ak-Ek+(1/muk)*Yk,I,J,m,n)+...
        D-Ak-Ek+(1/muk)*Yk-MtOmega(D-Ak-Ek+(1/muk)*Yk,I,J,m,n);
    Yk=Yk+muk*(D-Ak-Ek-Zk);
    muk=rho*muk;

     stopCriterion = norm(D-Ak-Ek-Zk, 'fro') / d_norm;
      disp([ ' r(F) ' num2str(rank(Ak))...
            ' |E|_0 ' num2str(length(find(abs(Ek)>0)))...
            ' |Z|_0 ' num2str(length(find(abs(Zk)>0)))...
            ' stopCriterion ' num2str(stopCriterion)  ' iter1 ' num2str(iter1) ' mu ' num2str(muk)]);
    if stopCriterion < tol
        converged1 = true;
    end
    if ~converged1&&iter1>=maxiter
        disp('Maximum iterations reached');
        converged1=true;
    end
end

A=Ak;
E=Ek;
Z=Zk;
end

author: andrewssobral


text:
% (this is just an interface)
function varargout = mydatejd_inv (varargin)
    [varargout{1:nargout}] = mydatejdi (varargin{:});
end

%!test
%! % just check whether it runs:
%! mydatejd_inv(1);


text:
function [u, info] = solvebvpLinear(L, rhs, Ninit, pref, displayInfo)
%SOLVEBVP  Solve a linear CHEBOP BVP system.
%
% [U, INFO] = SOLVEBVPLINEAR(N, L, RHS, PREF, DISPLAYINFO), where:
%   N is a CHEBOP
%   L is a linear CHEBOP
%   RHS is a CHEBMATRIX
%   PREF is a CHEBOPPREF
%   DISPLAYINFO is a function handle
%
% attempts to solve the linear BVP
%
%       L*U = RHS + boundary conditions specified by L
%
% The output U is a CHEBMATRIX, and INFO is a MATLAB struct with useful
% information. This method should generally not called directly by the user,
% which should rather call the CHEBOP/MLDIVIDE or CHEBOP/SOLVEBVP methods.
%
% Observe that when this method is called, any affine parts of the original
% CHEBOP have been absorbed into RHS.
%
% See also: chebop/solvebvp

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

% Get defaults:
if ( nargin < 3 )
    pref = cheboppref();
end

% All values of the LINOPCONSTRAINT stored in L will be of incorrect sign when
% returned from LINEARIZE(), if we want to use it for a LINOP backslash. This is
% because when problems are solved with LINOP backslash, the solution to the
% problem is the output itself, while in a Newton iteration, we have to add the
% output of the LINOP solution to the current guess. Thus, flip the signs of the
% values of L.constraint:
L.constraint = -L.constraint;

% Solve the linear problem:
del = linsolve(L, rhs, pref);

if ( ~isempty(Ninit) )
    % If Ninit is not empty, N will have been linearized around Ninit. In that
    % case, we need to regard the solution del obtained above as a Newton
    % correction to Ninit.
    u = Ninit + del;
    u = simplify(u);
else
    u = del;
end

% Norm of residual:
normRes = norm(L*u - rhs, 'fro');

if ( nargin > 3 )
    % Print information after linear problem has been solved:
    displayInfo('linear', u, normRes, pref)
end

if ( nargout > 1 )
    % Return the norm of the residual in the INFO struct:
    info.error = normRes;
end

end

author: chebfun

repo: chebfun


text:
%
% DEPRECATED
%

h = findobj('Tag','EEGLAB');
i = get(h,'Children');
k = get(i,'Tag');
tags = regexpi(k,'\w*win\w*','match');
tags = [tags{:}];
for a=1:length(tags)
      hh = findobj('Tag',tags{a});
      set(hh,'backgroundcolor','g')
end
hm = findobj('Tag', 'Frame1');
set(hm,'backgroundcolor','g')
hm = findobj('Tag', 'EEGLAB');
set(hm,'Color','g')

author: ucdavis

repo: erplab

save_path: github-repos/MATLAB/ucdavis-erplab


text:
function x = p20_start ( n )

%*****************************************************************************80
%
%% P20_START returns a starting point for optimization for problem 20.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    19 March 2000
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, integer N, the number of variables X.
%
%    Output, real X(N), a starting point for the optimization.
%
  x = zeros ( n, 1 );

  return
end

author: johannesgerer

repo: jburkardt-m