text:
classdef BackgroundSubtractorMOG2 < handle
    %BACKGROUNDSUBTRACTORMOG2  Gaussian Mixture-based Background/Foreground Segmentation Algorithm
    %
    % The class implements the Gaussian mixture model background subtraction
    % described in [Zivkovic2004] and [Zivkovic2006].
    %
    % The code is very fast and performs also shadow detection. Number of
    % Gausssian components is adapted per pixel.
    %
    % The algorithm similar to the standard Stauffer&Grimson algorithm with
    % additional selection of the number of the Gaussian components based on
    % [Zivkovic04recursiveunsupervised].
    %
    % ## References
    % [Zivkovic2004]:
    % > Zoran Zivkovic. "Improved adaptive gaussian mixture model for
    % > background subtraction". In Pattern Recognition, 2004. ICPR 2004.
    % > Proceedings of the 17th International Conference on, volume 2,
    % > pages 28-31. IEEE, 2004.
    % > [PDF](http://www.zoranz.net/Publications/zivkovic2004ICPR.pdf).
    %
    % [Zivkovic2006]:
    % > Zoran Zivkovic and Ferdinand van der Heijden. "Efficient adaptive
    % > density estimation per image pixel for the task of background
    % > subtraction". Pattern recognition letters, 27(7):773-780, 2006.
    % > [PDF](http://www.zoranz.net/Publications/zivkovicPRL2006.pdf).
    %
    % [Zivkovic04recursiveunsupervised]:
    % > Zoran Zivkovic and Ferdinand van der Heijden, "Recursive unsupervised
    % > learning of finite mixture models", IEEE Trans. on Pattern Analysis
    % > and Machine Intelligence, vol.26, no.5, pages 651-656, 2004.
    % > [PDF](http://www.zoranz.net/Publications/zivkovic2004PAMI.pdf).
    %
    % [Prati03detectingmoving]:
    % > Andrea Prati, Ivana Mikic, Mohan M. Trivedi, Rita Cucchiara.
    % > "Detecting Moving Shadows: Algorithms and Evaluation", IEEE PAMI, 2003.
    %
    % See also: cv.BackgroundSubtractorMOG2.BackgroundSubtractorMOG2,
    %  cv.BackgroundSubtractorMOG2.apply,
    %  cv.BackgroundSubtractorMOG2.getBackgroundImage,
    %  vision.ForegroundDetector
    %

    properties (SetAccess = private)
        % Object ID
        id
    end

    properties (Dependent)
        % The number of last frames that affect the background model.
        History
        % The number of gaussian components in the background model.
        % The model needs to be reinitalized to reserve memory.
        NMixtures
        % The "background ratio" parameter of the algorithm.
        % If a foreground pixel keeps semi-constant value for about
        % `BackgroundRatio * History` frames, it's considered background and
        % added to the model as a center of a new component. It corresponds to
        % `TB` parameter in the paper.
        BackgroundRatio
        % The variance threshold for the pixel-model match.
        % The main threshold on the squared Mahalanobis distance to decide if
        % the sample is well described by the background model or not. Related
        % to `Cthr` from the paper.
        VarThreshold
        % The variance threshold for the pixel-model match used for new
        % mixture component generation.
        % Threshold for the squared Mahalanobis distance that helps decide
        % when a sample is close to the existing components (corresponds to
        % `Tg` in the paper). If a pixel is not close to any component, it is
        % considered foreground or added as a new component.
        % `3 sigma => Tg=3*3=9` is default. A smaller `Tg` value generates
        % more components. A higher `Tg` value may result in a small number of
        % components but they can grow too large.
        VarThresholdGen
        % The initial variance of gaussian components. default 15
        VarInit
        % Minimmum variance. default 4
        VarMin
        % Maximum variance. default 5*15
        VarMax
        % The complexity reduction threshold.
        % This parameter defines the number of samples needed to accept to
        % prove the component exists. `CT=0.05` is a default value for all the
        % samples. By setting `CT=0` you get an algorithm very similar to the
        % standard Stauffer&Grimson algorithm.
        ComplexityReductionThreshold
        % The shadow detection flag.
        % If true, the algorithm detects shadows and marks them. See
        % cv.BackgroundSubtractorMOG2.BackgroundSubtractorMOG2 for details.
        DetectShadows
        % The shadow value.
        % Shadow value is the value used to mark shadows in the foreground
        % mask. Default value is 127. Value 0 in the mask always means
        % background, 255 means foreground.
        ShadowValue
        % The shadow threshold.
        % A shadow is detected if pixel is a darker version of the background.
        % The shadow threshold (Tau in the paper) is a threshold defining how
        % much darker the shadow can be. Tau=0.5 means that if a pixel is more
        % than twice darker then it is not shadow.
        % See [Prati03detectingmoving].
        ShadowThreshold
    end

    %% BackgroundSubtractor
    methods
        function this = BackgroundSubtractorMOG2(varargin)
            %BACKGROUNDSUBTRACTORMOG2  Creates MOG2 Background Subtractor
            %
            %     bs = cv.BackgroundSubtractorMOG2()
            %     bs = cv.BackgroundSubtractorMOG2('OptionName', optionValue, ...)
            %
            % ## Options
            % * __History__ Length of the history. default 500
            % * __VarThreshold__ Threshold on the squared Mahalanobis distance
            %   between the pixel and the model to decide whether a pixel is
            %   well described by the background model. This parameter does
            %   not affect the background update. default 16
            % * __DetectShadows__ If true, the algorithm will detect shadows
            %   and mark them. It decreases the speed a bit, so if you do not
            %   need this feature, set the parameter to false. default true
            %
            % See also: cv.BackgroundSubtractorMOG2
            %
            this.id = BackgroundSubtractorMOG2_(0, 'new', varargin{:});
        end

        function delete(this)
            %DELETE  Destructor
            %
            %     bs.delete()
            %
            % See also: cv.BackgroundSubtractorMOG2
            %
            if isempty(this.id), return; end
            BackgroundSubtractorMOG2_(this.id, 'delete');
        end

        function fgmask = apply(this, im, varargin)
            %APPLY  Updates the background model and computes the foreground mask
            %
            %     fgmask = bs.apply(im)
            %     fgmask = bs.apply(im, 'OptionName', optionValue, ...)
            %
            % ## Input
            % * __im__ Next video frame, uint8 or single. Floating point frame
            %   will be used without scaling and should be in range [0,255].
            %
            % ## Output
            % * __fgmask__ The output foreground mask as an 8-bit binary image
            %   (0 for background, 255 for foregound, and `ShadowValue` for
            %   shadows if `DetectShadows` is true).
            %
            % ## Options
            % * __LearningRate__ The value between 0 and 1 that indicates how
            %   fast the background model is learnt. Negative parameter value
            %   makes the algorithm to use some automatically chosen learning
            %   rate. 0 means that the background model is not updated at all,
            %   1 means that the background model is completely reinitialized
            %   from the last frame. default -1
            %
            % See also: cv.BackgroundSubtractorMOG2.getBackgroundImage
            %
            fgmask = BackgroundSubtractorMOG2_(this.id, 'apply', im, varargin{:});
        end

        function bgImg = getBackgroundImage(this)
            %GETBACKGROUNDIMAGE  Computes a background image
            %
            %     bgImg = bs.getBackgroundImage()
            %
            % ## Output
            % * __bgImg__ The output background image, which is the mean of
            %   all background Gaussians.
            %
            % ### Note
            % Sometimes the background image can be very blurry, as it contain
            % the average background statistics.
            %
            % See also: cv.BackgroundSubtractorMOG2.apply
            %
            bgImg = BackgroundSubtractorMOG2_(this.id, 'getBackgroundImage');
        end
    end

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

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

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

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

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

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

        function value = get.NMixtures(this)
            value = BackgroundSubtractorMOG2_(this.id, 'get', 'NMixtures');
        end
        function set.NMixtures(this, value)
            BackgroundSubtractorMOG2_(this.id, 'set', 'NMixtures', value);
        end

        function value = get.BackgroundRatio(this)
            value = BackgroundSubtractorMOG2_(this.id, 'get', 'BackgroundRatio');
        end
        function set.BackgroundRatio(this, value)
            BackgroundSubtractorMOG2_(this.id, 'set', 'BackgroundRatio', value);
        end

        function value = get.VarThreshold(this)
            value = BackgroundSubtractorMOG2_(this.id, 'get', 'VarThreshold');
        end
        function set.VarThreshold(this, value)
            BackgroundSubtractorMOG2_(this.id, 'set', 'VarThreshold', value);
        end

        function value = get.VarThresholdGen(this)
            value = BackgroundSubtractorMOG2_(this.id, 'get', 'VarThresholdGen');
        end
        function set.VarThresholdGen(this, value)
            BackgroundSubtractorMOG2_(this.id, 'set', 'VarThresholdGen', value);
        end

        function value = get.VarInit(this)
            value = BackgroundSubtractorMOG2_(this.id, 'get', 'VarInit');
        end
        function set.VarInit(this, value)
            BackgroundSubtractorMOG2_(this.id, 'set', 'VarInit', value);
        end

        function value = get.VarMin(this)
            value = BackgroundSubtractorMOG2_(this.id, 'get', 'VarMin');
        end
        function set.VarMin(this, value)
            BackgroundSubtractorMOG2_(this.id, 'set', 'VarMin', value);
        end

        function value = get.VarMax(this)
            value = BackgroundSubtractorMOG2_(this.id, 'get', 'VarMax');
        end
        function set.VarMax(this, value)
            BackgroundSubtractorMOG2_(this.id, 'set', 'VarMax', value);
        end

        function value = get.ComplexityReductionThreshold(this)
            value = BackgroundSubtractorMOG2_(this.id, 'get', 'ComplexityReductionThreshold');
        end
        function set.ComplexityReductionThreshold(this, value)
            BackgroundSubtractorMOG2_(this.id, 'set', 'ComplexityReductionThreshold', value);
        end

        function value = get.DetectShadows(this)
            value = BackgroundSubtractorMOG2_(this.id, 'get', 'DetectShadows');
        end
        function set.DetectShadows(this, value)
            BackgroundSubtractorMOG2_(this.id, 'set', 'DetectShadows', value);
        end

        function value = get.ShadowValue(this)
            value = BackgroundSubtractorMOG2_(this.id, 'get', 'ShadowValue');
        end
        function set.ShadowValue(this, value)
            BackgroundSubtractorMOG2_(this.id, 'set', 'ShadowValue', value);
        end

        function value = get.ShadowThreshold(this)
            value = BackgroundSubtractorMOG2_(this.id, 'get', 'ShadowThreshold');
        end
        function set.ShadowThreshold(this, value)
            BackgroundSubtractorMOG2_(this.id, 'set', 'ShadowThreshold', value);
        end
    end

end

author: kyamagu

repo: mexopencv

save_path: github-repos/MATLAB/kyamagu-mexopencv


text:
function [new_img,new_seg,rimgs,rsegs] = Create_ROIs(img, seg,patchSize)
% Create the ROIs needed, x*y*z
% Also randomly crop patches to increase background information
% -------------------------------------------------------------
    % the agumented size of surrounding context
    aug_r =20;

    img = single(img);

    % get bounding box
    [s1, s2, s3] = size(seg);
    positions = regionprops(seg>0, 'BoundingBox');
    position  = positions(1).BoundingBox;
    cx = max(ceil(position(1)) - aug_r/2, 1);
    cy = max(ceil(position(2)) - aug_r/2, 1);
    cz = max(ceil(position(3)) - aug_r/2, 1);
    x  = min(ceil(position(4)) + aug_r, s2 - cx + 1);
    y  = min(ceil(position(5)) + aug_r, s1 - cy + 1);
    z  = min(ceil(position(6)) + aug_r, s3 - cz + 1);

    % crop ROI
    new_img = img(cy:cy+y-1, cx:cx+x-1, cz:cz+z-1, :);
    new_seg = seg(cy:cy+y-1, cx:cx+x-1, cz:cz+z-1);

    % padding the ROI to patchSize*patchSize*patchSize
    pading_size = max(0, ceil((patchSize-size(img))/2));
    new_img = padarray(new_img,pading_size);
    new_seg = padarray(new_seg,pading_size,100);

    % randomly crop patches from the whole volume
    rimgs = [];
    rsegs  = [];
    for k = 0: 0
        [sy, sx, sz] = size(img);
        cy = floor(rand()*(sy -y)+1);
        cx = floor(rand()*(sx -x)+1);
        cz = max(floor(rand()*(sz - z+1)), 1);
        random_per_imgs = img(cy:cy+y-1, cx:cx+x-1, cz:cz+z-1,:);
        random_per_seg1 = seg(cy:cy+y-1, cx:cx+x-1, cz:cz+z-1);
        rimgs{end+1}  = random_per_imgs;
        rsegs{end+1} = random_per_seg1;
    end
end

author: yulequan

repo: HeartSeg

save_path: github-repos/MATLAB/yulequan-HeartSeg


text:
function Tbig = multiply_by_pot(Tbig, Tsmall)
% MULTIPLY_BY_POT Tbig *= Tsmall
% Tbig = multiply_by_pot(Tbig, Tsmall)
%
% Tsmall's domain must be a subset of Tbig's domain.

smallp = extend_domain_table(Tsmall.p, Tsmall.domain, Tsmall.sizes, Tbig.domain, Tbig.sizes);
Tbig.p = Tbig.p .* smallp;

smallu = extend_domain_table(Tsmall.u, Tsmall.domain, Tsmall.sizes, Tbig.domain, Tbig.sizes);
Tbig.u = Tbig.u + smallu;

author: bayesnet

repo: bnt

save_path: github-repos/MATLAB/bayesnet-bnt


text:
function L=ishappy(n)

N=n;
L=0;
NHistory=n;
while L==0
    N=sum(str2num(num2str(N).').^2); %Sum digits
    if N==1 %n is a happy number
        L=1;
    elseif ismember(N,NHistory); %Kill we are in a loop
        break
    end
    NHistory=[NHistory N];
end
NHistory=[NHistory N];

%%
% _*GIBBON footer text*_
%
% License: <https://github.com/gibbonCode/GIBBON/blob/master/LICENSE>
%
% GIBBON: The Geometry and Image-based Bioengineering add-On. A toolbox for
% image segmentation, image-based modeling, meshing, and finite element
% analysis.
%
% Copyright (C) 2006-2022 Kevin Mattheus Moerman and the GIBBON contributors
%
% This program is free software: you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation, either version 3 of the License, or
% (at your option) any later version.
%
% This program is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
% GNU General Public License for more details.
%
% You should have received a copy of the GNU General Public License
% along with this program.  If not, see <http://www.gnu.org/licenses/>.

author: gibbonCode

repo: GIBBON

save_path: github-repos/MATLAB/gibbonCode-GIBBON


text:
function [R, t] = read4PCSResults(filename)

 M = dlmread(filename);
 R = M(1:3, 1:3);
 t = M(1:3, 4);

end

author: intellhave

repo: SDRSAC

save_path: github-repos/MATLAB/intellhave-SDRSAC


text:
%Saving data into ASPAR type 5 format

report_this_filefun(mfilename('fullpath'));

str = [];
[newmatfile] = ['t1.sum'];
[ newpath] = [hodi '/aspar/'];

newt3 = [];
do = [ ' cd ' newpath ];
eval(do)

% lets add the mainshock as the fisrt and largest event...
l = newt2.Date > maepi(1,3);
newt3 =  newt2(l,1:9);
newt3 = [maepi(1,1:9) ; newt3 ];
lam = (newt3(:,2)-floor(newt3(:,2)))*100*6/10;
lom = (newt3(:,1)-floor(newt3(:,1)))*100*6/10;

s = [floor(newt3(:,3:5))  newt3(:,8:9) floor(newt3(:,2)) lam floor(abs((newt3(:,1))))  lom  newt3(:,7) newt3(:,6)];

fid = fopen([newpath newmatfile],'w') ;;
fprintf(fid,'%5.0f %5.0f %5.0f %5.0f %5.0f %5.0f %7.3f %5.0f %7.3f %7.3f %7.3f\n',s');
fclose(fid);
clear s
return

author: CelsoReyes


text:
function y = quickselect(x, k)
% Brings the largest k entries to the front of array
    y = mex_quickselect(x, k);
end

author: indigits

repo: sparse-plex


text:
function [ n_data, n, x, fx ] = w_polynomial_values ( n_data )

%*****************************************************************************80
%
%% W_POLYNOMIAL_VALUES returns values of Chebyshev polynomials W(n,x).
%
%  Discussion:
%
%    In Mathematica, the function can be evaluated by:
%
%      u = Sqrt[(x+1)/2],
%      ChebyshevU[2*n,u]
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    23 April 2012
%
%  Author:
%
%    John Burkardt
%
%  Reference:
%
%    Milton Abramowitz and Irene Stegun,
%    Handbook of Mathematical Functions,
%    US Department of Commerce, 1964.
%
%    Stephen Wolfram,
%    The Mathematica Book,
%    Fourth Edition,
%    Wolfram Media / Cambridge University Press, 1999.
%
%  Parameters:
%
%    Input/output, integer N_DATA.  The user sets N_DATA to 0 before the
%    first call.  On each call, the routine increments N_DATA by 1, and
%    returns the corresponding data; when there is no more data, the
%    output value of N_DATA will be 0 again.
%
%    Output, integer N, the order of the function.
%
%    Output, real X, the point where the function is evaluated.
%
%    Output, real FX, the value of the function.
%
  n_max = 13;

  fx_vec = [ ...
     1.000000000000000E+00, ...
     2.600000000000000E+00, ...
     3.160000000000000E+00, ...
     2.456000000000000E+00, ...
     0.769600000000000E+00, ...
    -1.224640000000000E+00, ...
    -2.729024000000000E+00, ...
    -3.141798400000000E+00, ...
    -2.297853440000000E+00, ...
    -0.534767104000000E+00, ...
     1.442226073600000E+00, ...
     2.842328821760000E+00, ...
     3.105500041216000E+00 ];

  n_vec = [ ...
     0,  1,  2, ...
     3,  4,  5, ...
     6,  7,  8, ...
     9, 10, 11, ...
    12 ];

  x_vec = [ ...
     0.8E+00, ...
     0.8E+00, ...
     0.8E+00, ...
     0.8E+00, ...
     0.8E+00, ...
     0.8E+00, ...
     0.8E+00, ...
     0.8E+00, ...
     0.8E+00, ...
     0.8E+00, ...
     0.8E+00, ...
     0.8E+00, ...
     0.8E+00 ];

  if ( n_data < 0 )
    n_data = 0;
  end

  n_data = n_data + 1;

  if ( n_max < n_data )
    n_data = 0;
    n = 0;
    x = 0.0;
    fx = 0.0;
  else
    n = n_vec(n_data);
    x = x_vec(n_data);
    fx = fx_vec(n_data);
  end

  return
end

author: johannesgerer

repo: jburkardt-m


text:
function update_visualization_func = show_video(img_files, video_path, resize_image)
%SHOW_VIDEO
%   Visualizes a tracker in an interactive figure, given a cell array of
%   image file names, their path, and whether to resize the images to
%   half size or not.
%
%   This function returns an UPDATE_VISUALIZATION function handle, that
%   can be called with a frame number and a bounding box [x, y, width,
%   height], as soon as the results for a new frame have been calculated.
%   This way, your results are shown in real-time, but they are also
%   remembered so you can navigate and inspect the video afterwards.
%   Press 'Esc' to send a stop signal (returned by UPDATE_VISUALIZATION).
%
%   Joao F. Henriques, 2014
%   http://www.isr.uc.pt/~henriques/
%
%   revised by: Yang Li, August, 2014
%   http://ihpdep.github.io

	%store one instance per frame
	num_frames = numel(img_files);
	boxes = cell(num_frames,1);

	%create window
	[fig_h, axes_h, unused, scroll] = videofig(num_frames, @redraw, [], [], @on_key_press);  %#ok, unused outputs
	set(fig_h, 'Number','off', 'Name', ['Tracker - ' video_path])
	axis off;

	%image and rectangle handles start empty, they are initialized later
	im_h = [];
	rect_h = [];

	update_visualization_func = @update_visualization;
	stop_tracker = false;

	function stop = update_visualization(frame, box)
		%store the tracker instance for one frame, and show it. returns
		%true if processing should stop (user pressed 'Esc').
		boxes{frame} = box;
		scroll(frame);
		stop = stop_tracker;
	end

	function redraw(frame)
		%render main image
		im = imread([video_path img_files{frame}]);
% 		if size(im,3) > 1,
% 			im = rgb2gray(im);
% 		end
		if resize_image,
			im = imresize(im, 0.5);
		end

		if isempty(im_h),  %create image
			im_h = imshow(im, 'Border','tight', 'InitialMag',200, 'Parent',axes_h);
		else  %just update it
			set(im_h, 'CData', im)
		end

		%render target bounding box for this frame
		if isempty(rect_h),  %create it for the first time
			rect_h = rectangle('Position',[0,0,1,1], 'EdgeColor','g', 'Parent',axes_h);
		end
		if ~isempty(boxes{frame}),
			set(rect_h, 'Visible', 'on', 'Position', boxes{frame});
		else
			set(rect_h, 'Visible', 'off');
		end
	end

	function on_key_press(key)
		if strcmp(key, 'escape'),  %stop on 'Esc'
			stop_tracker = true;
		end
	end

end

author: ihpdep

repo: samf

save_path: github-repos/MATLAB/ihpdep-samf


text:
function [out] = percolation_3(S,Smax)
%percolation_3

% Copyright (C) 2019, 2021 Wouter J.M. Knoben, Luca Trotter
% This file is part of the Modular Assessment of Rainfall-Runoff Models
% Toolbox (MARRMoT).
% MARRMoT is a free software (GNU GPL v3) and distributed WITHOUT ANY
% WARRANTY. See <https://www.gnu.org/licenses/> for details.

% Flux function
% ------------------
% Description:  Non-linear percolation (empirical)
% Constraints:  -
% @(Inputs):    S    - current storage [mm]
%               Smax - maximum contributing storage [mm]

out = Smax^(-4)./(4)*(4/9)^(4)*S^5;

end

author: wknoben

repo: MARRMoT

save_path: github-repos/MATLAB/wknoben-MARRMoT


text:
function [node,PSLG,part] = bfsgeo2(node,PSLG,seed)
%BFSGEO2 partition geometry about "seeds" via breadth-first
%search.
%   [NODE,EDGE,PART] = BFSGEO2(NODE,EDGE,SEED) returns a set
%   of 2-manifold geometry partitions by expanding about a
%   list of "seed" points. Partitions expand in a breadth-
%   first sense until a polygon bounday is encountered. SEED
%   is an S-by-2 array of XY-coordinates to expand around,
%   NODE is an N-by-2 array of polygon vertices, and EDGE is
%   an E-by-2 array of polygon edge indexing. Each row in
%   EDGE represents an edge of the polygon, such that
%   NODE(EDGE(JJ,1),:) and NODE(EDGE(JJ,2),:) are the XY-co-
%   ordinates of the endpoints of the JJ-TH edge. PART is an
%   S-by-1 cell array of geometry partitions, where each
%   PART{KK} is a list of edge indices into EDGE that define
%   the KK-TH partition.
%
%   This function may be useful when seeking to partition
%   complex, non-manifold geometry into a format that's app-
%   ropriate for the triangulation algorithms in REFINE2.
%
%   See also REFINE2, FIXGEO2, BFSTRI2

%-----------------------------------------------------------
%   Darren Engwirda : 2017 --
%   Email           : de2363@columbia.edu
%   Last updated    : 11/10/2017
%-----------------------------------------------------------

%---------------------------------------------- basic checks
    if ( ~isnumeric(node) || ...
         ~isnumeric(PSLG) || ...
         ~isnumeric(seed) )
        error('bfsgeo2:incorrectInputClass' , ...
            'Incorrect input class.') ;
    end

%---------------------------------------------- basic checks
    if (ndims(node) ~= +2 || ...
        ndims(PSLG) ~= +2 || ...
        ndims(seed) ~= +2 )
        error('bfsgeo2:incorrectDimensions' , ...
            'Incorrect input dimensions.');
    end

    if (size(node,2)~= +2 || ...
        size(PSLG,2)~= +2 || ...
        size(seed,2)~= +2 )
        error('bfsgeo2:incorrectDimensions' , ...
            'Incorrect input dimensions.');
    end

    nnod = size(node,1) ;
    nedg = size(PSLG,1) ;

%---------------------------------------------- basic checks
    if (min([PSLG(:)])<+1 || ...
        max([PSLG(:)])>nnod)
        error('bfsgeo2:invalidInputs', ...
            'Invalid EDGE input array.') ;
    end

%------------------------------------------ assemble full CDT
   [node,PSLG,tria] = deltri2 (node,PSLG) ;

%------------------------------------------ find seeds in CDT
   [sptr,stri] = findtria(node,tria,seed) ;

    okay = sptr(:,2) ...
        >= sptr(:,1) ;
    itri = stri(sptr(okay,1));

%------------------------------------------ PART for all seed
    part = {} ;

    for ipos = +1 : size (itri,1)

    %-------------- BFS about current tria.
       [mark] = ...
            bfstri2(PSLG,tria,itri(ipos)) ;

    %-------------- match tria./poly. edges
        edge = [
            tria(mark,[1,2]) ;
            tria(mark,[2,3]) ;
            tria(mark,[3,1]) ;
               ] ;

        edge = sort(edge,+2) ;
        PSLG = sort(PSLG,+2) ;

       [same,epos] = ...
            setset2(edge(:,1:2),PSLG) ;

    %-------------- find match multiplicity
        epos = epos(epos>+0) ;
        epos = sort(epos);

        eidx = ...
            find(diff(epos)) ;

        eptr = ...
            [1;eidx+1;length(epos)+1] ;
        enum = ...
            eptr(2:end)-eptr(1:end-1) ;

    %---------- select singly-matched edges
        part{ipos} = ...
            epos(eptr(enum == 1)) ;

    end

end

author: dengwirda


text:
function r8mat_print ( m, n, a, title )

%*****************************************************************************80
%
%% R8MAT_PRINT prints an R8MAT.
%
%  Discussion:
%
%    An R8MAT is an array of R8's.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    06 September 2004
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, integer M, the number of rows in A.
%
%    Input, integer N, the number of columns in A.
%
%    Input, real A(M,N), the matrix.
%
%    Input, string TITLE, a title to be printed.
%
  r8mat_print_some ( m, n, a, 1, 1, m, n, title );

  return
end


text:
% Codes for CVPR-15 work `Face Alignment by Coarse-to-Fine Shape Searching'
% Any question please contact Shizhan Zhu: zhshzhutah2@gmail.com
% Released on July 25, 2015

function [resultPose,weight] = poseVoting(resultPose_samples,voteIterTot)
% resultPose_samples : M * dim * samplings
% voteIterTot : scalar
%
% resultPose : M * dim
% weight : M * samplings

if ~exist('voteIterTot','var')
    voteIterTot = 100;
end;
M = size(resultPose_samples,1);
dim = size(resultPose_samples,2);
samplings = size(resultPose_samples,3);
% assert(dim <= 3);
resultPose = zeros(M,dim);
weight = zeros(M,samplings);
parfor k = 1:M
    candidate = reshape(resultPose_samples(k,:,:),dim,samplings);
    d = dist(candidate);
    d = d / mean(d(:));
    ed = exp(-d.^2);
    for i = 1:size(ed,1), ed(i,i) = 0;end;
    w = zeros(samplings,voteIterTot);
    w(:,1) = 1 / samplings;
    for i = 1:voteIterTot-1
        for j = 1:samplings, w(:,i+1) = w(:,i+1) + ed(:,j).*w(j,i);end;
        w(:,i+1) = w(:,i+1) / sum(w(:,i+1));
    end;
    weight(k,:) = w(:,voteIterTot)';
    resultPose(k,:) = w(:,voteIterTot)' * candidate';
end;

end

author: zhusz


text:
classdef (Abstract) MimModelCollection  < MimModel

    methods
        function addItem(obj, itemId)
            [currentCollection, ~] = obj.getOrRun();
            if ~ismember(itemId, currentCollection)
                currentCollection{end + 1} = itemId;
                obj.setValue(currentCollection);
            end
        end

        function removeItem(obj, itemId)
            [currentCollection, ~] = obj.getOrRun();
            if ismember(itemId, currentCollection)
                currentCollection = setdiff(currentCollection, itemId);
                obj.setValue(currentCollection);
            end
        end
    end
end

author: tomdoel

repo: pulmonarytoolkit


text:
function n = ndims(v,varargin)
% overloads ndims

n = ndims(v.x,varargin{:});

author: mtex-toolbox

repo: mtex

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


text:
function cobj2 = subset(cobj, varargin)
    %CATALOG.SUBSET Create a new catalogObject by subsetting based
    %on indices.
    % Usage:
    %   cobj.subset('indices', [])
    % or:
    %   cobj.subset('start_time', [], 'end_time', [])
    % can also just give one of start_time or end_time

    %if nargin>0
        % Parse required, optional and param-value pair arguments,
        % set default values, and add validation conditions
        p = inputParser;
        p.addParameter('start_time', [], @isnumeric) % positional
        p.addParameter('end_time', [], @isnumeric)
        p.addParameter('indices', [], @isnumeric)
        p.parse(varargin{:});
        fields = fieldnames(p.Results);
        for i=1:length(fields)
            field=fields{i};
            val = p.Results.(field);
            eval(sprintf('%s = val;',field));
        end
    %end

%     datestr(min(cobj.otime))
%     datestr(max(cobj.otime))
%     datestr(start_time)
%     datestr(end_time)

    if isempty(indices)
        i1 = 1:numel(cobj.otime);
        if start_time
            i1 = find(cobj.otime >= start_time);
            debug.print_debug(1,sprintf('Found %d events after %s',numel(i1), datestr(start_time)));
        end
        i2 = 1:numel(cobj.otime);
        if end_time
            i2 = find(cobj.otime <= end_time);
            debug.print_debug(1,sprintf('Found %d events before %s',numel(i1), datestr(end_time)));
        end
        indices = intersect(i1, i2);
    end

    cobj2 = cobj;
    N = numel(cobj.otime);
    debug.print_debug(1,sprintf('Subsetting from %d events to %d events',N, numel(indices)));
    cobj2.otime = cobj.otime(indices);
    if numel(cobj.lon)==N
        cobj2.lon = cobj.lon(indices);
    end
    if numel(cobj.lat)==N
        cobj2.lat = cobj.lat(indices);
    end
    if numel(cobj.depth)==N
        cobj2.depth = cobj.depth(indices);
    end
    if numel(cobj.mag)==N
        cobj2.mag = cobj.mag(indices);
    end
    if numel(cobj.magtype)==N
        dummy=[];
%         try % there is a limit on subsetting cell array
%         cobj2.magtype = cobj.magtype{indices};
%         end
        for c=1:numel(indices)
            dummy{c} = cobj.magtype{indices(c)};
        end
        cobj2.magtype = dummy;
        clear dummy
    end
    if numel(cobj.etype)==N
        dummy=[];
        for c=1:numel(indices)
            dummy{c} = cobj.etype{indices(c)};
        end
        cobj2.etype = dummy;
        clear dummy;
    end
    if numel(cobj.arrivals)==N
        cobj2.arrivals = cobj.arrivals(indices);
    end
    if numel(cobj.waveforms)==N
        cobj2.waveforms = cobj.waveforms(indices);
%         counter = 1;
%         indices
%         for thisindex=indices
%             cobj2.waveforms{counter} = cobj.waveforms{thisindex};
%             counter = counter + 1;
%         end
    end
    if numel(cobj.ontime)==N
        cobj2.ontime = cobj.ontime(indices);
    end
    if numel(cobj.offtime)==N
        cobj2.offtime = cobj.offtime(indices);
    end

    %% Added by Glenn Thompson 2018-05-01 to add a request start & end time
    if isempty(start_time)
        cobj2.request.start_time = min(cobj2.otime);
    else
        cobj2.request.start_time = start_time;
    end
    if isempty(end_time)
        cobj2.request.end_time = max(cobj2.otime);
    else
        cobj2.request.end_time = end_time;
    end

end

author: geoscience-community-codes

repo: GISMO


text:
classdef ConstraintEnum < matlab.mixin.SetGet
    %ConstraintEnum Summary of this class goes here
    %   Detailed explanation goes here

    enumeration
        UniversalTime('Universal Time','GenericMAConstraint','Universal Time')

        PosX('Position Vector (X)','GenericMAConstraint','Position Vector (X)')
        PosY('Position Vector (Y)','GenericMAConstraint','Position Vector (Y)')
        PosZ('Position Vector (Z)','GenericMAConstraint','Position Vector (Z)')
        VelX('Velocity Vector (X)','GenericMAConstraint','Velocity Vector (X)')
        VelY('Velocity Vector (Y)','GenericMAConstraint','Velocity Vector (Y)')
        VelZ('Velocity Vector (Z)','GenericMAConstraint','Velocity Vector (Z)')

        PosVectMag('Position Vector Magnitude','GenericMAConstraint','Position Vector Magnitude')
        VelVectMag('Velocity Vector Magnitude','GenericMAConstraint','Velocity Vector Magnitude')

        SMA('Semi-major Axis','GenericMAConstraint','Semi-major Axis')
        Ecc('Eccentricity','GenericMAConstraint','Eccentricity')
        Inc('Inclination','GenericMAConstraint','Inclination')
        RAAN('Right Asc. of the Asc. Node','GenericMAConstraint','Right Asc. of the Asc. Node')
        ArgPeri('Argument of Periapsis','GenericMAConstraint','Argument of Periapsis')
        TrueAnom('True Anomaly','GenericMAConstraint','True Anomaly')
        MeanAnom('Mean Anomaly','GenericMAConstraint','Mean Anomaly')
        OrbitPeriod('Orbital Period','GenericMAConstraint','Orbital Period')

        RadiusPeriapsis('Radius of Periapsis','GenericMAConstraint','Radius of Periapsis')
        RadiusApoapsis('Radius of Apoapsis','GenericMAConstraint','Radius of Apoapsis')
        AltitudeApoapsis('Altitude of Apoapsis','GenericMAConstraint','Altitude of Apoapsis')
        AltitudePeriapsis('Altitude of Periapsis','GenericMAConstraint','Altitude of Periapsis')

        H1('Equinoctial H1','GenericMAConstraint','Equinoctial H1')
        K1('Equinoctial K1','GenericMAConstraint','Equinoctial K1')
        H2('Equinoctial H2','GenericMAConstraint','Equinoctial H2')
        K2('Equinoctial K2','GenericMAConstraint','Equinoctial K2')

        FlightPathAngle('Flight Path Angle','GenericMAConstraint','Flight Path Angle');
        Longitude('Longitude (East)','GenericMAConstraint','Longitude (East)')
        Latitude('Latitude (North)','GenericMAConstraint','Latitude (North)')
        Altitude('Altitude','GenericMAConstraint','Altitude')
        VelAz('Velocity Azimuth','GenericMAConstraint','Velocity Azimuth')
        VelEl('Velocity Elevation','GenericMAConstraint','Velocity Elevation')
        LongitudeDriftRate('Longitudinal Drift Rate','GenericMAConstraint','Longitudinal Drift Rate')

        C3Energy('C3 Energy','GenericMAConstraint','C3 Energy');
        SecPastPeri('Seconds Past Periapsis','GenericMAConstraint','Seconds Past Periapsis');

        SurfaceVelocity('Surface Velocity','GenericMAConstraint','Surface Velocity')
        VerticalVel('Vertical Velocity','GenericMAConstraint','Vertical Velocity')

        HyperVelUnitX('Hyperbolic Velocity Unit Vector X','GenericMAConstraint','Hyperbolic Velocity Unit Vector X')
        HyperVelUnitY('Hyperbolic Velocity Unit Vector Y','GenericMAConstraint','Hyperbolic Velocity Unit Vector Y')
        HyperVelUnitZ('Hyperbolic Velocity Unit Vector Z','GenericMAConstraint','Hyperbolic Velocity Unit Vector Z')
        HyperVelRA('Hyperbolic Velocity Vector Right Ascension','GenericMAConstraint','Hyperbolic Velocity Vector Right Ascension')
        HyperVelDec('Hyperbolic Velocity Vector Declination','GenericMAConstraint','Hyperbolic Velocity Vector Declination')
        HyperVelMag('Hyperbolic Velocity Magnitude','GenericMAConstraint','Hyperbolic Velocity Magnitude')

        SolarBetaAngle('Solar Beta Angle','GenericMAConstraint','Solar Beta Angle')
        TotalScMass('Total Spacecraft Mass','GenericMAConstraint','Total Spacecraft Mass')

        PitchAngle('Pitch Angle','PitchAngleConstraint',[]);
        RollAngle('Roll Angle','RollAngleConstraint',[]);
        YawAngle('Yaw Angle','YawAngleConstraint',[]);
        Throttle('Throttle','ThrottleConstraint',[]);

        BankAngle('Bank Angle','BankAngleConstraint',[]);
        AngleOfAttack('Angle of Attack','AngleOfAttackConstraint',[]);
        SideSlipAngle('Side Slip Angle','SideSlipAngleConstraint',[]);

        InertialAngleOfAttack('Inertial Angle of Attack','InertialAngleOfAttackConstraint',[]);
        InertialBankAngle('Inertial Bank Angle','InertialBankAngleConstraint',[]);
        InertialSideSlipAngle('Inertial Side Slip Angle','InertialSideSlipAngleConstraint',[]);

        TotalThrust('Total Thrust','TotalThrustConstraint',[]);
        Thrust2Weight('Sea Level Thrust to Weight','ThrustToWeightConstraint',[]);

        StopwatchValue('Stopwatch Value','StopwatchValueConstraint',[]);
        ExtremumValue('Extremum Value','ExtremumValueConstraint',[]);

        TwoBodyImpactTime('Two-Body Time To Impact','TwoBodyImpactPointTime',[]);
        TwoBodyImpactLong('Two-Body Impact Longitude','TwoBodyImpactPointLongitude',[]);
        TwoBodyImpactLat('Two-Body Impact Latitude','TwoBodyImpactPointLatitude',[]);

        EventDeltaVExpended('Event Delta-V Expended','EventDeltaVExpendedConstraint',[]);

        CalculusCalculation('Calculus Calculation','CalculusCalculationValueConstraint',[]);

        CumulativePwrStorageStateOfCharge('Electrical Power Cumulative Storage State of Charge', 'CumPwrStorageStateOfChargeConstraint', []);

        GroundObjAz('Ground Object Azimuth', 'GroundObjAzConstraint', []);
        GroundObjEl('Ground Object Elevation', 'GroundObjElConstraint', []);
        GroundObjRng('Ground Object Range', 'GroundObjRangeConstraint', []);

        GeoVectorMagntide('Geometric Vector Magnitude', 'GeometricVectorMagConstraint', []);
        GeoVectorXComponent('Geometric Vector X', 'GeometricVectorXConstraint', []);
        GeoVectorYComponent('Geometric Vector Y', 'GeometricVectorYConstraint', []);
        GeoVectorZComponent('Geometric Vector Z', 'GeometricVectorZConstraint', []);

        GeoAngleMagnitude('Geometric Angle Magnitude', 'GeometricAngleMagConstraint', []);

        TankMass('Tank Mass', 'TankMassConstraint', []);
        TankMassFlowRate('Tank Mass Flow Rate', 'TankMassFlowRateConstraint', []);

        Plugin('Plugin Value', 'PluginConstraint', []);
    end

    enumeration(Hidden)
        BodyCentricPositionX('Body-centric Position (X)','GenericMAConstraint','Body-centric Position (X)')
        BodyCentricPositionY('Body-centric Position (Y)','GenericMAConstraint','Body-centric Position (Y)')
        BodyCentricPositionZ('Body-centric Position (Z)','GenericMAConstraint','Body-centric Position (Z)')
        BodyCentricVelocityX('Body-centric Velocity (X)','GenericMAConstraint','Body-centric Velocity (X)')
        BodyCentricVelocityY('Body-centric Velocity (Y)','GenericMAConstraint','Body-centric Velocity (Y)')
        BodyCentricVelocityZ('Body-centric Velocity (Z)','GenericMAConstraint','Body-centric Velocity (Z)')
        SunCentricPositionX('Sun-centric Position (X)','GenericMAConstraint','Sun-centric Position (X)')
        SunCentricPositionY('Sun-centric Position (Y)','GenericMAConstraint','Sun-centric Position (Y)')
        SunCentricpositionZ('Sun-centric Position (Z)','GenericMAConstraint','Sun-centric Position (Z)')

        RadiusSpacecraft('Radius of Spacecraft','GenericMAConstraint','Radius of Spacecraft')
        SpeedOfSpacecraft('Speed of Spacecraft','GenericMAConstraint','Speed of Spacecraft')

        DistToRefCelBody('Distance to Ref. Celestial Body','GenericMAConstraint','Distance to Ref. Celestial Body')

        CentralBodyId('Central Body ID','GenericMAConstraint','Central Body ID')
        FuelOxMass('Liquid Fuel/Ox Mass','GenericMAConstraint','Liquid Fuel/Ox Mass')
        MonopropMass('Monopropellant Mass','GenericMAConstraint','Monopropellant Mass')
        XenonMass('Xenon Mass','GenericMAConstraint','Xenon Mass')

        TimeContinuity('Time Continuity','TimeContinuityConstraint',[]);
        PositionContinuityX('Position Continuity (X)','PositionContinuityConstraintX',[]);
        PositionContinuityY('Position Continuity (Y)','PositionContinuityConstraintY',[]);
        PositionContinuityZ('Position Continuity (Z)','PositionContinuityConstraintZ',[]);
        VelocityContinuityX('Velocity Continuity (X)','VelocityContinuityConstraintX',[]);
        VelocityContinuityY('Velocity Continuity (Y)','VelocityContinuityConstraintY',[]);
        VelocityContinuityZ('Velocity Continuity (Z)','VelocityContinuityConstraintZ',[]);
    end

    properties
        name char = '';
        class char = '';
        constructorInput1
    end

    methods
        function obj = ConstraintEnum(name, class, constructorInput1)
            obj.name = name;
            obj.class = class;
            obj.constructorInput1 = constructorInput1;
        end
    end

    methods(Static)
        function listBoxStr = getListBoxStr()
            m = enumeration('ConstraintEnum');
            [~,I] = sort({m.name});
            listBoxStr = {m(I).name};
        end

        function [ind, enum] = getIndForName(name)
            m = enumeration('ConstraintEnum');
            [~,I] = sort({m.name});
            m = m(I);
            ind = find(ismember({m.name},name),1,'first');
            enum = m(ind);
        end
    end
end

author: Arrowstar

repo: ksptot

save_path: github-repos/MATLAB/Arrowstar-ksptot


text:
function jed = datenum_to_jed ( dn )

%*****************************************************************************80
%
%% DATENUM_TO_JED converts a MATLAB date number to a JED.
%
%  Discussion:
%
%    The MATLAB "datenum" function accepts a string defining
%    a date and returns a datenumber:
%
%      dn = datenum ( 'Aug 17 1939' )
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    15 January 2009
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, real DN, a MATLAB date number.
%
%    Output, real JED, the Julian Ephemeris Date.
%
  jed = dn + 1721058.5;

  return
end


text:
function pos = getComPosition(obj)
    % Returns the symbolic representation of the robot manipulator
    %
    % Return values:
    % pos: the 3-Dimensional SO(3) position vector of the CoM of the system
    % @type SymExpression

    n_link = length(obj.Links);
    links = cell(1,n_link);
    for i=1:n_link
        links{i}.Mass = obj.Links(i).Mass;
        % links{i}.Inertia = obj.Links(i).Inertia;
        links{i}.gst0 = obj.Links(i).gst0;
        frame = obj.Links(i).Reference;
        while isempty(frame.TwistPairs)
            frame = frame.Reference;
            if isempty(frame)
                error('The coordinate system is not fully defined.');
            end
        end

        links{i}.TwistPairs = frame.TwistPairs;
        % links{i}.ChainIndices = frame.ChainIndices;
    end

    pos = eval_math_fun('ComputeComPosition',links);

end

author: ayonga

repo: frost-dev

save_path: github-repos/MATLAB/ayonga-frost-dev


text:
function Pvalue=myfisher33(x)
%P=MYFISHER33(X)- Fisher's Exact Probability Test on 3x3 matrix.
%Fisher's exact test of 3x3 contingency tables permits calculation of
%precise probabilities in situation where, as a consequence of small cell
%frequencies, the much more rapid normal approximation and chi-square
%calculations are liable to be inaccurate. The Fisher's exact test involves
%the computations of several factorials to obtain the probability of the
%observed and each of the more extreme tables. Factorials growth quickly,
%so it's necessary use logarithms of factorials. This computations is very
%easy in Matlab because x!=gamma(x+1) and log(x!)=gammaln(x+1). This
%function is now fully vectorized to speed up the computation.
% Syntax: 	myfisher33(x)
%
%     Inputs:
%           X - 3x3 data matrix
%     Outputs:
%           P - 2-tailed p-value
%
%   Example:
%
%                A   B   C
%           -------------------
%      X         4   1   1
%           -------------------
%      Y         1   3   1
%           -------------------
%      Z         0   1   5
%           -------------------
%
%   x=[4 1 1; 1 3 1; 0 1 5];
%
%   Calling on Matlab the function:
%             myfisher33(x)
%
% 3x3 matrix Fisher's exact test: 301 tables were evaluated
% -------------------------------------------------------
% 2-tail p-value = 0.0322177822
% -------------------------------------------------------
%
%           Created by Giuseppe Cardillo
%           giuseppe.cardillo-edta@poste.it
%
% To cite this file, this would be an appropriate format:
% Cardillo G. (2007) MyFisher33: a very compact routine for Fisher's exact
% test on 3x3 matrix
% http://www.mathworks.com/matlabcentral/fileexchange/15482

%Input Error handling
if ~isequal(size(x),[3 3])
    error('Input matrix must be a 3x3 matrix')
end
if ~all(isfinite(x(:))) || ~all(isnumeric(x(:)))
    error('Warning: all X values must be numeric and finite')
end
if ~isequal(x(:),round(x(:)))
    error('Warning: X data matrix values must be whole numbers')
end

Rs=sum(x,2); %rows sum
Cs=sum(x); %columns sum
N=sum(Rs); %Number of observations

%If necessed, rearrange matrix
if ~issorted(Cs)
    [Cs,ind]=sort(Cs);
    x=x(:,ind);
end
if ~issorted(Rs)
    [Rs,ind]=sort(Rs);
    x=x(ind,:);
end

Kf=sum(gammaln([Rs' Cs]+1))-gammaln(N+1); %The costant factor K=log(prod(R!)*prod(C!)/N!)
zf=gammaln(x+1); %compute log(x!)
op=exp(Kf-sum(zf(:))); %compute the p-value of the observed matrix

% First step: Compute all possible primary tables (play on the first two
% degrees of freedom)
I=0:1:min(Rs(1),Cs(1)); %all possible values of X(1,1)
J=min(Rs(1)-I,Cs(2)); %max value of X(1,2) given X(1,1)
et=sum(J+1); %primary tables to evaluate
m1=zeros(et,3); %matrix preallocation
%set the arrays of indices
idxAstop=cumsum(J+1);
idxAstart=[1 idxAstop(1:end-1)+1];
%Build-up the matrix
for K=1:length(I)
    m1(idxAstart(K):idxAstop(K),1)=I(K);
    m1(idxAstart(K):idxAstop(K),2)=0:1:J(K);
end
%Put all the possible values of X(1,3) given X(1,1) and X(1,2)
m1(:,3)=Rs(1)-sum(m1(:,1:2),2);

%Second step:
%Compute all possible secondary tables (play on the second two
% degrees of freedom)
m2=repmat(Cs,et,1);
Cs2=m2-m1; %Residual sum of rows given the first row
L=min(Rs(2),Cs2(:,1)); %find the max values of X(2,1) given X(1,1) and X(1,2)
KK=max(L)+1; %Remember that zero is a possible value...
%matrix preallocation
UpB=zeros(et,KK);
for K=1:KK
%find the max values of X(2,2) given X(1,1), X(1,2) and X(2,2)
    UpB(:,K)=min(Cs2(:,2),Rs(2)-(K-1));
end
%matrix preallocation
LoB=zeros(et,KK); InT=repmat(NaN,et,KK);
sectab=zeros(1,et); %array preallocation
for K=1:et
    %find the min values of X(2,2) given X(1,1), X(1,2) and X(2,2)
    z=L(K)+1;
    a=0:1:L(K);
    LoB(K,1:z)=max(0,Rs(2)-a-Cs2(K,3));
    %Compute the range of variation of X(2,2)
    InT(K,1:z)=UpB(K,1:z)-LoB(K,1:z)+1;
    %Compute how many secondary tables are generated from each primary table
    sectab(K)=sum(InT(K,1:z),2);
end
et2=sum(sectab); %total tables to evaluate
%Matrix and vector preallocation
M=zeros(et2,9);
%Set indices
InT=InT'; InT(isnan(InT))=[];
idxAstop=cumsum(InT(:))';
idxAstart=[1 idxAstop(1:end-1)+1];
idxBstop=cumsum(sectab); idxBstart=[1 idxBstop(1:end-1)+1];
%Build-up the secondary tables
z=1;
for K=1:et
    %Expand the primary table
    M(idxBstart(K):idxBstop(K),1:3)=repmat(m1(K,:),sectab(K),1);
    I=0:1:L(K); %values of X(2,1) given X(1,1) and X(1,2)
    %Put them in the matrix in the correct position
    for J=1:length(I)
        M(idxAstart(z):idxAstop(z),4)=I(J);
        M(idxAstart(z):idxAstop(z),5)=LoB(K,J):1:UpB(K,J); %values of X(2,2) given X(1,1), X(1,2) and X(2,2)
        z=z+1;
    end
end
%Complete the table
M(:,6)=Rs(2)-sum(M(:,4:5),2);
M(:,7:9)=repmat(Cs,et2,1)-M(:,1:3)-M(:,4:6);
zf=gammaln(M+1); %compute log(x!)
np=exp(Kf-sum(zf,2)); %compute the p-value of each possible matrix
P=sum(np(np<op));

%display results
tr=repmat('-',1,55); %Set up the divisor
fprintf('3x3 matrix Fisher''s exact test: %0.0f tables were evaluated\n',et2)
disp(tr)
fprintf('2-tail p-value = %0.10f\n',P);
disp(tr)
if nargout
    Pvalue=P;
end

author: Sable

repo: mcbench-benchmarks


text:
%
% [A,b,c,K]=convertf(A,b,c,K)
%
% converts free variables in a SeDuMi problem into nonnegative LP variables.
%
function [A,b,c,K]=convertf(A,b,c,K)
%
% Get the number of constraints.
%
m=length(b);
%
%  Deal with the following special case.  If A is transposed, transpose
%  it again so that it is of the right size.
%
[Am,An]=size(A);
if (Am ~= m)
  if (An == m)
    fprintf('Transposing A to match b \n');
    A=A';
  else
    fprintf('A is not of the correct size to match b \n');
    return
  end
end
%
%  Deal with the following special case:  if c==0, then c should really
%  be a zero vector of the appropriate size.
%
if (c == 0)
  fprintf('Expanding c to the appropriate size\n');
  [Am,An]=size(A);
  c=zeros(An,1);
end
%
% If c is empty, then act as if it was zero.
%
if (isempty(c))
  fprintf('Expanding empty c to zeros of the appropriate size\n');
  [Am,An]=size(A);
  c=zeros(An,1);
end
%
% If c is a row vector, make it a column vector.
%
[cm,cn]=size(c);
if (cn > cm)
  c=c';
end
%
% Check for any free LP variables and rewrite them as the differences of
% regular LP variables.
%
if (isfield(K,'f'))
  nfree=K.f
  fprintf('Converting %d free variables to LP variables\n',nfree);
  if (isfield(K,'l'))
    nlin=K.l;
  else
    nlin=0;
  end
  [Am,An]=size(A);
  Anew=[A(:,1:nfree) -A(:,1:nfree) A(:,nfree+1:An)];
  A=Anew;
  cnew=[c(1:nfree); -c(1:nfree); c(nfree+1:An)];
  c=cnew;

  K.l=nlin+2*nfree;
  K.f=0;
end

author: vigente

repo: gerardus

save_path: github-repos/MATLAB/vigente-gerardus


text:
function mm_io_test10 ( )

%*****************************************************************************80
%
%% MM_IO_TEST10 tests MM_FILE_READ.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    01 May 2004
%
%  Author:
%
%    John Burkardt
%
  fprintf ( 1, '\n' );
  fprintf ( 1, 'MM_IO_TEST10:\n' );
  fprintf ( 1, '  MM_VALUE_READ reads the value records \n' );
  fprintf ( 1, '  of a Matrix Market file.\n' );

  input_filename = 'mm_io_test09.mm';

  input_unit = fopen ( input_filename );

  if ( input_unit < 0 )
    fprintf ( 1, '\n' );
    fprintf ( 1, 'MM_IO_TEST10 - Fatal error!\n' );
    fprintf ( 1, '  Could not open the input file.\n' );
    error ( 'MM_IO_TEST10 - Fatal error!' );
  end

  [ id, type, rep, field, symm, nrow, ncol, nnz, indx, jndx, value ) ...
    = mm_file_read ( input_unit );

  fclose ( input_unit );

  fprintf ( 1, '\n' );
  fprintf ( 1, '  Here is the information read from the file:\n' );
  fprintf ( 1, '\n' );

  mm_header_print ( input_file, id, type, rep, field, symm );

  mm_size_print ( input_filename, rep, symm, nrow, ncol, nnz );

  mm_values_print ( rep, field, nnz, indx, jndx, value );

  fprintf ( 1, '\n' );
  fprintf ( 1, 'MM_IO_TEST10:\n' );
  fprintf ( 1, '  The MM file "%s" containing \n', input_filename );
  fprintf ( 1, '  the header, comments, size and values, has been read.\n' );
  fprintf ( 1, '\n' );

  return
end


text:
% Demo1: Dimensionality reduction or feature extraction using KPCA
% X: training samples
%
%%%%%%%%%%%%%%%%%%%%%%%%%%     KPCA       %%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Improve the performance by adjusting the following parameters:
% 1. sigma  % kernel width
% 2. dims   % Dimensionality
% ---------------------------------------------------------------------%

clc
clear all
close all

addpath(genpath(pwd))

% 4 circles
load circledata
X = circledata;
figure
for i = 1:4
    scatter(X(1+250*(i-1):250*i,1),X(1+250*(i-1):250*i,2))
    hold on
end

% Train KPCA model
model = kpca_train(X,'type',0,'sigma',5,'dims',2);

% Visualize the result of dimensionality reduction
figure
for i = 1:4
    scatter(model.mappedX(1+250*(i-1):250*i,1), ...
        model.mappedX(1+250*(i-1):250*i,2))
    hold on
end

author: zhangqianqianQQ

repo: MachineVisionAlgorithm


text:
function ds_info=gen_dataset_info(ds_config)

if ds_config.use_custom_data

    ds_info=gen_ds_info_custom_data(ds_config, ds_config.class_info);

else

    cache_dir=ds_config.ds_info_cache_dir;
    cache_file=fullfile(cache_dir, 'my_dataset_info.mat');

    if ~my_check_file(cache_file)
        ds_info=ds_config.gen_ds_info_fn(ds_config);
        my_save_file(cache_file, ds_info);
    else
        ds_info=my_load_file(cache_file);
        ds_info=ds_info.data_obj;
    end

end

end

function ds_info=gen_ds_info_custom_data(ds_config, class_info)

    ds_info=[];

    img_dir=ds_config.img_data_dir;
    img_files=my_list_file(img_dir)';
    img_num=length(img_files);
    ds_info.img_idxes=uint32(1:img_num);
    ds_info.img_num=img_num;
    ds_info.train_idxes=uint32([]);
    ds_info.test_idxes=ds_info.img_idxes;

    img_names=cell(img_num, 1);
    for img_idx=1:img_num
        one_name=img_files{img_idx};
        [~, one_name]=fileparts(one_name);
        img_names{img_idx}=one_name;
    end

    ds_info.img_files=img_files;
    ds_info.img_names=img_names;
    ds_info.mask_files=[];
    ds_info.train_idxes=[];
    ds_info.test_idxes=uint32(1:img_num)';
    ds_info.name=ds_config.ds_name;

    ds_info.data_dir_idxes_img=repmat(uint8(1), [img_num, 1]);
    ds_info.data_dir_idxes_mask=[];
    ds_info.data_dirs={img_dir};

    ds_info.class_info=class_info;

end

author: guosheng

repo: refinenet

save_path: github-repos/MATLAB/guosheng-refinenet


text:
%Copyright 2008-2009  Marius Muja (mariusm@cs.ubc.ca). All rights reserved.
%Copyright 2008-2009  David G. Lowe (lowe@cs.ubc.ca). All rights reserved.
%
%THE BSD LICENSE
%
%Redistribution and use in source and binary forms, with or without
%modification, are permitted provided that the following conditions
%are met:
%
%1. Redistributions of source code must retain the above copyright
%   notice, this list of conditions and the following disclaimer.
%2. Redistributions in binary form must reproduce the above copyright
%   notice, this list of conditions and the following disclaimer in the
%   documentation and/or other materials provided with the distribution.
%
%THIS SOFTWARE IS PROVIDED BY THE AUTHOR ``AS IS'' AND ANY EXPRESS OR
%IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
%OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
%IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY DIRECT, INDIRECT,
%INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT
%NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
%DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
%THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
%(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF
%THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.

function [indices, dists] = flann_search(data, testset, n, search_params)
%NN_SEARCH  Fast approximate nearest neighbors search
%
% Performs a fast approximate nearest neighbor search using an
% index constructed using flann_build_index or directly a
% dataset.

% Marius Muja, January 2008

    algos = struct( 'linear', 0, 'kdtree', 1, 'kmeans', 2, 'composite', 3, 'saved', 254, 'autotuned', 255 );
    center_algos = struct('random', 0, 'gonzales', 1, 'kmeanspp', 2 );
    log_levels = struct('none', 0, 'fatal', 1, 'error', 2, 'warning', 3, 'info', 4);

    default_params = struct('algorithm', 'kdtree' ,'checks', 32, 'eps', 0.0, 'sorted', 1, 'max_neighbors', -1, 'cores', 1, 'trees', 4, 'branching', 32, 'iterations', 5, 'centers_init', 'random', 'cb_index', 0.4, 'target_precision', 0.9,'build_weight', 0.01, 'memory_weight', 0, 'sample_fraction', 0.1, 'log_level', 'warning', 'random_seed', 0);

    if ~isstruct(search_params)
        error('The "search_params" argument must be a structure');
    end

    params = default_params;
    fn = fieldnames(search_params);
    for ii = 1:length(fn),
        name = cell2mat(fn(ii));
        params.(name) = search_params.(name);
    end
    if ~isnumeric(params.algorithm),
        params.algorithm = value2id(algos,params.algorithm);
    end
    if ~isnumeric(params.centers_init),
        params.centers_init = value2id(center_algos,params.centers_init);
    end
    if ~isnumeric(params.log_level),
        params.log_level = value2id(log_levels,params.log_level);
    end

    if (size(data,1)==1 && size(data,2)==1)
        % we already have an index
        [indices,dists] = nearest_neighbors('index_find_nn', data, testset, n, params);
    else
        % create the index and search
        [indices,dists] = nearest_neighbors('find_nn', data, testset, n, params);
    end
end

function id = value2id(map,value)
    id = map.(value);
end

function value = id2value(map, id)
    fields = fieldnames(map);
    for ii = 1:length(fields),
        val = cell2mat(fields(ii));
        if map.(val) == id
            value = val;
            break;
        end
    end
end

repo: gspbox


text:
function [p] = gui_eeg_open(p,command,parent)

% gui_eeg_open - Load an EEG data file into matlab workspace
%
% Usage: [p] = gui_eeg_open(p,[command],[parent])
%
% p is a structure, generated by 'eeg_toolbox_defaults'
% command is either 'init' or 'load'
% parent is a handle to the gui that calls this gui, useful
% for updating the UserData field of the parent from this gui.
% The p structure is returned to the parent when the parent
% handle is given.
%

% $Revision: 1.1 $ $Date: 2009-04-28 22:13:56 $

% Licence:  GNU GPL, no express or implied warranties
% History:  01/2002, Darren.Weber_at_radiology.ucsf.edu
%           08/2002, Darren.Weber_at_radiology.ucsf.edu
%                    added EMSE avg handling
%                    added interpolation of zero point option
%           10/2003, Darren.Weber_at_radiology.ucsf.edu
%                    created gui_avg_open and adapted this
%                    function to handle CNT/EEG data formats
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

if ~exist('p','var'),
 [p] = eeg_toolbox_defaults;
elseif isempty(p),
 [p] = eeg_toolbox_defaults;
end

if ~exist('command','var'), command = 'init'; end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Paint the GUI
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

switch command,

  case 'init',

    if exist('parent','var'),
      EEGOpen = INIT(p,parent);
    else
      EEGOpen = INIT(p,'');
    end

  case 'plot',

    EEGOpen = get(gcbf,'Userdata');

    EEGOpen.p = cnt_open(EEGOpen.p,EEGOpen.gui);
    set(EEGOpen.gui,'Userdata',EEGOpen); % update GUI with CNT data

   [p] = gui_updateparent(EEGOpen,0);

    if isequal(get(EEGOpen.handles.Bhold,'Value'),0),
      close gcbf;
      if isfield(EEGOpen,'parent'),
        parent = EEGOpen.parent.gui;
      else
        parent = [];
      end
    else
      parent = EEGOpen.gui;
    end

    plotfig = figure('Name',EEGOpen.p.cnt.file,...
      'NumberTitle','off',...
      'UserData',EEGOpen.p);
    movegui(plotfig,'center');

    plot(EEGOpen.p.cnt.timeArray,EEGOpen.p.cnt.data); axis tight;
    eeg_plot_metric;

    if isempty(parent), [Xpoint,Ypoint] = eeg_crosshair('init',EEGOpen.p);
    else                [Xpoint,Ypoint] = eeg_crosshair('init',EEGOpen.p,parent);
    end

  case 'save',

    fprintf('\nGUI_EEG_OPEN: Save As not implemented yet.\n');

  case 'return',

    EEGOpen = get(gcbf,'Userdata');

    set(EEGOpen.gui,'Pointer','watch');

    EEGOpen.p = cnt_open(EEGOpen.p);

    set(EEGOpen.gui,'Pointer','arrow');

    if isequal(get(EEGOpen.handles.Bhold,'Value'),1),
     [p] = gui_updateparent(EEGOpen,0);
    else
     [p] = gui_updateparent(EEGOpen);
      close gcbf;
    end

  otherwise,

    EEGOpen = get(gcbf,'Userdata');
    GUI.parent = EEGOpen.parent;
    gui_updateparent(GUI);
    close gcbf;

end

return

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [EEGOpen] = INIT(p,parent)
% GUI General Parameters

GUIwidth  = 500;
GUIheight = 120;

version = '$Revision: 1.1 $';
name = sprintf('EEG File Open [v %s]\n',version(11:15));

GUI = figure('Name',name,'Tag','EEG_OPEN',...
  'NumberTitle','off',...
  'MenuBar','none','Position',[1 1 GUIwidth GUIheight]);
movegui(GUI,'center');

Font.FontName   = 'Helvetica';
Font.FontUnits  = 'Pixels';
Font.FontSize   = 12;
Font.FontWeight = 'normal';
Font.FontAngle  = 'normal';

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% CNT Data Selection and Parameters

G.Title_data = uicontrol('Parent',GUI,'Style','text','Units','Normalized',Font, ...
  'Position',[.01 .75 .17 .2],...
  'String','Data Type:','HorizontalAlignment','left');

switch lower(p.cnt.type),
  case 'scan4x_cnt',  cntType = 1;
  case 'scan3x_cnt',  cntType = 2;
  otherwise,          cntType = 1;
end

%  'String',{'ASCII_CNT' 'Scan4x_CNT' 'Scan3x_CNT' 'EEGLAB'},'Value',cntType,...

G.PcntType = uicontrol('Tag','PcntType','Parent',GUI,'Style','popupmenu',...
  'Units','Normalized',Font,  ...
  'Position',[.20 .75 .2 .2],...
  'String',{'scan4x_cnt' 'scan3x_cnt'},'Value',cntType,...
  'Callback',strcat('EEGOpen = get(gcbf,''Userdata'');',...
  'EEGOpen.p.cnt.type = popupstr(EEGOpen.handles.PcntType);',...
  'set(gcbf,''Userdata'',EEGOpen); clear EEGOpen;'));

G.Title_path = uicontrol('Parent',GUI,'Style','text','Units','Normalized',Font, ...
  'Position',[.01 .50 .17 .2],...
  'String','Path','HorizontalAlignment','left');
G.EcntPath = uicontrol('Parent',GUI,'Style','edit','Units','Normalized',Font,  ...
  'Position',[.20 .50 .58 .2], 'String',p.cnt.path,...
  'Callback',strcat('EEGOpen = get(gcbf,''Userdata'');',...
  'EEGOpen.p.cnt.path = get(EEGOpen.handles.EcntPath,''String'');',...
  'set(gcbf,''Userdata'',EEGOpen); clear EEGOpen;'));

G.Title_file = uicontrol('Parent',GUI,'Style','text','Units','Normalized',Font, ...
  'Position',[.01 .25 .17 .2],...
  'String','File','HorizontalAlignment','left');
G.EcntFile = uicontrol('Parent',GUI,'Style','edit','Units','Normalized',Font,  ...
  'Position',[.20 .25 .58 .2], 'String',p.cnt.file,...
  'Callback',strcat('EEGOpen = get(gcbf,''Userdata'');',...
  'EEGOpen.p.cnt.file = get(EEGOpen.handles.EcntFile,''String'');',...
  'set(gcbf,''Userdata'',EEGOpen); clear EEGOpen;'));

Font.FontWeight = 'bold';

% BROWSE: Look for the data
browsecommand = strcat('EEGOpen = get(gcbf,''Userdata'');',...
  'cd(EEGOpen.p.cnt.path);',...
  '[file, path] = uigetfile(',...
  '{''*.cnt'', ''NeuroScan (*.cnt)'';', ...
  ' ''*.*'',   ''All Files (*.*)''},', ...
  '''Select EEG File'');',...
  'if ~isequal(path,0), EEGOpen.p.cnt.path = path; end;',...
  'if ~isequal(file,0), EEGOpen.p.cnt.file = file; end;',...
  'set(EEGOpen.handles.EcntPath,''String'',EEGOpen.p.cnt.path);',...
  'set(EEGOpen.handles.EcntFile,''String'',EEGOpen.p.cnt.file);',...
  'set(gcbf,''Userdata'',EEGOpen); clear EEGOpen file path;');
G.BcntFile = uicontrol('Parent',GUI,'Style','pushbutton','Units','Normalized',Font, ...
  'Position',[.01 .01 .17 .2], 'String','BROWSE',...
  'BackgroundColor',[0.8 0.8 0.0],...
  'ForegroundColor', [1 1 1], 'HorizontalAlignment', 'center',...
  'Callback', browsecommand );

% PLOT: Load & plot the data!
G.Bplot = uicontrol('Parent',GUI,'Style','pushbutton','Units','Normalized', Font, ...
  'Position',[.20 .01 .18 .2],...
  'String','PLOT','BusyAction','queue',...
  'TooltipString','Plot the EEG data and return p struct.',...
  'BackgroundColor',[0.0 0.5 0.0],...
  'ForegroundColor', [1 1 1], 'HorizontalAlignment', 'center',...
  'Callback',strcat('EEGOpen = get(gcbf,''Userdata'');',...
  'p = gui_eeg_open(EEGOpen.p,''plot'');',...
  'clear EEGOpen;'));

% Save As
G.Bsave = uicontrol('Parent',GUI,'Style','pushbutton','Units','Normalized', Font, ...
  'Position',[.40 .01 .18 .2],'HorizontalAlignment', 'center',...
  'String','SAVE AS','TooltipString','EEG File Conversion Tool (not implemented yet)',...
  'BusyAction','queue',...
  'Visible','off',...
  'BackgroundColor',[0.0 0.0 0.75],...
  'ForegroundColor', [1 1 1], 'HorizontalAlignment', 'center',...
  'Callback',strcat('EEGOpen = get(gcbf,''Userdata'');',...
  'p = gui_eeg_open(EEGOpen.p,''save'');',...
  'clear EEGOpen;'));

% Quit, return file parameters
G.Breturn = uicontrol('Parent',GUI,'Style','pushbutton','Units','Normalized', Font, ...
  'Position',[.60 .01 .18 .2],...
  'String','RETURN','BusyAction','queue',...
  'TooltipString','Return p struct to workspace and parent GUI.',...
  'BackgroundColor',[0.75 0.0 0.0],...
  'ForegroundColor', [1 1 1], 'HorizontalAlignment', 'center',...
  'Callback',strcat('EEGOpen = get(gcbf,''Userdata'');',...
  'p = gui_eeg_open(EEGOpen.p,''return'');',...
  'clear EEGOpen;'));

% Cancel
G.Bcancel = uicontrol('Parent',GUI,'Style','pushbutton','Units','Normalized', Font, ...
  'Position',[.80 .01 .18 .2],...
  'String','CANCEL','BusyAction','queue',...
  'TooltipString','Close, do not return parameters.',...
  'BackgroundColor',[0.75 0.0 0.0],...
  'ForegroundColor', [1 1 1], 'HorizontalAlignment', 'center',...
  'Callback',strcat('EEGOpen = get(gcbf,''Userdata'');',...
  'p = gui_eeg_open(EEGOpen.p,''cancel'');',...
  'clear EEGOpen;'));

% Hold GUI Open checkbox
G.Bhold = uicontrol('Parent',GUI,'Style','checkbox','Units','Normalized', Font, ...
  'Position',[.80 .25 .18 .2],'String','Hold GUI','BusyAction','queue',...
  'TooltipString','EEG File Load GUI remains open after ''Plot'' or ''Return'' commands.',...
  'Value',p.hold,'HorizontalAlignment', 'center');

% Store userdata
if exist('parent','var'), EEGOpen.parent.gui = parent; end
EEGOpen.gui = GUI;
EEGOpen.handles = G;
EEGOpen.p = p;
set(GUI,'Userdata',EEGOpen);
set(GUI,'HandleVisibility','callback');

return

author: PatternRecognition

repo: OpenBMI


text:
function exact = p06_exact ( )

%*****************************************************************************80
%
%% P06_EXACT returns the exact integral for problem 6.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    28 July 2007
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Output, real EXACT, the estimated value of the integral.
%
  exact = 0.00056103711148387120640;

  return
end


text:
function [ ap, det, inert ] = zhpdi ( ap, n, ipvt, job )

%*****************************************************************************80
%
%% ZHPDI: determinant, inertia and inverse of a complex hermitian matrix.
%
%  Discussion:
%
%    The routine uses the factors from ZHPFA.
%
%    The matrix is stored in packed form.
%
%    A division by zero will occur if the inverse is requested and ZHPCO has
%    set RCOND == 0.0 or ZHPFA has set INFO ~= 0.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    13 May 2007
%
%  Author:
%
%    MATLAB version by John Burkardt
%
%  Reference:
%
%    Jack Dongarra, Jim Bunch, Cleve Moler, Pete Stewart,
%    LINPACK User's Guide,
%    SIAM, 1979,
%    ISBN13: 978-0-898711-72-1,
%    LC: QA214.L56.
%
%  Parameters:
%
%    Input, complex AP(N*(N+1)/2); the factored matrix
%    from ZHPFA.
%
%    Input, integer N, the order of the matrix.
%
%    Input, integer IPVT(N), the pivot vector from ZHPFA.
%
%    Input, integer JOB, has the decimal expansion ABC where:
%    if C ~= 0, the inverse is computed,
%    if B ~= 0, the determinant is computed,
%    if A ~= 0, the inertia is computed.
%    For example, JOB = 111 gives all three.
%
%    Output, complex AP(N*(N+1)/2); if the inverse was requested, then
%    the upper triangle of the inverse of the original matrix, stored in packed
%    form.  The columns of the upper triangle are stored sequentially in a
%    one-dimensional array.
%
%    Output, real DET(2), if requested, the determinant of the original matrix.
%    Determinant = DET(1) * 10.0**DET(2) with 1.0 <= abs ( DET(1) ) < 10.0
%    or DET(1) = 0.0.
%
%    Output, integer INERT(3), if requested, the inertia of the original matrix.
%    INERT(1) = number of positive eigenvalues.
%    INERT(2) = number of negative eigenvalues.
%    INERT(3) = number of zero eigenvalues.
%
  det = [];
  inert = [];

  noinv = mod ( job, 10 ) == 0;
  nodet = floor ( mod ( job, 100 ) / 10 ) == 0;
  noert = floor ( mod ( job, 1000 ) / 100 ) == 0;

  if ( ~nodet | ~noert )

    if ( ~noert )
      inert(1:3) = 0;
    end

    if ( ~nodet )
      det(1) = 1.0;
      det(2) = 0.0;
    end

    t = 0.0;
    ik = 0;

    for k = 1 : n

      kk = ik + k;
      d = real ( ap(kk) );
%
%  Check if 1 by 1
%
      if ( ipvt(k) <= 0 )
%
%  2 by 2 block
%  Use DET (D  S; S  C)  =  ( D / T * C - T ) * T, T = abs ( S )
%  to avoid underflow/overflow troubles.
%  Take two passes through scaling.  Use T for flag.
%
        if ( t == 0.0 )
          ikp1 = ik + k;
          kkp1 = ikp1 + k;
          t = abs ( ap(kkp1) );
          d = ( d / t ) * real ( ap(kkp1+1) ) - t;
        else
          d = t;
          t = 0.0;
        end

      end

      if ( ~noert )

        if ( 0.0 < d )
          inert(1) = inert(1) + 1;
        elseif ( d < 0.0 )
          inert(2) = inert(2) + 1;
        elseif ( d == 0.0 )
          inert(3) = inert(3) + 1;
        end

      end

      if ( ~nodet )

        det(1) = det(1) * d;

        if ( det(1) ~= 0.0 )

          while ( abs ( det(1) ) < 1.0 )
            det(1) = det(1) * 10.0;
            det(2) = det(2) - 1.0;
          end

          while ( 10.0 <= abs ( det(1) ) )
            det(1) = det(1) / 10.0;
            det(2) = det(2) + 1.0;
          end

        end

      end

      ik = ik + k;

    end

  end
%
%  Compute inverse(A).
%
  if ( ~noinv )

    k = 1;
    ik = 0;

    while ( k <= n )

      km1 = k - 1;
      kk = ik + k;
      ikp1 = ik + k;
      kkp1 = ikp1 + k;
%
%  1 by 1
%
      if ( 0 <= ipvt(k) )

        ap(kk) = 1.0 / real ( ap(kk) );

        if ( 1 <= km1 )

          work(1:km1) = ap(ik+1:ik+km1);

          ij = 0;
          for j = 1 : km1
            jk = ik + j;
            ap(jk) = conj ( ap(ij+1:ij+j) ) * tranpose ( work(1:j) );
            ap(ik+1:ik+j-1) = ap(ik+1:ik+j-1) + work(j) * ap(ij+1:ij+j-1);
            ij = ij + j;
          end

          ap(kk) = ap(kk) + real ( conj ( work(1:km1) ) * transpose ( ap(ik+1:ik+km1) ) );

        end

        kstep = 1;
%
%  2 by 2
%
      else

        t = abs ( ap(kkp1) );
        ak = real ( ap(kk) ) / t;
        akp1 = real ( ap(kkp1+1) ) / t;
        akkp1 = ap(kkp1) / t;
        d = t * ( ak * akp1 - 1.0 );
        ap(kk) = akp1 / d;
        ap(kkp1+1) = ak / d;
        ap(kkp1) = -akkp1 / d;

        if ( 1 <= km1 )

          work(1:km1) = ap(ikp1+1:ikp1+km1);

          ij = 0;
          for j = 1 : km1
            jkp1 = ikp1 + j;
            ap(jkp1) = conj ( ap(ij+1:ij+j) ) * transpose ( work(1:j) );
            ap(ikp1+1:ikp1+j-1) = ap(ikp1+1:ikp1+j-1) ...
              + work(j) * ap(ij+1:ij+j-1);
            ij = ij + j;
          end

          ap(kkp1+1) = ap(kkp1+1) ...
            + real ( conj ( work(1:km1) ) * transpose ( ap(ikp1+1) ) );

          ap(kkp1) = ap(kkp1) ...
            + conj ( ap(ik+1:ik+km1) ) * transpose ( ap(ikp1+1:ikp1+km1) );

          work(1:km1) = ap(ik+1:ik+km1);

          ij = 0;

          for j = 1 : km1
            jk = ik + j;
            ap(jk) = conj ( ap(ij+1:ij+j) ) * transpose ( work(1:j) );
            ap(ik+1:ik+j-1) = ap(ik+1:ik+j-1) + work(j) * ap(ij+1:ij+j-1);
            ij = ij + j;
          end

          ap(kk) = ap(kk) ...
            + real ( conj ( work(1:km1) ) * transpose ( ap(ik+1:ik+km1) ) );

        end

        kstep = 2;

      end
%
%  Swap
%
      ks = abs ( ipvt(k) );

      if ( ks ~= k )

        iks = ( ks * ( ks - 1 ) ) / 2;

        temp             = ap(iks+1:iks+ks);
        ap(iks+1:iks+ks) = ap(ik+1:ik+ks);
        ap(ik+1:ik+ks)   = temp;

        ksj = ik + ks;

        for jb = ks : k

          j = k + ks - jb;
          jk = ik + j;

          temp    = conj ( ap(jk) );
          ap(jk)  = conj ( ap(ksj) );
          ap(ksj) = temp;

          ksj = ksj - ( j - 1 );

        end

        if ( kstep ~= 1 )

          kskp1 = ikp1 + ks;

          temp      = ap(kskp1);
          ap(kskp1) = ap(kkp1);
          ap(kkp1)  = temp;

        end

      end

      ik = ik + k;

      if ( kstep == 2 )
        ik = ik + k + 1;
      end

      k = k + kstep;

    end

  end

  return
end


text:
function copyOrbitToClipboardFromText(hEpoch, hSMA, hEcc, hInc, hRAAN, hArg, hAnomaly, isTrueAnomaly, bodyID)
%copyOrbitToClipboardFromText Summary of this function goes here
%   Detailed explanation goes here
    global GLOBAL_OrbitClipboard;

    if(~isempty(hEpoch))
        GLOBAL_OrbitClipboard(1) = str2double(get(hEpoch, 'String'));
    else
        GLOBAL_OrbitClipboard(1) = 0;
    end

    GLOBAL_OrbitClipboard(2) = str2double(get(hSMA, 'String'));
    GLOBAL_OrbitClipboard(3) = str2double(get(hEcc, 'String'));
    GLOBAL_OrbitClipboard(4) = deg2rad(str2double(get(hInc, 'String')));
    GLOBAL_OrbitClipboard(5) = deg2rad(str2double(get(hRAAN, 'String')));
    GLOBAL_OrbitClipboard(6) = deg2rad(str2double(get(hArg, 'String')));

    if(~isempty(hAnomaly))
        anomaly = deg2rad(str2double(get(hAnomaly, 'String')));
        if(isTrueAnomaly == true)
            tru = anomaly;
        else
            mean = anomaly;
            ecc = str2double(get(hEcc, 'String'));
            tru = computeTrueAnomFromMean(mean, ecc);
        end
        GLOBAL_OrbitClipboard(7) = tru;
    else
        GLOBAL_OrbitClipboard(7) = 0;
    end

    GLOBAL_OrbitClipboard(8) = bodyID;
end


text:
classdef Ephem < handle
% Class for opening JPL ephemeris files and calculating raw positions
% Copyright(c) 2014 by Jonathan Kipling Knight (drkipknight@aol.com)
% DEPENDENCIES: This class requires that ephemeris files are available online by
% ftp or have been previously downloaded and placed in the Matlab path.
  properties (Constant)
    % body numbers for use with planet_ephemeris

    Mercury               = 1; % Mercury body number
    Venus                 = 2; % Venus body number
    Earth                 = 3; % Earth body number
    Mars                  = 4; % Mars body number
    Jupiter               = 5; % Jupiter body number
    Saturn                = 6; % Saturn body number
    Uranus                = 7; % Uranus body number
    Neptune               = 8; % Neptune body number
    Pluto                 = 9; % Pluto body number
    Moon                  = 10; % Geocentric Moon body number
    Sun                   = 11; % Sun body number
    SolarSystemBarycenter = 12; % Solarsystem Barycenter body number
    EarthMoonBarycenter   = 13; % Earth-Moon Barycenter body number
    Nutations             = 14; % nutations body number
    Librations            = 15; % librations body number
    LunarEulerRates       = 16; % lunar euler angel rates body number
    TTmTDB                = 17; % TT-TDB body number
  end
  properties (Constant)
    % state numbers for use with state or interpolate

    MercuryState               = 1; % Mercury state number
    VenusState                 = 2; % Venus state number
    EarthMoonBarycenterState   = 3; % Earth-Moon Barycenter state number
    MarsState                  = 4; % Mars state number
    JupiterState               = 5; % Jupiter state number
    SaturnState                = 6; % Saturn state number
    UranusState                = 7; % Uranus state number
    NeptuneState               = 8; % Neptune state number
    PlutoState                 = 9; % Pluto state number
    MoonState                  = 10; % Geocentric Moon state number
    SunState                   = 11; % Sun state number
    NutationsState             = 12; % nutations state number
    LibrationsState            = 13; % librations state number
    LunarEulerState            = 14; % lunar euler angle rates state number
    TTState                    = 15; % TT-TDB state number
  end
  properties (Constant)
    % list of available ASCII Development Ephemerides
    de_numbers = {'102','200','202','403','405','406','410','413','414','418','421','422','423','424','430','430t','431'};
    % list of minimum years for ASCII Development Ephemerides
    minyears = [-1400, 1600, 1900, 1600, 1600,-3000, 1960, 1900, 1600, 1900, 1900,-3000, 1800,-3000, 1550, 1550,-13000];
    % list of increments between files for ASCII Development Ephemerides
    incyears = [  300,   20,   50,  100,   20,  100,   20,   25,  100,  150,  150,  100,   50,  100,  100,  100,  1000];
    % list of maximum years for ASCII Development Ephemerides
    maxyears = [ 2800, 2160, 2050, 2100, 2200, 2900, 2000, 2025, 2100, 1900, 2050, 2900, 2150, 2900, 2550, 2550, 16000];
  end
  properties (Constant)
    TTLsize = 84; % max character length of the 3 TTL lines in the binary ephemeris
    CNAMsize = 6; % max character length of the Constant Name in the binary ephemeris
    OLDCONMAX = 400; % old maximum for number of constants
    JPLftpSite = 'ftp://ssd.jpl.nasa.gov/pub/eph/planets/ascii'; % site for ASCII files
  end
  properties
    de_number = 405; % DE number of data stored here
    % flag defining physical units of the output states.
    %          = 1, km and km/sec
    %          = 0, au and au/day
    % files are native km and km/sec
    KM = 1;

    SS = [0.0,0.0,0.0]; % JD of first record, JD of last record, days between records
    JPLAU = 0.0; % definition of the astronomical unit (km)
    EM_RATIO = 0.0; % ratio of Earth's mass to Moon's mass
    RECORD_LENGTH = 0; % number of bytes in record
    IPT = zeros(3,Ephem.TTState); % table of indices into record
    scan = []; % linear array of values in records
    indexInc = 0; % increment between records in scan
    indexOffset = 1; % offset for buffer
    numRecs = 0; % number of records
    jd_begin = 0.0; % Julian Date of first record
    jd_end = 0.0; % Julian Date of last record
    jd_inc = 0.0; % days between records
    masses = zeros(Ephem.Mercury,Ephem.EarthMoonBarycenter); % G*mass (AU^3/Day^2) of body number 1-13
    output = 1; % message output
  end
  methods
    function s = Ephem( denum )
      % constructor
      if nargin < 1
        denum = s.de_number;
      end
      if isnumeric(denum)
        denum = sprintf('%d',denum);
      end
      s = s.setDEConstants( denum );
    end
    [ position, velocity, acceleration, error ] = planet_ephemeris( s, tjd, target, center );
    [ target_pos, target_vel, target_acc, error ] = state( s, jed, target );
    [ position, velocity, acceleration ] = interpolate( s, nrl, target, t );
    [ s, error ] = openDEasc( s, denum, fn, dontSave );
    [ s, error ] = openDEeph( s, denum, fn, dontSave );
    [ s, error ] = setDEConstants( s, denum );
    function [ rtn ] = hasNutations( s )
      % check if ephemeris file has nutations
      rtn = (s.IPT(2,Ephem.NutationsState) > 0);
    end
    function [ rtn ] = hasLibrations( s )
      % check if ephemeris file has librations
      [~,nc] = size(s.IPT);
      if nc >= Ephem.LibrationsState && s.IPT(2,Ephem.LibrationsState) > 0
        rtn = true;
      else
        rtn = false;
      end
    end
    function [ rtn ] = hasTTmTDB( s )
      % check if ephemeris file has TT-TDB
      [~,nc] = size(s.IPT);
      if nc >= Ephem.TTState && s.IPT(2,Ephem.TTState) > 0
        rtn = true;
      else
        rtn = false;
      end
    end
    function [ rtn ] = hasLunarEulerAngleRates( s )
      % check if ephemeris file has lunar Euler angle rates
      [~,nc] = size(s.IPT);
      if nc >= Ephem.LunarEulerState && s.IPT(2,Ephem.LunarEulerState) > 0
        rtn = true;
      else
        rtn = false;
      end
    end
  end
  methods (Static)
    function vector = newVector()
      % create a new position (row vector)
      vector = zeros(1,3);
    end
    [ fr1, fr2 ] = splitDouble( tt );
    [ error ] = asc2eph( denum, fn, OUTFILE );
    [ st, error ] = readHeader( denum );
    [ scan, error ] = read_asc( denum, fn, dontSave );
    [ fn ] = asc_name( denum, year );
    [ fn ] = eph_name( denum, year );
    function [ NCOEFF ] = numCoeff( IPT )
      % calculate number of coefficients in record based on IPT matrix
      KMX = 0;
      KHI = 0;
      [~,nc] = size(IPT);
      for I = 1:nc
        if (IPT(1,I) > KMX) && (IPT(2,I) > 0)
          KMX = IPT(1,I);
          KHI = I;
        end
      end
      if KHI == 0
        NCOEFF = 0;
        return;
      end

      if (KHI == Ephem.TTState)
        ND = 1;
      elseif (KHI == Ephem.NutationsState)
        ND = 2;
      else
        ND = 3;
      end

      NCOEFF = IPT(1,KHI)+ND*IPT(2,KHI)*IPT(3,KHI)-1;
    end
    function [rtn] = createAllEphFiles( denum, regen )
      % download and convert all ASCII ephemeris files to binary
      if nargin < 2
        regen = false;
      end
      if isnumeric(denum)
        denum = sprintf('%d',denum);
      end
      index = {};
      for i=1:length(Ephem.de_numbers)
        if strcmp(Ephem.de_numbers{i},denum)
          index = i;
          break;
        end
      end
      if isempty(index)
        rtn = 1;
        return;
      end
      for yr=Ephem.minyears(index):Ephem.incyears(index):Ephem.maxyears(index)
        infile = Ephem.asc_name(denum,yr);
        outfile = Ephem.eph_name(denum,yr);
        if regen || exist(outfile,'file') ~= 2
          rtn = Ephem.asc2eph(denum,infile,outfile);
        end
      end
    end
    rtn = TESTZERODAY( denum, testCase );
    rtn = TESTPOINTS( denum, testCase );
  end
%   methods (Test)
%     function test102(testCase)
%       % test DE102
%       Ephem.TESTZERODAY(102, testCase);
%       Ephem.TESTPOINTS(102, testCase);
%     end
%     function test200(testCase)
%       % test DE200
%       %Ephem.TESTZERODAY(200, testCase);
%       Ephem.TESTPOINTS(200, testCase);
%     end
%     function test405(testCase)
%       % test DE405
%       Ephem.TESTZERODAY(405, testCase);
%       Ephem.TESTPOINTS(405, testCase);
%     end
%     function test430(testCase)
%       % test DE430
%       Ephem.TESTZERODAY(430, testCase);
%       Ephem.TESTPOINTS(430, testCase);
%     end
%     function test430t(testCase)
%       % test DE430t
%       Ephem.TESTZERODAY('430t', testCase);
%       Ephem.TESTPOINTS('430t', testCase);
%     end
%     function test431(testCase)
%       % test DE431
%       Ephem.TESTZERODAY(431, testCase);
%       Ephem.TESTPOINTS(431, testCase);
%     end
%   end

end

author: goGPS-Project


text:
%{
MIT License

Copyright (c) 2018 ZjxRS

Permission is hereby granted, free of charge, to any person obtaining a copy
of this software and associated documentation files (the "Software"), to deal
in the Software without restriction, including without limitation the rights
to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
copies of the Software, and to permit persons to whom the Software is
furnished to do so, subject to the following conditions:

The above copyright notice and this permission notice shall be included in all
copies or substantial portions of the Software.

THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
SOFTWARE.
%}
% This is lifted from https://github.com/ZjxRS/ImageProcess/tree/master/Matlab/ENVI_ReadAndWrite

function D=envidataread(datafile,info)
%ENVIDATAREAD Reads data of ENVI image file according to ENVI image file header information.
%D=ENVIDATAREAD(DATAFILE,INFO) provides image data (D).
%D will be in whatever number format (double, int, etc.) as in the ENVI
%file.

%Original version by Ian Howat, Ohio State Universtiy, ihowat@gmail.com
%Thanks to Yushin Ahn and Ray Jung
%Heavily modified by Felix Totir.
%

if nargin < 2
    error('You must specify the information about data (output of envihdrread)');
end

%-improved- assert(info.header_offset==0,'Non-zero header offset not supported');

%% Set binary format parameters
switch info.byte_order
    case {0}
        machine = 'ieee-le';
    case {1}
        machine = 'ieee-be';
    otherwise
        machine = 'n';
end

iscx=false; %if it is complex
switch info.data_type
    case {1}
        format = 'uint8';
    case {2}
        format= 'int16';
    case{3}
        format= 'int32';
    case {4}
        format= 'float';
    case {5}
        format= 'double';
    case {6}
        iscx=true;
        format= 'float';
    case {9}
        iscx=true;
        format= 'double';
    case {12}
        format= 'uint16';
    case {13}
        format= 'uint32';
    case {14}
        format= 'int64';
    case {15}
        format= 'uint64';
    otherwise
        error(['File type number: ',num2str(dtype),' not supported']);
end

% *** BSQ Format ***
% [Band 1]
% R1: C1, C2, C3, ...
% R2: C1, C2, C3, ...
%  ...
% RN: C1, C2, C3, ...
%
% [Band 2]
%  ...
% [Band N]

% *** BIL Format ***
% [Row 1]
% B1: C1, C2, C3, ...
% B2: C1, C2, C3, ...
%
%  ...
% [Row N]

% *** BIP Format ***
% Row 1
% C1: B1 B2 B3, ...
% C2: B1 B2 B3, ...
% ...
% Row N

n = info.lines*info.samples*info.bands;

fid=fopen(datafile,'r');
fread(fid,info.header_offset,'uint8',0,machine); %we skip the header
if ~iscx
    if strcmp(info.file_type, 'ENVI Classification')
        D = fread(fid,n,'uint8',0,machine);
    else
        D = fread(fid,n,format,0,machine); %alternatively, multibandreader (matlab function) might be used
    end
else
    D = fread(fid,2*n,format,0,machine); %alternatively, multibandreader (matlab function) might be used
    D = complex(D(1:2:end),D(2:2:end));
end
fclose(fid);

switch lower(info.interleave)
    case {'bsq'}
        D = reshape(D,[info.samples,info.lines,info.bands]);
        D = permute(D,[2,1,3]);
    case {'bil'}
        D = reshape(D,[info.samples,info.bands,info.lines]);
        D = permute(D,[3,1,2]);
    case {'bip'}
        D = reshape(D,[info.bands,info.samples,info.lines]);
        D=permute(D,[3,2,1]);
end

author: Bobholamovic

repo: ChangeDetectionToolbox


text:
function plotData(x, y)
%PLOTDATA Plots the data points x and y into a new figure
%   PLOTDATA(x,y) plots the data points and gives the figure axes labels of
%   population and profit.

% ====================== YOUR CODE HERE ======================
% Instructions: Plot the training data into a figure using the
%               "figure" and "plot" commands. Set the axes labels using
%               the "xlabel" and "ylabel" commands. Assume the
%               population and revenue data have been passed in
%               as the x and y arguments of this function.
%
% Hint: You can use the 'rx' option with plot to have the markers
%       appear as red crosses. Furthermore, you can make the
%       markers larger by using plot(..., 'rx', 'MarkerSize', 10);

figure; % open a new figure window

% ============================================================

end

author: zsiciarz

repo: ml-coursera


text:
function m = modes()
% modes of operation
m.OPEN = 1;
m.CLOSE = 2;
m.KSEARCH = 3;
m.PRISEARCH=4;
m.FRSEARCH=5;

author: CHLNDDEV


text:
function C = mrdivide(A, B)
%/    CHEBOP right division.
%   C = A/B, where ones of a A is a CHEBOP and is B a scalar is equivalent
%   to A*(1/B), repectively.
%
%   All other instances returns an error.
%
% See also CHEBOP/MTIMES, CHEBOP/RDIVIDE.

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

if ( isnumeric(B) )
    C = mtimes(A, 1./B);
else
    error('CHEBOP:MRDIVIDE:NotSupported', '%s/%s is not supported.', ...
        class(A), class(B));
end

end

author: chebfun

repo: chebfun

save_path: github-repos/MATLAB/chebfun-chebfun


text:
function [h,a] = PlotDistribution2(var1,var2,varargin)

%PlotDistribution2 - Plot bivariate data along with respective distributions.
%
%  USAGE
%
%    [h,a] = PlotDistribution2(var1,var2,<options>)
%
%    Using cell arrays will overlay variable pairs.
%
%    var1           variable 1
%    var2           variable 2
%    <options>      optional list of property-value pairs (see table below)
%
%    =========================================================================
%     Properties    Values
%    -------------------------------------------------------------------------
%     'nBins'       number of bins for distributions (default = 100)
%     'xbins'       [m M n], lower and upper bounds, and number of bins,
%                   respectively, for x axis (default [min max 100])
%     'ybins'       [m M n], lower and upper bounds, and number of bins,
%                   respectively, for y axis (default [min max 100])
%     'smooth'      standard deviation of Gaussian kernel (default = 5)
%    =========================================================================
%
%  OUTPUT
%
%    h              handles to the plots and histograms (one line per set of
%                   variables)
%    a              axes for the plots and histograms: [main top right]

% Copyright (C) 2004-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.

colors = [1 0.2 0;0 0.4 1;0 1 0.2;1 0.6 0;0.8 0 1;0.8 1 0];
% Default values
nBins = 100;
smooth = 5;
xbins = [];
ybins = [];

% Check number of parameters
if nargin < 2 | mod(length(varargin),2) ~= 0,
  error('Incorrect number of parameters (type ''help PlotDistribution2'' for details).');
end

if isa(var1,'cell'),
	if ~isa(var2,'cell'),
		error('One parameter is a cell array and the other is not (type ''help PlotDistribution2'' for details).');
	end
	if length(var1) ~= length(var2),
		error('Variable lists have different lengths (type ''help PlotDistribution2'' for details).');
	end
else
	var1 = {var1};
	var2 = {var2};
end

% Check parameter sizes
for i = 1:length(var1),
	if ~isdvector(var1{i}) | ~isdvector(var2{i}) | any(size(var1{i})~=size(var2{i})),
		error('Individual variables must be vectors of identical length (type ''help PlotDistribution2'' for details).');
	end
end

% Parse parameter list
for i = 1:2:length(varargin),
	if ~ischar(varargin{i}),
		error(['Parameter ' num2str(i+2) ' is not a property (type ''help PlotDistribution2'' for details).']);
	end
	switch(lower(varargin{i})),
		case 'nbins',
		nBins = varargin{i+1};
		if ~isiscalar(nBins,'>0'),
			error('Incorrect value for property ''nBins'' (type ''help PlotDistribution2'' for details).');
		end
	case 'xbins',
		xbins = varargin{i+1};
		if ~isdvector(xbins,'#3') | xbins(1) > xbins(2) | ~isiscalar(xbins(3),'>0'),
			error('Incorrect value for property ''xbins'' (type ''help PlotDistribution2'' for details).');
		end
	case 'ybins',
		ybins = varargin{i+1};
		if ~isdvector(ybins,'#3') | ybins(1) > ybins(2) | ~isiscalar(ybins(3),'>0'),
			error('Incorrect value for property ''ybins'' (type ''help PlotDistribution2'' for details).');
		end
	case 'smooth',
		smooth = varargin{i+1};
		if ~isdscalar(smooth,'>=0'),
			error('Incorrect value for property ''smooth'' (type ''help PlotDistribution2'' for details).');
		end
	otherwise,
		error(['Unknown property ''' num2str(varargin{i}) ''' (type ''help PlotDistribution2'' for details).']);
  end
end

% Resize current axes to make room for distribution plots
main = gca;
position = get(main,'position');
x = position(1);
h1 = position(3); % current width
hg = h1/10; % horizontal gap
h2 = h1 - hg; % new total width
y = position(2);
v1 = position(4); % current height
vg = v1/10; % vertical gap
v2 = v1 - vg; % new total height
%  set(main,'position',[x+h2/3+hg y+v2/3+vg h2*2/3 v2*2/3]);
set(main,'position',[x y h2*2/3 v2*2/3]);
right = axes('position',[x+h2*2/3+hg y h2/3 v2*2/3]);
top = axes('position',[x y+v2*2/3+vg h2*2/3 v2/3]);

% Axis limits
if isempty(xbins),
	m1 = min(var1{1});
	M1 = max(var1{1});
	for i = 2:length(var1),
		m1 = min([m1;var1{i}]);
		M1 = max([M1;var1{i}]);
	end
	nBinsX = nBins;
else
	m1 = xbins(1);
	M1 = xbins(2);
	nBinsX = xbins(3);
end
if isempty(ybins),
	m2 = min(var2{1});
	M2 = max(var2{1});
	for i = 2:length(var1),
		m2 = min([m2;var2{i}]);
		M2 = max([M2;var2{i}]);
	end
	nBinsY = nBins;
else
	m2 = ybins(1);
	M2 = ybins(2);
	nBinsY = ybins(3);
end
d1 = M1-m1;
m1 = m1-d1/10;
M1 = M1+d1/10;
d2 = M2-m2;
m2 = m2-d2/10;
M2 = M2+d2/10;

h = [];
a = [main top right];

for i = 1:length(var1),
	% Plot bivariate data
	axes(main);hold on;
	p = plot(var1{i},var2{i},'o');
	k = mod(i,size(colors,1)+1);
	set(p,'MarkerFaceColor',colors(k,:),'MarkerEdgeColor',colors(i,:));
	if i > 1,
		m1 = min([m1 pm1]);
		M1 = max([M1 pM1]);
		m2 = min([m2 pm2]);
		M2 = max([M2 pM2]);
	end
	pm1 = m1;
	pM1 = M1;
	pm2 = m2;
	pM2 = M2;
	xlim([m1 M1]);ylim([m2 M2]);

	% Distribution of var1
	axes(top);hold on;
	[dist1,x1] = hist(var1{i},m1:(M1-m1)/nBinsX:M1);
	dist1 = Smooth(dist1/sum(dist1),smooth);
	bv = bar(x1,dist1,1);
	k = mod(i,size(colors,1)+1);
	set(bv,'FaceColor',colors(k,:),'EdgeColor',colors(i,:));
	set(top,'xtick',[]);
	xlim([m1 M1]);

	% Distribution of var2
	axes(right);hold on;
	[dist2,x2] = hist(var2{i},m2:(M2-m2)/nBinsY:M2);
	dist2 = Smooth(dist2/sum(dist2),smooth);
	bh = barh(x2,dist2,1);
	k = mod(i,size(colors,1)+1);
	set(bh,'FaceColor',colors(k,:),'EdgeColor',colors(i,:));
	set(right,'ytick',[]);
	ylim([m2 M2]);

	% store
	h = [h;p bv bh];
end

% Adjust axis locations
set(main,'xaxislocation','bottom','yaxislocation','left');
set(right,'xaxislocation','bottom','yaxislocation','right');
set(top,'xaxislocation','top','yaxislocation','left');
axes(main);

author: buzsakilab

repo: buzcode

save_path: github-repos/MATLAB/buzsakilab-buzcode


text:
function x = fibonacci_lattice_q_nodes ( k )

%*****************************************************************************80
%
%% FIBONACCI_LATTICE_Q_NODES returns Fibonacci lattice nodes in 2D.
%
%  Discussion:
%
%    Because this is a standard lattice rule, it is really only suited
%    for functions which are periodic, of period 1, in both X and Y.
%
%    The number of nodes returned is
%
%      M = fibonacci ( k ).
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    19 November 2008
%
%  Author:
%
%    John Burkardt
%
%  Reference:
%
%    Ian Sloan, Stephen Joe,
%    Lattice Methods for Multiple Integration,
%    Oxford, 1994,
%    ISBN: 0198534728,
%    LC: QA311.S56
%
%  Parameters:
%
%    Input, integer K, the index of the Fibonacci number to be used.
%    K must be at least 3.
%
%    Output, real X(2,M), the nodes.
%
  dim_num = 2;

  m = fibonacci ( k );

  z(1) = 1;
  z(2) = fibonacci ( k - 1 );

  for j = 0 : m - 1
    x(1:dim_num,j+1) = mod ( j * z(1:dim_num) / m, 1.0 );
  end

  return
end

author: johannesgerer

repo: jburkardt-m


text:
function [ bb, bx ] = quadbf ( x, it, in, xc, node )

%*****************************************************************************80
%
%% QUADBF evaluates the quadratic basis functions and derivatives.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    26 April 2006
%
%  Parameters:
%
%    Input, real X, the evaluation point.
%
%    Input, integer IT, the element in which X lies.
%
%    Input, integer IN, the node associated with the basis function.
%
%    Input, real XC(*), the coordinates of the nodes.
%
%    Input, integer NODE(ELEMENT_NUM,3), the nodes that make up
%    each element.
%
%    Output, real BB, BX, the values of the basis function and its
%    derivative with respect to X.
%
  in1 = in;
  in2 = mod ( in, 3 ) + 1;
  in3 = mod ( in + 1, 3 ) + 1;

  i1 = node(it,in1);
  i2 = node(it,in2);
  i3 = node(it,in3);

  den = ( xc(i3) - xc(i1) ) * ( xc(i2) - xc(i1) );

  bb = ( xc(i3) - x ) * ( xc(i2) - x ) / den;

  bx = ( 2.0 * x - xc(i2) - xc(i3) ) / den;

  return
end