text:
function metadata = GetMetaData(filenames,fields)
%GETMETADATA Gets metadata from SAR files for specified fields (when
%available)

metadata = repmat(cell2struct(cell(numel(fields),1),fields),0,0);  % Empty struct

waith = waitbar(0,'Processing Metadata');

for i=1:length(filenames)

    Temp = [];

    a = sprintf('Processing file %d of %d',i,length(filenames));
    waitbar(double(i)/double(length(filenames)),waith,a);
    PHD = 0;
    try
        reader_obj = open_reader(filenames{i});
        if iscell(reader_obj)
            reader_obj = reader_obj{1};
        end
        meta = reader_obj.get_meta();
        SamplePixel = reader_obj.read_chip([1 1],[1 1]);
        reader_obj.close();
    catch
        try
            reader_obj = open_ph_reader(filenames{i});
            if iscell(reader_obj)
                reader_obj = reader_obj{1};
            end
            meta = reader_obj.get_meta();
            PHD = 1;
            if ~isfield(meta,'SCPCOA')
                % Get info for a few specific pulses
                if isfield(reader_obj,'read_raw')
                    read_fun = reader_obj.read_raw;
                elseif isfield(reader_obj,'read_cphd')
                    read_fun = reader_obj.read_cphd;
                end
                [~, vbmeta] = read_fun([1 ...
                    round(meta.Data.Channel(1).NumVectors/2) ...
                    round(meta.Data.Channel(1).NumVectors/2) + 1 ...
                    meta.Data.Channel(1).NumVectors],[]); %first, last and two center pulses (to get velocity)
            end
            meta = derived_phd_fields(meta,vbmeta);
        catch
            continue;
        end
    end

    %loop through fieldnames and populate as metadat supports
    for j=1:length(fields)
        switch fields{j}
            case 'Corename'
                if isfield(meta,'CollectionInfo') && isfield(meta.CollectionInfo,'CoreName')
                    Temp.Corename = meta.CollectionInfo.CoreName;
                elseif isfield(meta,'CollectionID') && isfield(meta.CollectionID,'CoreName')
                    Temp.Corename = meta.CollectionID.CoreName;
                else
                    Temp.Corename = '';
                end
            case 'DataType'  %todo: if PHD then set to PHD
                if PHD
                    Temp.DataType = 'PHD';
                else
                    if isreal(SamplePixel)
                        Temp.DataType = 'Detected';
                    else
                        Temp.DataType = 'Complex';
                    end
                end
            case 'PixelType'
                if isfield(meta,'ImageData') && isfield(meta.ImageData,'PixelType')
                    Temp.PixelType = meta.ImageData.PixelType;
                else
                    Temp.PixelType = '';
                end
            case 'Lat'
                if isfield(meta,'GeoData') && isfield(meta.GeoData,'SCP') && ...
                   isfield(meta.GeoData.SCP,'LLH') && isfield(meta.GeoData.SCP.LLH,'Lat')
                    Temp.Lat = meta.GeoData.SCP.LLH.Lat;
                else
                    Temp.Lat = '';
                end
            case 'Lon'
                if isfield(meta,'GeoData') && isfield(meta.GeoData,'SCP') && ...
                   isfield(meta.GeoData.SCP,'LLH') && isfield(meta.GeoData.SCP.LLH,'Lon')
                    Temp.Lon = meta.GeoData.SCP.LLH.Lon;
                else
                    Temp.Lon = '';
                end
            case 'CollectTime'
                if isfield(meta,'Timeline') && isfield(meta.Timeline,'CollectStart')
                    Temp.CollectTime = datestr(meta.Timeline.CollectStart);
                elseif isfield(meta,'Global') && isfield(meta.Global,'CollectStart')
                    Temp.CollectTime = datestr(meta.Global.CollectStart);
                elseif isfield(meta,'Global') && isfield(meta.Global,'Timeline') && ...
                        isfield(meta.Global.Timeline,'CollectionStart')
                    Temp.CollectTime = datestr(meta.Global.Timeline.CollectionStart);
                else
                    Temp.CollectTime = '';
                end
            case 'Azimuth'
                if isfield(meta,'SCPCOA') && isfield(meta.SCPCOA,'AzimAng')
                    Temp.Azimuth = meta.SCPCOA.AzimAng;
                else
                    Temp.Azimuth = '';
                end
            case 'Graze'
                if isfield(meta,'SCPCOA') && isfield(meta.SCPCOA,'GrazeAng')
                    Temp.Graze = meta.SCPCOA.GrazeAng;
                else
                    Temp.Graze = '';
                end
            case 'NumRows'
                if isfield(meta,'ImageData') && isfield(meta.ImageData,'NumRows')
                    Temp.NumRows = meta.ImageData.NumRows;
                else
                    Temp.NumRows = '';
                end
            case 'NumColumns'
                if isfield(meta,'ImageData') && isfield(meta.ImageData,'NumCols')
                    Temp.NumColumns = meta.ImageData.NumCols;
                else
                    Temp.NumColumns = '';
                end
            case 'AzSS'
                if isfield(meta,'Grid') && isfield(meta.Grid,'Col') && ...
                   isfield(meta.Grid.Col,'SS')
                    Temp.AzSS = meta.Grid.Col.SS;
                else
                    Temp.AzSS = '';
                end
            case 'RgSS'
                if isfield(meta,'Grid') && isfield(meta.Grid,'Row') && ...
                   isfield(meta.Grid.Row,'SS')
                    Temp.RgSS = meta.Grid.Row.SS;
                else
                    Temp.RgSS = '';
                end
            case 'Duration' %todo: check for TStart and TEnd
                if isfield(meta,'Timeline') && isfield(meta.Timeline,'CollectDuration')
                    Temp.Duration = meta.Timeline.CollectDuration;
                else
                    Temp.Duration = '';
                end
            case 'RadarMode'
                if isfield(meta,'CollectionInfo') && isfield(meta.CollectionInfo,'RadarMode') && ...
                   isfield(meta.CollectionInfo.RadarMode,'ModeType')
                    Temp.RadarMode = meta.CollectionInfo.RadarMode.ModeType;
                elseif isfield(meta,'CollectionID') && isfield(meta.CollectionID,'RadarMode') && ...
                   isfield(meta.CollectionID.RadarMode,'ModeType')
                    Temp.RadarMode = meta.CollectionID.RadarMode.ModeType;
                else
                    Temp.RadarMode = '';
                end
            case 'Polarization'
                if isfield(meta,'ImageFormation') && isfield(meta.ImageFormation,'TxRcvPolarizationProc')
                    Temp.Polarization = meta.ImageFormation.TxRcvPolarizationProc;
                elseif isfield(meta,'RadarCollection') && isfield(meta.RadarCollection,'RcvChannels') && ...
                       isfield(meta.RadarCollection.RcvChannels,'ChanParameters') && ...
                       isfield(meta.RadarCollection.RcvChannels.ChanParameters(1),'TxRcvPolarization')
                    Temp.Polarization = meta.RadarCollection.RcvChannels.ChanParameters(1).TxRcvPolarization;
                elseif isfield(meta,'Channel') && isfield(meta.Channel,'Parameters') && ...
                        isfield(meta.Channel.Parameters,'Polarization') && ...
                        all(isfield(meta.Channel.Parameters(1).Polarization,{'TxPol','RcvPol'}))
                    Temp.Polarization = [meta.Channel.Parameters(1).Polarization.TxPol ...
                        ':' meta.Channel.Parameters(1).Polarization.RcvPol];
                elseif isfield(meta,'Channel') && isfield(meta.Channel,'Parameters') && ...
                        isfield(meta.Channel.Parameters,'RcvPol')
                    Temp.Polarization = meta.Channel.Parameters(1).RcvPol;
                else
                    Temp.Polarization = '';
                end
            case 'ImagePlane'
                if isfield(meta,'Grid') && isfield(meta.Grid,'ImagePlane')
                    Temp.ImagePlane = meta.Grid.ImagePlane;
                else
                    Temp.ImagePlane = '';
                end
            case 'HAE'
                if isfield(meta,'GeoData') && isfield(meta.GeoData,'SCP') && ...
                   isfield(meta.GeoData.SCP,'LLH') && isfield(meta.GeoData.SCP.LLH,'HAE')
                    Temp.HAE = meta.GeoData.SCP.LLH.HAE;
                else
                    Temp.HAE = '';
                end
            case 'SCPTime'
                if isfield(meta,'SCPCOA') && isfield(meta.SCPCOA,'SCPTime')
                    Temp.SCPTime = meta.SCPCOA.SCPTime;
                else
                    Temp.SCPTime = '';
                end
            case 'Side'
                if isfield(meta,'SCPCOA') && isfield(meta.SCPCOA,'SideOfTrack')
                    Temp.Side = meta.SCPCOA.SideOfTrack;
                else
                    Temp.Side = '';
                end
            case 'SlantRange'
                if isfield(meta,'SCPCOA') && isfield(meta.SCPCOA,'SlantRange')
                    Temp.SlantRange = meta.SCPCOA.SlantRange;
                else
                    Temp.SlantRange = '';
                end
            case 'GroundRange'
                if isfield(meta,'SCPCOA') && isfield(meta.SCPCOA,'GroundRange')
                    Temp.GroundRange = meta.SCPCOA.GroundRange;
                else
                    Temp.GroundRange = '';
                end
            case 'DCA'
                if isfield(meta,'SCPCOA') && isfield(meta.SCPCOA,'DopplerConeAng')
                    Temp.DCA = meta.SCPCOA.DopplerConeAng;
                else
                    Temp.DCA = '';
                end
            case 'Twist'
                if isfield(meta,'SCPCOA') && isfield(meta.SCPCOA,'TwistAng')
                    Temp.Twist = meta.SCPCOA.TwistAng;
                else
                    Temp.Twist = '';
                end
            case 'Slope'
                if isfield(meta,'SCPCOA') && isfield(meta.SCPCOA,'SlopeAng')
                    Temp.Slope = meta.SCPCOA.SlopeAng;
                else
                    Temp.Slope = '';
                end
            case 'Layover'
                if isfield(meta,'SCPCOA') && isfield(meta.SCPCOA,'LayoverAng')
                    Temp.Layover = meta.SCPCOA.LayoverAng;
                else
                    Temp.Layover = '';
                end
            case 'Filename'
                Temp.Filename = filenames{i};
            case 'ROV'
                if isfield(meta,'SCPCOA') && isfield(meta.SCPCOA,'SlantRange') && ...
                   isfield(meta.SCPCOA,'ARPVel')
                   Vel = norm([meta.SCPCOA.ARPVel.X meta.SCPCOA.ARPVel.Y meta.SCPCOA.ARPVel.Z]);
                   Temp.ROV = meta.SCPCOA.SlantRange./Vel;
                else
                   Temp.ROV = '';
                end
            case 'PRF'
                if isfield(meta,'Timeline') && isfield(meta.Timeline,'IPP') && isfield(meta.Timeline.IPP,'Set') && ...
                   isfield(meta.Timeline.IPP.Set,'IPPStart') && isfield(meta.Timeline.IPP.Set,'IPPEnd')  && ...
                   isfield(meta.Timeline.IPP.Set,'TStart') && isfield(meta.Timeline.IPP.Set,'TEnd')
                   Temp.PRF = (meta.Timeline.IPP.Set.IPPEnd-meta.Timeline.IPP.Set.IPPStart)./ ...
                                          (meta.Timeline.IPP.Set.TEnd-meta.Timeline.IPP.Set.TStart);
                else
                   Temp.PRF = '';
                end
            case 'SARAlt'
                if isfield(meta,'SCPCOA') && isfield(meta.SCPCOA,'ARPPos')
                    R = norm([meta.SCPCOA.ARPPos.X meta.SCPCOA.ARPPos.Y meta.SCPCOA.ARPPos.Z]);
                    Temp.SARAlt = (R-6378135)./1000; %km
                elseif isfield(meta,'Position') && isfield(meta.Position,'ARPPoly')
                    R = norm([meta.Position.ARPPoly.X(1) meta.Position.ARPPoly.Y(1) meta.Position.ARPPoly.Z(1)]);
                    Temp.SARAlt = (R-6378135)./1000; %km
                else
                    Temp.SARAlt = '';
                end
            case 'CenterFrequency'
                if isfield(meta,'RadarCollection') && isfield(meta.RadarCollection,'TxFrequency') && ...
                   all(isfield(meta.RadarCollection.TxFrequency,{'Min','Max'}))
                    Temp.CenterFrequency = (meta.RadarCollection.TxFrequency.Max + ...
                                            meta.RadarCollection.TxFrequency.Min)/2;
                else
                    Temp.CenterFrequency = '';
                end
            case 'TxBandwidth'
                if isfield(meta,'RadarCollection') && isfield(meta.RadarCollection,'TxFrequency') && ...
                   all(isfield(meta.RadarCollection.TxFrequency,{'Min','Max'}))
                    Temp.TxBandwidth = meta.RadarCollection.TxFrequency.Max - ...
                                            meta.RadarCollection.TxFrequency.Min;
                else
                    Temp.TxBandwidth = '';
                end
            case 'TxPulseLength'
                if isfield(meta,'RadarCollection') && isfield(meta.RadarCollection,'Waveform') && ...
                   isfield(meta.RadarCollection.Waveform,'WFParameters') && ...
                   isfield(meta.RadarCollection.Waveform.WFParameters,'TxPulseLength')
                    Temp.TxPulseLength = meta.RadarCollection.Waveform.WFParameters(1).TxPulseLength;
                else
                    Temp.TxPulseLength = '';
                end
            case 'TxChirpRate'
                if isfield(meta,'RadarCollection') && isfield(meta.RadarCollection,'Waveform') && ...
                   isfield(meta.RadarCollection.Waveform,'WFParameters') && ...
                   isfield(meta.RadarCollection.Waveform.WFParameters,'TxFMRate')
                    Temp.TxChirpRate = meta.RadarCollection.Waveform.WFParameters(1).TxFMRate;
                else
                    Temp.TxChirpRate = '';
                end
            case 'SampleRate'
                if isfield(meta,'RadarCollection') && isfield(meta.RadarCollection,'Waveform') && ...
                   isfield(meta.RadarCollection.Waveform,'WFParameters') && ...
                   isfield(meta.RadarCollection.Waveform.WFParameters,'ADCSampleRate')
                    Temp.SampleRate = meta.RadarCollection.Waveform.WFParameters(1).ADCSampleRate;
                else
                    Temp.SampleRate = '';
                end
            case 'RxWindowLength'
                if isfield(meta,'RadarCollection') && isfield(meta.RadarCollection,'Waveform') && ...
                   isfield(meta.RadarCollection.Waveform,'WFParameters') && ...
                   isfield(meta.RadarCollection.Waveform.WFParameters,'RcvWindowLength')
                    Temp.RxWindowLength = meta.RadarCollection.Waveform.WFParameters(1).RcvWindowLength;
                else
                    Temp.RxWindowLength = '';
                end
            case 'NumPulses'
                if isfield(meta,'Data') && isfield(meta.Data,'Channel') && ...
                   isfield(meta.Data.Channel,'NumVectors')
                    Temp.NumPulses = meta.Data.Channel.NumVectors;
                elseif isfield(meta,'Timeline') && isfield(meta.Timeline,'IPP') && ...
                        isfield(meta.Timeline.IPP,'Set') && ...
                        isfield(meta.Timeline.IPP.Set,'IPPEnd')
                    Temp.NumPulses = max([meta.Timeline.IPP.Set.IPPEnd]);
                else
                    Temp.NumPulses = '';
                end
            case 'NumSamples'
                if isfield(meta,'Data') && isfield(meta.Data,'Channel') && ...
                   isfield(meta.Data.Channel,'NumSamples')
                    Temp.NumSamples = meta.Data.Channel.NumSamples;
                elseif isfield(meta,'RadarCollection') && isfield(meta.RadarCollection,'Waveform') && ...
                   isfield(meta.RadarCollection.Waveform,'WFParameters') && ...
                   all(isfield(meta.RadarCollection.Waveform.WFParameters,{'ADCSampleRate','RcvWindowLength'}))
                    Temp.NumSamples = meta.RadarCollection.Waveform.WFParameters(1).ADCSampleRate * ...
                       meta.RadarCollection.Waveform.WFParameters(1).RcvWindowLength;
                else
                    Temp.NumSamples = '';
                end
            otherwise
                Temp.(fields{j}) = '';
        end
    end
    metadata(end+1) = Temp;
end

close(waith);

end

author: ngageoint


text:
function [fBValue, fAValue, fStdDev, fMc, fMeanMag] = plot_FMD(mCatalog, varargin)
        % Creates and plots a frequency magnitude distribution including the b-value
    %  [fBValue, fAValue, fStdDev, fMc, fMeanMag] = PLOT_FMD(mCatalog,showCumulative, hAxes, sSymbol, sColor, bPlotB, nCalculateMC, binWidth)
    % -------------------------------------------------------------------------------------------

    %
    % PLOT_FMD(mCatalog) opens a figure and plots the frequency magnitude distribution
    %   with standard parameters
    %
    % Input parameters:
    %   mCatalog        Earthquake catalog
    %
    % Named Input Parameters:  use as PLOT_FMD(... , Name, value)
    %   ShowCumulative  Plot cumulative frequency magnitude distribution (true) or non-cumulative (false)
    %   Axes           Handle of axes to plot the frequency magnitude distribution
    %   Symbol         Type of symbol for plotting (refer to Matlab 'plot')
    %   Color          Color of plot (refer to Matlab 'plot')
    %   ShowBval          Also plot the b-value line with markers
    %   CalcMethod    Method to determine the magnitude of completeness
    %                   1: Maximum curvature
    %                   2: Fixed Mc = minimum magnitude (Mmin)
    %                   3: Mc90 (90% probability)
    %                   4: Mc95 (95% probability)
    %                   5: Best combination (Mc95 - Mc90 - maximum curvature)
    %   BinWidth        Magnitude binning of the catalog (default 0.1)
    %
    % Output parameters:
    %   fBValue         Calculated b-value
    %   fAValue         Calculateda-value
    %   fStdDev         Standard deviation of b-value
    %   fMc             Magnitude of completeness
    %   fMeanMag        Determined mean magnitude
    %
    % Danijel Schorlemmer
    % June 16, 2003

    report_this_filefun();

    p = inputParser();
    p.addRequired('mCatalog');
    p.addParameter('ShowCumulative', true);
    p.addParameter('Axes',          []);
    P.addParameter('Symbol',        's');
    p.addParameter('Color',         'k');
    p.addParameter('ShowBval',      true)
    p.addParameter('CalcMethod',    2);
    p.addParameter('BinWidth',      0.1);
    p.addParameter('MarkerSize',    12);
    p.KeepUnmatched = true;
    p.parse(varargin{:});

    if isempty(p.Results.Axes)
        figure;
        hAxes=newplot;
    else
        hAxes = p.Results.Axes;
    end

    binWidth = p.Results.BinWidth;
    sColor = p.Results.Color;
    % Activate given axes
    axes(hAxes);

    % Create the frequency magnitude distribution vector
    [nEvsCum, nEvsNonCum, magX] = calc_FMD(mCatalog.Magnitude);

    % Plot the frequency magnitude distribution
    if p.Results.showCumulative
        hPlot = semilogy(magX, nEvsCum);
    else
        hPlot = semilogy(magX, nEvsNonCum);
    end

    hPlot.Symbol = p.Results.Symbol;
    hPlot.Color = p.Results.Color;
    hPlot.MarkerSize = p.Results.MarkerSize;

    if p.Results.ShowBval
        % Add further plots to the axes
        hAxes.NextPlot = 'add';

        % Calculate magnitude of completeness
        fMc = calc_Mc(mCatalog, p.Results.CalcMethod, binWidth);

        % Determine the positions of 'x'-markers
        nIndexLo = find((magX < fMc + 0.05) & (magX > fMc - 0.05));

        % Plot the 'x'-marker
        hPlot = semilogy(magX(nIndexLo), nEvsCum(nIndexLo));
        hPlot.Marker = 'x';
        hPlot.Color = sColor;
        hPlot.LineWidth = 2.5;
        hPlot.MarkerSize = 12;

        hPlot = semilogy(magX(1), nEvsCum(1));
        hPlot.Marker = 'x';
        hPlot.Color = sColor;
        hPlot.LineWidth = 2.5;
        hPlot.MarkerSize = 12;

        % Calculate the b-value etc. for M > Mc
        vSel = mCatalog.Magnitude >= fMc-(binWidth/2);
        fMeanMag = mean(mCatalog.Magnitude(vSel));
        [ fBValue, fStdDev, fAValue] =  calc_bmemag(mCatalog.Magnitude(vSel), binWidth);

        % Plot the line representing the b-value
        vPoly = [-1*fBValue fAValue];

        fMagHi = magX(1);
        vMagnitudes = magX(fMc - 0.0001 <= magX & magX <= fMagHi);

        fBFunc = 10.^(polyval(vPoly, vMagnitudes));
        hPlot = semilogy(vMagnitudes, fBFunc, sColor);
        hPlot.LineWidth = 2.0;
    end
end

author: CelsoReyes


text:
function fig_hdl = img2curve
% img2curve
%-------------------------------------------------------------------------------
% File name   : img2curve.m
% Generated on: 11-Mar-2013 10:54:18
% Description :
%
% img2curve helps import data from a graph in a photo file, received from a
% scan, catalog, pdf etc. by marking points on the curve.
% Output will include the marked points and a polynomial fit.
% Work flow:
% 1. Load image file.
% 2. Mark scale size for X and Y axis.
% 3. Set a value for a point on the image (0,0 or other) to locate CS.
% 4. Mark curve points.
% 5. Export to mat file, excel or work space.
%
% Good Luck,
% Amihay Blau

%-------------------------------------------------------------------------------

% Initialize handles structure
handles = struct();

% Create all UI controls
build_gui();

% Assign function output
fig_hdl = handles.img2curve_figure;

%% ---------------------------------------------------------------------------
	function build_gui()
% Creation of all uicontrols

		% --- FIGURE -------------------------------------
		handles.img2curve_figure = figure( ...
			'Tag', 'img2curve_figure', ...
			'Units', 'centimeters', ...
			'Position', [4.9968253125 8.46023333333333 8.38091864583333 12.2673383333333], ...
			'Name', 'img2curve', ...
			'MenuBar', 'none', ...
			'NumberTitle', 'off', ...
			'Color', [0.804 0.878 0.969], ...
			'Resize', 'on');

		% --- PANELS -------------------------------------
		handles.import_uipanel = uipanel( ...
			'Parent', handles.img2curve_figure, ...
			'Tag', 'import_uipanel', ...
			'UserData', zeros(1,0), ...
			'Units', 'centimeters', ...
			'Position', [0.264382291666667 7.34982770833334 8.0107834375 4.70600479166667], ...
			'FontSize', 12, ...
			'BackgroundColor', [0.839 0.91 0.851], ...
			'Title', 'Import Image');

		handles.curve_uipanel = uipanel( ...
			'Parent', handles.img2curve_figure, ...
			'Tag', 'curve_uipanel', ...
			'UserData', zeros(1,0), ...
			'Units', 'centimeters', ...
			'Position', [0.264382291666667 3.7277903125 8.0107834375 3.48984625], ...
			'FontSize', 12, ...
			'BackgroundColor', [0.839 0.91 0.851], ...
			'Title', 'Define Curve');

		handles.save_uipanel = uipanel( ...
			'Parent', handles.img2curve_figure, ...
			'Tag', 'save_uipanel', ...
			'UserData', zeros(1,0), ...
			'Units', 'centimeters', ...
			'Position', [0.264382291666667 0.211505833333333 8.0107834375 3.35765510416667], ...
			'FontSize', 12, ...
			'BackgroundColor', [0.839 0.91 0.851], ...
			'Title', 'Save Curve Data');

		% --- STATIC TEXTS -------------------------------------
		handles.import_xScale_text = uicontrol( ...
			'Parent', handles.import_uipanel, ...
			'Tag', 'import_xScale_text', ...
			'Style', 'text', ...
			'Units', 'normalized', ...
			'Position', [0.381270903010033 0.33116883116883 0.535117056856187 0.246753246753247], ...
			'FontSize', 10, ...
			'BackgroundColor', [0.839 0.91 0.851], ...
			'String', {'X Scale','Y Scale'}, ...
			'HorizontalAlignment', 'left');

		handles.curve_PolOrder_text = uicontrol( ...
			'Parent', handles.curve_uipanel, ...
			'Tag', 'curve_PolOrder_text', ...
			'Style', 'text', ...
			'Units', 'normalized', ...
			'Position', [0.468227424749164 0.722222222222222 0.337792642140468 0.175925925925926], ...
			'FontSize', 10, ...
			'BackgroundColor', [0.839 0.91 0.851], ...
			'String', 'Polynom Order');

		handles.curve_Pol_text = uicontrol( ...
			'Parent', handles.curve_uipanel, ...
			'Tag', 'curve_Pol_text', ...
			'Style', 'text', ...
			'Units', 'normalized', ...
			'Position', [0.0602006688963211 0.0740740740740741 0.82943143812709 0.527777777777778], ...
			'FontSize', 10, ...
			'BackgroundColor', [0.839 0.91 0.851], ...
			'String', 'Curve Equation');

		% --- PUSHBUTTONS -------------------------------------
		handles.import_loadImg_pushbutton = uicontrol( ...
			'Parent', handles.import_uipanel, ...
			'Tag', 'import_loadImg_pushbutton', ...
			'UserData', zeros(1,0), ...
			'Style', 'pushbutton', ...
			'Units', 'normalized', ...
			'Position', [0.0735785953177258 0.669194186256364 0.247491638795987 0.25974025974026], ...
			'FontSize', 10, ...
			'String', 'Load Image', ...
			'TooltipString', 'Load curve image', ...
			'CData', zeros(1,0), ...
			'Callback', @import_loadImg_pushbutton_Callback);

		handles.import_getScaleFile_pushbutton = uicontrol( ...
			'Parent', handles.import_uipanel, ...
			'Tag', 'import_getScaleFile_pushbutton', ...
			'Style', 'pushbutton', ...
			'Units', 'normalized', ...
			'Position', [0.377926421404682 0.668831168831169 0.324414715719064 0.25974025974026], ...
			'FontSize', 10, ...
			'String', 'Get Scale File', ...
			'TooltipString', 'Get a previous scale file', ...
			'Enable', 'off', ...
			'Callback', @import_getScaleFile_pushbutton_Callback);

		handles.import_scale_pushbutton = uicontrol( ...
			'Parent', handles.import_uipanel, ...
			'Tag', 'import_scale_pushbutton', ...
			'Style', 'pushbutton', ...
			'Units', 'normalized', ...
			'Position', [0.0735785953177258 0.363999381061558 0.240802675585284 0.201298701298701], ...
			'FontSize', 10, ...
			'String', 'Set Scale', ...
			'TooltipString', 'Set image scale by marking points on axes and writing the difference between these points', ...
			'Enable', 'off', ...
			'Callback', @import_scale_pushbutton_Callback);

		handles.import_MarkRef_pushbutton = uicontrol( ...
			'Parent', handles.import_uipanel, ...
			'Tag', 'import_MarkRef_pushbutton', ...
			'Style', 'pushbutton', ...
			'Units', 'normalized', ...
			'Position', [0.0769230769230769 0.0977656148277906 0.471571906354515 0.201298701298701], ...
			'FontSize', 10, ...
			'String', 'Mark Reference Point', ...
			'TooltipString', 'set the C.S. by marking a point in the image. point value is set on the right', ...
			'Enable', 'off', ...
			'Callback', @import_MarkRef_pushbutton_Callback);

		handles.curve_markPoints_pushbutton = uicontrol( ...
			'Parent', handles.curve_uipanel, ...
			'Tag', 'curve_markPoints_pushbutton', ...
			'Style', 'pushbutton', ...
			'Units', 'normalized', ...
			'Position', [0.0635451505016723 0.675925925925926 0.404682274247492 0.305555555555555], ...
			'FontSize', 10, ...
			'String', 'Mark Curve Points', ...
			'TooltipString', 'Mark points on the curve', ...
			'Enable', 'off', ...
			'Callback', @curve_markPoints_pushbutton_Callback);

		handles.save_pushbutton = uicontrol( ...
			'Parent', handles.save_uipanel, ...
			'Tag', 'save_pushbutton', ...
			'Style', 'pushbutton', ...
			'Units', 'centimeters', ...
			'Position', [5.49915166666667 2.08862010416666 2.1414965625 0.846023333333334], ...
			'String', 'Save to mat', ...
			'TooltipString', 'Save results', ...
			'Enable', 'off', ...
			'Callback', @save_pushbutton_Callback);

		handles.save_excel_pushbutton = uicontrol( ...
			'Parent', handles.save_uipanel, ...
			'Tag', 'save_excel_pushbutton', ...
			'Style', 'pushbutton', ...
			'Units', 'centimeters', ...
			'Position', [5.49915166666667 1.16328208333333 2.1414965625 0.846023333333334], ...
			'String', 'Save to excel', ...
			'TooltipString', 'Save results', ...
			'Enable', 'off', ...
			'Callback', @save_excel_pushbutton_Callback);

		handles.save_ws_pushbutton = uicontrol( ...
			'Parent', handles.save_uipanel, ...
			'Tag', 'save_ws_pushbutton', ...
			'Style', 'pushbutton', ...
			'Units', 'centimeters', ...
			'Position', [5.49915166666667 0.237944062499997 2.1414965625 0.846023333333334], ...
			'String', 'Save to WS', ...
			'TooltipString', 'Save results', ...
			'Enable', 'off', ...
			'Callback', @save_ws_pushbutton_Callback);

		% --- EDIT TEXTS -------------------------------------
		handles.import_MarkRefX_edit = uicontrol( ...
			'Parent', handles.import_uipanel, ...
			'Tag', 'import_MarkRefX_edit', ...
			'Style', 'edit', ...
			'Units', 'normalized', ...
			'Position', [0.648829431438127 0.123739640801817 0.0869565217391304 0.155844155844156], ...
			'FontSize', 10, ...
			'BackgroundColor', [1 1 1], ...
			'String', '0', ...
			'TooltipString', 'Reference point X value', ...
			'Callback', @import_MarkRefX_edit_Callback, ...
			'CreateFcn', @import_MarkRefX_edit_CreateFcn);

		handles.import_MarkRefY_edit = uicontrol( ...
			'Parent', handles.import_uipanel, ...
			'Tag', 'import_MarkRefY_edit', ...
			'Style', 'edit', ...
			'Units', 'normalized', ...
			'Position', [0.749163879598662 0.123739640801817 0.0869565217391304 0.155844155844156], ...
			'FontSize', 10, ...
			'BackgroundColor', [1 1 1], ...
			'String', '0', ...
			'TooltipString', 'Reference point Y value', ...
			'Callback', @import_MarkRefY_edit_Callback, ...
			'CreateFcn', @import_MarkRefY_edit_CreateFcn);

		handles.curve_PolOrder_edit = uicontrol( ...
			'Parent', handles.curve_uipanel, ...
			'Tag', 'curve_PolOrder_edit', ...
			'Style', 'edit', ...
			'Units', 'normalized', ...
			'Position', [0.832775919732442 0.703703703703704 0.130434782608696 0.203703703703704], ...
			'FontSize', 10, ...
			'BackgroundColor', [1 1 1], ...
			'String', '3', ...
			'TooltipString', 'Polynom order', ...
			'Callback', @curve_PolOrder_edit_Callback, ...
			'CreateFcn', @curve_PolOrder_edit_CreateFcn);

		handles.save_fName_edit = uicontrol( ...
			'Parent', handles.save_uipanel, ...
			'Tag', 'save_fName_edit', ...
			'Style', 'edit', ...
			'Units', 'centimeters', ...
			'Position', [0.502326354166667 1.00465270833333 4.78531947916666 1.08396739583333], ...
			'FontSize', 10, ...
			'BackgroundColor', [1 1 1], ...
			'String', 'Curve', ...
			'TooltipString', 'Save file name', ...
			'Callback', @save_fName_edit_Callback, ...
			'CreateFcn', @save_fName_edit_CreateFcn);

	end

%% ---------------------------------------------------------------------------
	function import_loadImg_pushbutton_Callback(hObject,evendata) %#ok<INUSD>

        [picName, dirNm] = uigetfile('*');
        if picName == 0
            return
        end

        if ~strcmp(dirNm, cd)
            picName = [dirNm picName];
        end
        img = imread(picName);

        handles.disp_fig = figure;

        imagesc([], [], flipdim(img,1));
        set(gca,'ydir','normal', 'XTick', [], 'YTick', []);

        % set(handles.disp_fig, 'Pointer', 'fullcrosshair')

        handles.disp_ax = gca;

        setappdata(handles.img2curve_figure, 'img', img);

        pushbutton_activate('import');

	end

%% ---------------------------------------------------------------------------
	function import_getScaleFile_pushbutton_Callback(hObject,evendata) %#ok<INUSD>

        [fName, dirNm] = uigetfile('scale.mat');
        if fName == 0
            return
        end

        S = load([dirNm fName]);
        try
            setappdata(handles.img2curve_figure, 'scale', S.scale);
            setappdata(handles.img2curve_figure, 'oo', S.oo);
            setappdata(handles.img2curve_figure, 'ooVal', S.ooVal);

            plotImBackground;

            pushbutton_activate;
        catch
            errordlg('Previous scale file loading failed!', 'Load error');
        end

	end

%% ---------------------------------------------------------------------------
    function import_scale_pushbutton_Callback(hObject,evendata) %#ok<INUSD>

        pushbutton_deactivate;

        scale = getAppD('scale', [1 1]);

        scale(1) = getAxisScale(1, scale(1));
        scale(2) = getAxisScale(2, scale(2));
        setappdata(handles.img2curve_figure, 'scale', scale);

        set(handles.import_xScale_text, 'String', {...
            ['X Scale 1:' num2str(scale(1), 3)],...
            ['Y Scale 1:' num2str(scale(2), 3)]});

        checkSaveScale;

        pushbutton_activate;

    end

%% ---------------------------------------------------------------------------
    function import_MarkRef_pushbutton_Callback(hObject,evendata) %#ok<INUSD>

        pushbutton_deactivate;

        scale = getAppD('scale', [1  1]);
        oo = getAppD('oo',[]);
        ooVal = getAppD('ooVal', [0  0]);

        axes(handles.disp_ax);
        t = title(handles.disp_ax, ['Mark (' num2str(ooVal) ') point'],...
            'FontSize', 16, 'Color', 'r');

        if isempty(oo)
            oo = ginput(1);
        else
            oo = ( ginput(1) ) ./scale + oo - ooVal ./ scale;
        end

        delete(t);

        pushbutton_activate;

        setappdata(handles.img2curve_figure, 'oo', oo);
        setOoVal;

        checkSaveScale;

    end

%% ---------------------------------------------------------------------------
    function curve_markPoints_pushbutton_Callback(hObject,evendata) %#ok<INUSD>

        scale = [];
        oo = [];
        ooVal = [];
        names = {'scale' 'oo' 'ooVal'};
        for i = 1:length(names)
            if isappdata(handles.img2curve_figure, names{i})
                eval([names{i} ' = getappdata(handles.img2curve_figure, ''' names{i} ''');']);
            else
                errordlg('Set scale and oo first!', 'Error');
                return
            end
        end

        pushbutton_deactivate;

        axes(handles.disp_ax);

        tit = title(handles.disp_ax, {'Mark as many curve points as you want' 'Press Enter when done'},...
            'FontSize', 16, 'Color', 'r');
        po = ginput;
        if isempty(po)
            try
                load('curvePoints.mat');
            catch %#ok<*CTCH>
                return
            end
        else
            save('curvePoints.mat', 'po');
        end

        setappdata(handles.img2curve_figure, 'po', po);

        calcPol;

        pushbutton_activate;

    end

%% ---------------------------------------------------------------------------
    function save_pushbutton_Callback(hObject,evendata) %#ok<INUSD>

        po = getappdata(handles.img2curve_figure, 'po');
        oo = getappdata(handles.img2curve_figure, 'oo');
        ooVal = getappdata(handles.img2curve_figure, 'ooVal');
        scale = getappdata(handles.img2curve_figure, 'scale');
        pol = getappdata(handles.img2curve_figure, 'pol');
        img = getappdata(handles.img2curve_figure, 'img');

        fName = [get(handles.save_fName_edit, 'String') '.mat'];

        save(fName, 'po', 'oo', 'ooVal', 'scale', 'pol', 'img');

        msgbox([fName ' saved successfully!'], 'Save Fit', 'help');

    end

%% ---------------------------------------------------------------------------
    function save_excel_pushbutton_Callback(hObject,evendata) %#ok<INUSD>

        po = getAppD('po', []);
        pol = getAppD('pol', []);

        if isempty(po) || isempty(pol)
            errordlg('Mark curve point first!', 'Excel export error');
        else
            fName = [get(handles.save_fName_edit, 'String') '.xls'];
            headStr = {'Polynomial fit:' '';
                poly2str(pol, 'x') '';
                'Polynomial coefficients:' '';
                '' '';
                '' '';
                'X' 'Y'};
            try
                xlswrite(fName, headStr, 'Sheet1', 'A1');
                xlswrite(fName, reshape(pol', 1, []) , 'Sheet1', 'A4');
                xlswrite(fName, po, 'Sheet1', 'A7');

                msgbox([fName ' saved successfully!'], 'Save Fit', 'help');
            catch
                errordlg({'Could not save to excel' 'Is excel file open?'}, 'Excel export error')
            end
        end

    end

%% ---------------------------------------------------------------------------
    function save_ws_pushbutton_Callback(hObject,evendata) %#ok<INUSD>

        names = {'po', 'oo', 'ooVal', 'scale', 'pol', 'img'};
        for i = 1:length(names)
            assignin( 'base', names{i}, getappdata(handles.img2curve_figure, names{i}) );
        end

        msgbox('Results exported to WS successfully!', 'WS export Fit', 'help');

    end

%% ---------------------------------------------------------------------------
	function import_MarkRefX_edit_Callback(hObject,evendata) %#ok<INUSD>

	end

%% ---------------------------------------------------------------------------
	function import_MarkRefX_edit_CreateFcn(hObject,evendata) %#ok<INUSD>

	end

%% ---------------------------------------------------------------------------
	function import_MarkRefY_edit_Callback(hObject,evendata) %#ok<INUSD>

	end

%% ---------------------------------------------------------------------------
	function import_MarkRefY_edit_CreateFcn(hObject,evendata) %#ok<INUSD>

	end

%% ---------------------------------------------------------------------------
	function curve_PolOrder_edit_Callback(hObject,evendata) %#ok<INUSD>

	end

%% ---------------------------------------------------------------------------
	function curve_PolOrder_edit_CreateFcn(hObject,evendata) %#ok<INUSD>

	end

%% ---------------------------------------------------------------------------
	function save_fName_edit_Callback(hObject,evendata) %#ok<INUSD>

	end

%% ---------------------------------------------------------------------------
	function save_fName_edit_CreateFcn(hObject,evendata) %#ok<INUSD>

    end

%% Service Functions

    function scale = getAxisScale(ax, scale)

        directionStr = {'x \rightarrow'   'y \uparrow'};

        axes(handles.disp_ax);
        t = title(handles.disp_ax,{ ['Get ' directionStr{ax} ' scale'], 'Mark 2 points, from small to large'}, ...
            'FontSize', 16, 'Color', 'r');
        Po = ginput(2);
        if isempty(Po)
            %     load('scale.mat')
            %     scale = scale(ax);
            return
        else
            PoSize = ( Po(2,ax)-Po(1,ax) ) / scale;
        end

        Answer = inputdlg('Difference between the 2 points:', 'Get Scale');
        Me = str2num(Answer{1}); %#ok<*ST2NM>
        scale = Me/PoSize;
        delete(t);
    end

    function setOoVal

        ooVal(1) = str2num(get(handles.import_MarkRefX_edit, 'String'));
        ooVal(2) = str2num(get(handles.import_MarkRefY_edit, 'String'));
        setappdata(handles.img2curve_figure, 'ooVal', ooVal);
    end

    function checkSaveScale

        if isappdata(handles.img2curve_figure, 'scale') && isappdata(handles.img2curve_figure, 'oo')
            scale = getappdata(handles.img2curve_figure, 'scale'); %#ok<*NASGU>
            oo = getappdata(handles.img2curve_figure, 'oo');
            ooVal = getappdata(handles.img2curve_figure, 'ooVal');

            save('scale.mat', 'scale', 'oo', 'ooVal');

            plotImBackground;
        end

    end

    function pol = calcPol

        po = getappdata(handles.img2curve_figure, 'po');
        oo = getappdata(handles.img2curve_figure, 'oo');
        ooVal = getappdata(handles.img2curve_figure, 'ooVal');
        scale = getappdata(handles.img2curve_figure, 'scale');

        polOrder = str2num(get(handles.curve_PolOrder_edit, 'String'));
        pol = Fitting(po(:,1), po(:,2), polOrder+1);

        setappdata(handles.img2curve_figure, 'pol', pol);

        set(handles.curve_Pol_text, 'String', poly2str(pol, 'x'));
        showRes;

    end

    function showRes

        po = getappdata(handles.img2curve_figure, 'po');
        oo = getappdata(handles.img2curve_figure, 'oo');
        ooVal = getappdata(handles.img2curve_figure, 'ooVal');
        scale = getappdata(handles.img2curve_figure, 'scale');
        pol = getappdata(handles.img2curve_figure, 'pol');

        if isappdata(handles.img2curve_figure, 'lin')
            try %#ok<TRYNC>
                lin = getappdata(handles.img2curve_figure, 'lin');
                delete(lin);
            end
        end

        x = po(:,1);
        y = po(:,2);
        yP = polyval( pol, po(:,1) );

        axes(handles.disp_ax);

        hold on
        % imagesc(scale(1), scale(2), img);
        lin = plot(x,y, 'o', x,yP,'--', 'linewidth', 2);
        hold off
        tit = title( poly2str(pol, 'x'), 'FontSize', 12, 'Color', 'k' );

        setappdata(handles.img2curve_figure, 'lin', lin);

    end

    function plotImBackground

        po = getappdata(handles.img2curve_figure, 'po');
        oo = getappdata(handles.img2curve_figure, 'oo');
        ooVal = getappdata(handles.img2curve_figure, 'ooVal');
        scale = getappdata(handles.img2curve_figure, 'scale');
        pol = getappdata(handles.img2curve_figure, 'pol');
        img = getappdata(handles.img2curve_figure, 'img');

        if isempty(ooVal)
            ooVal = [0 0];
        end

        oo = oo - ooVal./scale;

        minX = -1 * scale(1) * oo(1);
        maxX = scale(1) * (length(img(1,:,1)) - oo(1));

        minY = -1 * scale(2) * oo(2);
        maxY = scale(2) * (length(img(:,1,1)) - oo(2));

        imagesc([minX maxX], [minY maxY], flipdim(img,1), 'Parent', handles.disp_ax);
        set(handles.disp_ax,'ydir','normal');

    end

    function pushbutton_activate(indicat)

        pb = findobj(handles.img2curve_figure, 'Style', 'pushbutton');

        if nargin<1   % activate all
            set(pb, 'Enable', 'on');
        else
            pbb = [];
            for i = 1:length(pb)
                if strfind(get(pb(i), 'Tag'), indicat)
                    pbb(end+1) = pb(i); %#ok<AGROW>
                end
            end
            set(pbb, 'Enable', 'on');
        end

    end

    function pushbutton_deactivate

        pb = findobj(handles.img2curve_figure, 'Style', 'pushbutton');

        set(pb, 'Enable', 'off');

    end

    function val = getAppD(name, def)

        if isappdata(handles.img2curve_figure, name)
            val = getappdata(handles.img2curve_figure, name);
        else
            val =def;
        end

    end

end

function Pol = Fitting(x, y, m)

if nargin<3
    m = 2;
end

x = MakeColomVec(x);
y = MakeColomVec(y);

for q=1:m
    A=[];
    for p=1:q
        A=[A x.^(p-1)]; %#ok<AGROW>
    end
    d(1:q,q) =A\y; %#ok<AGROW>
    % plot(x, y, x,A*d(1:q,q),'--')
    % hold on
end

Pol = d(1:q,q);
Pol = Pol(end:-1:1);

    function v = MakeColomVec(v)

        v = reshape(v, length(v), 1);
    end
end

author: Sable

repo: mcbench-benchmarks


text:
function pass = test_pchip(pref)

% Test a scalar function:
x = 0:10;
y = sin(x);
f = chebfun.pchip(x, y);
tol = 10*eps;
pass(1) = norm(feval(f, x) - y) < tol;
pass(2) = numel(f.funs) == 10;
pass(3) = length(f) == 40;

% Test an array-valued function:
x = (0:10)';
y = [sin(x), cos(x)];
f = chebfun.pchip(x, y);
tol = 10*eps;
pass(4) = norm(feval(f, x) - y) < tol;
pass(5) = numel(f.funs) == 10;
pass(6) = length(f) == 40;

% Test a different domain:
dom = [.01, 10.01];
x = (0:10)';
y = [sin(x), cos(x)];
f = chebfun.pchip(x, y, dom);
tol = 10*eps;
pass(7) = all(f.domain == [dom(1), 1:10, dom(2)]);
pass(8) = norm(feval(f, x(2:end)) - y(2:end,:)) < tol;
pass(9) = numel(f.funs) == 11;
pass(10) = length(f) == 44;

% Test the example from the m-file:
x = -3:3;
y = [-1 -1 -1 0 1 1 1];
f = chebfun.pchip(x, y);
tol = 10*eps;
pass(11) = norm(feval(f, x) - y) < tol;
pass(12) = numel(f.funs) == 6;

end

author: chebfun

repo: chebfun

save_path: github-repos/MATLAB/chebfun-chebfun


text:
function [ Th ] = Threshold( image, t )
h=size(image,1);
w=size(image,2);
Th=zeros(h,w);
parfor i=1:h
    for j=1:w
        if image(i,j)>t
            Th(i,j)=1;
        else
            Th(i,j)=0;
        end
    end
end
end

author: zhangqianqianQQ

repo: MachineVisionAlgorithm


%   getBox:  Given the current state, returns a number from 1 to 162
%             designating the region of the state space encompassing the current state.
%             Returns a value of -1 if a failure state is encountered.

function box = getBox3(theta,thetaDot,x,xDot)
theta = rad2deg(theta);
thetaDot = rad2deg(thetaDot);
if (x < -1 || x > 1  || theta < -8 || theta > 8)
    box = -1;
else

if (theta<-1&&theta>=-8)
	thetaBucket = 1;
elseif (theta<0&&theta>=-1)
	thetaBucket = 2;
elseif (theta<1&&theta>=0)	% zero included
	thetaBucket = 3;
elseif (theta<8&&theta>=1)
	thetaBucket = 4;
end

if (x<-0.8&&x>=-1)
	xBucket = 1;
elseif (x<=0.8&&x>=-0.8)
	xBucket = 2;
elseif (x<=1&&x>0.8)
	xBucket = 3;
end

if (xDot<-0.5)
	xDotBucket = 1;
elseif (xDot>=-0.5&&xDot<=0.5)
	xDotBucket = 2;
else
	xDotBucket = 3;
end

if (thetaDot<-50)
	thetaDotBucket = 1;
elseif (thetaDot>=-50&&thetaDot<=50)
	thetaDotBucket = 2;
else
	thetaDotBucket = 3;
end

box = sub2ind([4,3,3,3],thetaBucket, thetaDotBucket,xBucket,xDotBucket);
end
return;


text:
classdef(Abstract) AbstractIntegrationTerminationCause < matlab.mixin.SetGet & matlab.mixin.Heterogeneous
    %AbstractIntegrationTerminationCause Summary of this class goes here
    %   Detailed explanation goes here

    properties

    end

    methods
        tf = shouldRestartIntegration(obj)

        newStateLogEntry = getRestartInitialState(obj, stateLogEntry)
    end

    methods(Static, Access=protected)
        function defaultObject = getDefaultScalarElement()
            defaultObject = SoITransitionDownIntTermCause();
        end
    end
end

author: Arrowstar

repo: ksptot

save_path: github-repos/MATLAB/Arrowstar-ksptot


text:
% Propagation of uncertainty as in Atlas
% Propagating uncertainty from node to node using relative motion.
%
%
function propagation_of_uncertainty_example
close all

poses = [0 0 0;
         10 0 0;
         10 10 pi/2;
         10 20 pi/2;
         0 20 pi;
         0 10 -pi/2;
         -10 -10 0];
n = 7
% figure, hold on
% plot(poses(:, 1), poses(:, 2), 'o')

S00 = diag([0.0 0.0 0.0])
sigmas = {}
sigmas{1} = S00
S = S00
for i=1:n-1
    S = propagate_uncertainty(poses(i, :), poses(i+1, :) , S)
    sigmas{i+1} = S;
end

% axis equal, hold on
% for i=1:n
%     plot_error_ellipse(poses(i,1:2), sigmas{i}(1:2, 1:2), 0.99)
% end

M = 500;
samples = zeros(M, 3);
sets = {};
sets{1} = samples;
for i=1:n-1
    samples = propagate_uncertainty_sampled(poses(i, :), poses(i+1, :) , samples);
    sets{i+1} = samples;
end

% plot everything
figure, hold on
plot(poses(:, 1), poses(:, 2), 'o')
axis equal, hold on
for i=1:n
   plot_error_ellipse(poses(i,1:2), sigmas{i}(1:2, 1:2), 0.99)
   plot(sets{i}(:,1), sets{i}(:,2), '.')
end

function Sres=propagate_uncertainty(posea, poseb, Si)

Sij = diag([0.1, 0.01, 0.001]);
Ti = T(posea(1), posea(2), posea(3));
Tj = T(poseb(1), poseb(2), poseb(3));
Tij = inv(Ti)*Tj;
xij = Tij(1,4);
yij = Tij(2,4);
thij = atan2(Tij(2,1), Tij(1,1));

[J1, J2] = Jacobians(posea(1), posea(2), posea(3), xij, yij, thij);

Sres = J1*Si*J1' + J2*Sij*J2';

function Set=propagate_uncertainty_sampled(ti, tj, Set)

Sjk = diag([0.1, 0.01, 0.001]);
Sjk = sqrt(Sjk);
Ti = T(ti(1), ti(2), ti(3));
Tj = T(tj(1), tj(2), tj(3));
Tij = inv(Ti)*Tj;
tij = t2v(Tij);
xij = tij(1);
yij = tij(2);
thij = tij(3);

for i=1:length(Set)
    xis = Set(i,1);
    yis = Set(i,2);
    this = Set(i, 3);
    Tis = T(xis, yis, this);
    % sample-based error propagation, add noise
    tx = xij + normrnd(0, Sjk(1,1));
    ty = yij + normrnd(0, Sjk(2,2));
    th = thij + normrnd(0, Sjk(3,3));
    Tijs = T(tx, ty, th);
    Tjs = Tis*Tijs;
    Set(i, :)= t2v(Tjs);
end

function A = T(x, y, th)

cth = cos(th);
sth = sin(th);

A=[cth  -sth  0   x;
   sth  cth   0   y;
    0   0     1   0;
    0   0     0   1];

function t = t2v(T)
x = T(1,4);
y = T(2,4);
th = atan2(T(2,1), T(1,1));
t = [x, y, th];

function [J1, J2] = Jacobians(xi, yi, thi, xij, yij, thij)

ci = cos(thi);
si = sin(thi);

J1=[1  0 -si*xij-ci*yij;
   0  1  ci*xij-si*yij;
   0   0   1];

J2=[ci -si   0;
   si  ci   0;
   0    0   1];

function plot_error_ellipse(mu, Sigma, p)
s = -2 * log(1 - p);
[V, D] = eig(Sigma * s);
t = linspace(0, 2 * pi, 50);
a = (V * sqrt(D)) * [cos(t(:))'; sin(t(:))'];
plot(a(1, :) + mu(1), a(2, :) + mu(2));

repo: ARTE


text:
function xi = cvx_cleanup_structure( xi )

%CVX_CLEANUP_STRUCTURE   Canonicalizes a matrix structure basis.
%    CVX_CLEANUP_STRUCTURE(X), where X is an m x n matrix, converts X
%    to row reduced echelon form with the zero rows removed. The matrices
%    X are assumed to come from CVX's matrix structure facility, and as
%    such certain assumptions are made about both X and its RREF: that it
%    is sparse and its nonzero elements are ratios of small integers. As a
%    result this routine should not be used for general RREF computations
%    without modifying it to remove some of the cleanup code.

% Reduce using an LU factorization
[LL,xi,PP] = lu(xi); %#ok
[m,n] = size(xi);

% Remove the entries that are close to zero
tol = 16 * eps;
xi  = xi .* ( abs(xi) > tol * norm(xi,'inf') );

% Find the locations of the leading element in each row. To do this we first
% find the first element in each row. Transposing xi insures that the
% indices are sorted properly to accomplish this.
[jj,ii] = find(xi');
if isempty(jj),
    xi = sparse(0,n);
    return
end
dd = [true;diff(ii)~=0];
ii = ii(dd);
jj = jj(dd);

% Sort the rows so that the leftmost nonzero is first (the LU factorization
% does this already much of the time, but in rank-degenerate cases further
% sorting is needed.) From that select a unique set of columns to use.
[jj,jndx] = sort(jj);
dd = [true;diff(jj)~=0];
ii = ii(jndx(dd));
jj = jj(dd);

% Divide through the rows ii by this full-rank triangle. The other rows
% must already be zero or be dependent upon these rows, so we remove them.
% Q is the result except that its columns are scrambled.
rr = length(ii);
j2 = (1:n)'; j2(jj) = [];
Q  = xi(ii,jj)\xi(ii,j2);

% Reduce roundoff error by converting the values to ratios of integers.
[i3,j3,vv] = find(Q);
[vn,vd] = rat(vv,tol);
xi = sparse([(1:rr)';i3],[jj;j2(j3)],[ones(rr,1);vn./vd],rr,n);

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

repo: REST


%==============================================================================
% This code is part of the Matlab-based toolbox
% FAIR - Flexible Algorithms for Image Registration.
% For details see
% - https://github.com/C4IR and
% - http://www.siam.org/books/fa06/
% ##2
%==============================================================================
%
% Tutorial for FAIR: regulariztion
%
% given a certain force field f on a 2D domain, this tutorial
% computes and visualizes the curvature regularized displacement u, such that
%
%  B'*B u = f
%
% where B is the discrete curvature operator
% see also getCurvatureMatrix
%==============================================================================

clear, close all, help(mfilename);

                    % 2D example
omega  = [0,1,0,1]; % physical domain
m      = [16,12];   % number of discretization points
n      = prod(m);   % number of cells, and staggered dimensions
hd     = prod((omega(2:2:end)-omega(1:2:end))./m);
alpha  = 10;

% generate cell centered and nodal grids
% note: computation is cell centered, others used for visualization
xc = getCellCenteredGrid(omega,m);
xn = getNodalGrid(omega,m);

% build elasticity operator on a staggered grid and visualize
B = getCurvatureMatrix(omega,m);
FAIRfigure(1); clf; subplot(2,2,1); spy(B); title('B curcature on cell-centered grid')

% create a force field o cell centered grid and visualize
fc  = [0*xc(1:n);-10*sin(pi*xc(1:n))];

% compute the displacement from (B'*B)*us = fs
% note B is derivative operator and has null space (constants)

% remove constant parts: E'*(f-E*w) = 0
E = [ones(n,1)*[1,0];ones(n,1)*[0,1]];
fc = fc - E*((E'*E)\(E'*fc));

% solve for displacements
warning off % matrix is singular, but MATLAB doesn't care
uc = hd*alpha*(B'*B)\fc;
warning on

% visualization on cell centered and nodal grids
un = grid2grid(uc,m,'centered','nodal');

% shortcuts for visualization,
vecNorm   = @(v) sqrt(sum(reshape(v,[],2).^2,2));
viewField = @(v) viewImage2Dsc(vecNorm(v),omega,m);
PG        = @(x) plotGrid(x,omega,m);

% visualize force field
FAIRfigure(1); subplot(2,2,3);
viewField(fc); hold on; colormap(gray); colorbar; PG(xn);
qh = quiver(xc(1:n),xc(n+1:end),fc(1:n),fc(n+1:end),1);
set(qh,'color','r','linewidth',1)
title('a force field')

% visualize displacement field
FAIRfigure(1); subplot(2,2,4);
viewField(uc); hold on; colormap(gray); colorbar; PG(xn);
qh = quiver(xc(1:n),xc(n+1:end),uc(1:n),uc(n+1:end),1);
set(qh,'color','g','linewidth',1)
title('the displacement field')

% visuaize the displaced grid
FAIRfigure(1); subplot(2,2,2);
xh = PG(xn); hold on; axis equal image; yh = PG(xn+un);
set(xh,'linewidth',2); set(xh,'linewidth',2,'color','g');
title('initial and displaced grids')

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


text:
function [T,degree,t,s] = tform2myT(tform)
T = tform.T';
T = inv(T);

t = T([1,2],3);

% assume s1 and s2 appear on the first and the second row of T
s = sqrt(T(1:2,1).^2 + T(1:2,2).^2);
R = T(1:2,1:2) ./ repmat(s,[1,2]);
theta = atan2(R(2,1), R(1,1));
degree = theta * 180/pi;

author: zhouyuanzxcv

repo: Hyperspectral


text:
function [MU,k]=vonMisesStat(T)

% function [MU,k]=vonMisesStat(T)
% ------------------------------------------------------------------------
%
%
% ------------------------------------------------------------------------

%%
MU=angle(mean(exp(1i.*T)));
Rsq=mean(cos(T)).^2+mean(sin(T)).^2;
R=sqrt(Rsq);

p=2;
ki=R.*(p-Rsq)./(1-Rsq);
diffTol=ki/1000;
qIter=1;
while 1
    kn=ki;
    A=besseli(p/2,ki)./ besseli((p/2)-1,ki);
    ki=ki-((A-R)./(1-A^2-((p-1)/ki)*A));
    if abs(kn-ki)<=diffTol
        break
    end
    qIter=qIter+1;
end
k=ki;

%%
% _*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:
classdef unimodalComponent
% This class is obsolet since MTEX 5.9. Use the class @SO3FunRBF instead.
% Anyway the class is preserved, so that saved @ODFs can be loaded.

methods (Static = true, Hidden=true)
  function odf = loadobj(s)
    psi = s.psi;
    if isempty(psi)
      psi = SO3DeLaValleePoussinKernel;
      warning(['MTEX is not able to load the SO3Kernel of the unimodalODF, ' ...
        'because its class has been deleted. A standard kernel is used instead.'])
    end
    odf = SO3FunRBF(s.center,psi,s.weights);
  end
end

end

author: mtex-toolbox

repo: mtex

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


text:
%IMATCH Template matching
%
% XM = IMATCH(IM1, IM2, U, V, H, S) is the position of the matching subimage of
% IM1 (template) within the image IM2.  The template in IM1 is centred at (U,V)
% and its half-width is H.
%
% The template is searched for within IM2 inside a rectangular region, centred
% at (U,V) and whose size is a function of S.  If S is a scalar the search
% region is [-S, S, -S, S] relative to (U,V).  More generally S is a 4-vector
% S=[umin, umax, vmin, vmax] relative to (U,V).
%
% The return value is XM=[DU,DV,CC] where (DU,DV) are the u- and v-offsets
% relative to (U,V) and CC is the similarity score for the best match in the
% search region.
%
% [XM,SCORE] = IMATCH(IM1, IM2, U, V, H, S) as above but also returns a matrix
% of matching score values for each template position tested. The rows
% correspond to horizontal positions of the template, and columns the vertical
% position.  The centre element corresponds to (U,V).
%
% Example::
%  Consider a sequence of images im(:,:,N) and we find corner points in the
%  k'th image
%        corners = icorner(im(:,:,k), 'nfeat', 20);
%  Now, for each corner we look for the 11x11 patch of surrounding pixels
%  in the next image, by searching within a 21x21 region
%        for corner=corners
%           xm = imatch(im(:,:,k), im(:,:,k+1), 5, 10);
%           if xm(3) > 0.8
%               fprintf('feature (%f,%f) moved by (%f,%f) pixels)\n', ...
%                   corner.u, corner.v, xm(1), xm(2) );
%           end
%        end
%
% Notes::
% - Useful for tracking a template in an image sequence where IM1 and IM2
%   are consecutive images in a template and (U,V) is the coordinate of
%   a corner point in IM1.
% - Is a MEX file.
% - IM1 and IM2 must be the same size.
% - ZNCC (zero-mean normalized cross correlation) matching is used as the
%   similarity measure. A perfect match score is 1.0 but anything above 0.8
%   is typically considered to be a good match.
%
% See also ISIMILARITY.

% Copyright (C) 1993-2011, by Peter I. Corke
%
% This file is part of The Machine Vision Toolbox for Matlab (MVTB).
%
% MVTB is free software: you can redistribute it and/or modify
% it under the terms of the GNU Lesser General Public License as published by
% the Free Software Foundation, either version 3 of the License, or
% (at your option) any later version.
%
% MVTB is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
% GNU Lesser General Public License for more details.
%
% You should have received a copy of the GNU Leser General Public License
% along with MVTB.  If not, see <http://www.gnu.org/licenses/>.

if ~exist('imatch', 'file')
    error('you need to build the MEX version of imatch, see vision/mex/README');
end

author: petercorke

repo: machinevision-toolbox-matlab


text:
classdef Preferences < wt.model.BaseModel
    % Base class for managing app preferences
    %
    % Notes:
    %   The preferences may be set to automatically load upon creating an
    %   instance, and save upon destruction of an instance of this class.

    % Copyright 2020-2021 The MathWorks Inc.

    %% Preference Properties
    properties (AbortSet)

        % Starting window position
        Position (1,4) double = [100 100 1000 700]

        % Last folder that was used when opening a file
        LastFolder (1,1) string {mustBeFolder} = pwd

        % List of recent session files
        RecentSessionPaths (:,1) string {mustBeFile}

    end %properties

    %% Saving and Loading Preferences
    methods (Sealed)

        function load(obj,prefGroup)
            % Load stored preferences

            % Get a list of preference properties
            propList = getPreferenceProperties(obj);

            % Get the preference values
            if ispref(prefGroup)
                prefValues = getpref(prefGroup);
            else
                prefValues = struct();
            end

            % Loop on each preference field
            for idx = 1:numel(propList)

                % Get the current preference
                thisProp = propList(idx).Name;
                thisDefault = propList(idx).DefaultValue;
                if isfield(prefValues,thisProp)
                    thisValue = prefValues.(thisProp);
                else
                    thisValue = thisDefault;
                end

                % Get the preference
                try
                    obj.(thisProp) = thisValue;
                catch err
                    warning('wt:Preferences:LoadFailure',...
                        ['Unable to load stored preference for %s:%s. '...
                        'Reverting to default value.\n%s'],...
                        prefGroup,thisProp,err.message);
                    obj.(thisProp) = thisDefault;
                end

            end %for idx = 1:numel(propList)

        end %function

        function save(obj,prefGroup)
            % Save preferences

            % Get a list of preference properties
            propList = getPreferenceProperties(obj);

            % Remove existing preferences first (in case the props changed)
            if ispref(prefGroup)
                rmpref(prefGroup);
            end

            % Loop on each preference field
            if ~isempty(prefGroup)
                for idx = 1:numel(propList)

                    % Get the current preference
                    thisProp = propList(idx).Name;
                    thisValue = obj.(thisProp);

                    % Set the preference
                    setpref(prefGroup,thisProp,thisValue);

                end %for idx = 1:numel(propList)
            end %if ~isempty(prefGroup)

        end %function

    end %methods

    %% Private methods
    methods (Access=private)

        function propList = getPreferenceProperties(obj)
            % Get the info for preference properties

            %mc = metaclass(obj);
            %propList = mc.PropertyList;
            propNames = properties(obj);
            for idx = numel(propNames):-1:1
                propList(idx) = findprop(obj,propNames{idx});
            end
            propsToKeep = ~[propList.Constant] & ...
                ~[propList.Transient] & ...
                ~[propList.Hidden] & ...
                [propList.HasDefault] & ...
                strcmp({propList.GetAccess},'public') & ...
                strcmp({propList.SetAccess},'public');
            propList(~propsToKeep) = [];

        end %function

    end %methods

end % classdef


text:
% This is the code of Hybrid_MSD algorithm:
%
% The code of Hybrid_MSD is provided by the authors of Hybrid_MSD.
% The interface is created by the authors of VIFB.

function img = run_Hybrid_MSD(imgVI, imgIR, visualization)

    %    The following is an implementation of an infrared and visible image
    %    fusion algorithm proposed in the paper:
    %
    %      "Perceptual fusion of infrared and visible images through a hybrid
    %      multi-scale decomposition with Gaussian and bilateral filters"
    %      Information Fusion, 2016.
    %
    %    This code is for testing purpose only.
    %    Some of the test images were obtained at
    %      http://www.imagefusion.org
    %      http://www.ece.lehigh.edu/SPCRL/IF/image_fusion.htm
    %
    %    Zhiqiang Zhou, Beijing Institute of Technology
    %    May, 2015

    path_Vis = imgVI.img;
    path_IR  = imgIR.img;

    % IR image
    img1 = imread(path_IR);

    % VI image
    img2 = imread(path_Vis);

    img1 = double(img1);
    img2 = double(img2);

    if visualization == 1
        paraShow1.fig = 'Visible image';
        paraShow2.fig = 'Infrared image';
        ShowImageGrad(img2, paraShow2);
        ShowImageGrad(img1, paraShow1);
    end

    tic;
    if size(img2, 3) == 1
        fuseimage = Hybrid_MSD(img2, img1);
    elseif size(img1,3) == 1
        fuseimage = zeros(size(img2));
        for i=1:3
            fuseimage(:,:,i) = Hybrid_MSD(img2(:,:,i),img1);
        end
    else
        fuseimage = zeros(size(img2));
        for i=1:3
           fuseimage(:,:,i) = Hybrid_MSD(img2(:,:,i),img1(:,:,i));
        end
    end
    toc;
    if visualization == 1
        figure,imshow(fuseimage, []);
    end
    img = uint8(fuseimage);

end

function res = Hybrid_MSD(img1, img2)

    nLevel = 4;
    lambda = 30;
    % lambda = 3000;

    %% ---------- Hybrid Multi-scale Decomposition --------------
    sigma = 2.0;
    sigma_r = 0.05;
    k = 2;

    M1 = cell(1, nLevel+1);
    M1L = cell(1, nLevel+1);
    M1{1} = img1;
    M1L{1} = img1;
    M1D = cell(1, nLevel+1);
    M1E = cell(1, nLevel+1);
    sigma0 = sigma;
    for j = 2:nLevel+1,
        w = floor(3*sigma0);
        h = fspecial('gaussian', [2*w+1, 2*w+1], sigma0);
        M1{j} = imfilter(M1{j-1}, h, 'symmetric');
        %M1L{j} = 255*bfilter2(M1L{j-1}/255,w,[sigma0, sigma_r/(k^(j-2))]);
        M1L{j} = 255*fast_bfilter2(M1L{j-1}/255,[sigma0, sigma_r/(k^(j-2))]);

        M1D{j} = M1{j-1} - M1L{j};
        M1E{j} = M1L{j} - M1{j};

        sigma0 = k*sigma0;
    end

    M2 = cell(1, nLevel+1);
    M2L = cell(1, nLevel+1);
    M2{1} = img2;
    M2L{1} = img2;
    M2D = cell(1, nLevel+1);
    M2E = cell(1, nLevel+1);
    sigma0 = sigma;
    for j = 2:nLevel+1,
        w = floor(3*sigma0);
        h = fspecial('gaussian', [2*w+1, 2*w+1], sigma0);
        M2{j} = imfilter(M2{j-1}, h, 'symmetric');
        %M2L{j} = 255*bfilter2(M2L{j-1}/255,w,[sigma0, sigma_r/(k^(j-2))]);
        M2L{j} = 255*fast_bfilter2(M2L{j-1}/255,[sigma0, sigma_r/(k^(j-2))]);

        M2D{j} = M2{j-1} - M2L{j};
        M2E{j} = M2L{j} - M2{j};

        sigma0 = k*sigma0;
    end

    %% ---------- Multi-scale Combination --------------
    for j = nLevel+1:-1:3
    b2 = abs(M2E{j});
    b1 = abs(M1E{j});
    R_j = max(b2-b1, 0);
    Emax = max(R_j(:));
    P_j = R_j/Emax;

    C_j = atan(lambda*P_j)/atan(lambda);

    % Base level combination
    sigma0 = 2*sigma0;
    if j == nLevel+1
        w = floor(3*sigma0);
        h = fspecial('gaussian', [2*w+1, 2*w+1], sigma0);
        lambda_Base = lambda;
        %lambda_Base = 30;
        C_N = atan(lambda_Base*P_j)/atan(lambda_Base);
        C_N = imfilter(C_N, h, 'symmetric');
        MF = C_N.*M2{j} + (1-C_N).*M1{j};
    end

    % Large-scale combination
    sigma0 = 1.0;
    w = floor(3*sigma0);
    h = fspecial('gaussian', [2*w+1, 2*w+1], sigma0);
    C_j = imfilter(C_j, h, 'symmetric');

    D_F = C_j.*M2E{j}+ (1-C_j).*M1E{j};
    MF = MF + D_F;
    D_F = C_j.*M2D{j}+ (1-C_j).*M1D{j};
    MF = MF + D_F;
    end

    % Small-scale combination
    sigma0 = 0.2;
    w = floor(3*sigma0);
    h = fspecial('gaussian', [2*w+1, 2*w+1], sigma0);
    C_0 = double(abs(M1E{2}) < abs(M2E{2}));
    C_0 = imfilter(C_0, h, 'symmetric');
    D_F = C_0.*M2E{2} + (1-C_0).*M1E{2};
    MF = MF + D_F;
    C_0 = abs(M1D{2}) < abs(M2D{2});
    C_0 = imfilter(C_0, h, 'symmetric');
    D_F = C_0.*M2D{2} + (1-C_0).*M1D{2};
    MF = MF + D_F;

    %% ---------- Fusion Result --------------
    % FI = ImRegular(MF);   % The intensities are regulated into [0, 255]
    FI = max(min(MF,255), 0);
    res = FI;
end

author: xingchenzhang

repo: VIFB

save_path: github-repos/MATLAB/xingchenzhang-VIFB


text:
% SP_TO_VTK: Export multipatch results to VTK format for plotting.
%
%  sp_to_vtk (u, space, geometry, npts, filename, fieldnames, [option], [lambda_lame, mu_lame])
%
% INPUT:
%
%     u:           vector of dof weights
%     space:       object representing the space of discrete functions (see sp_multipatch)
%     geometry:    geometry structure (see geo_load)
%     npts:        number of points along each parametric direction where to evaluate
%     filename:    name of the output file.
%     fieldnames:  how to name the saved variables in the vtk file
%     options:     cell array with the fields to plot
%                   accepted options are 'value' (default), 'gradient',
%                   and for vectors also 'curl', 'divergence', 'stress'
%     lambda_lame: function handle to the first Lame coefficient (only needed to compute 'stress')
%     mu_lame:     function handle for the second Lame coefficient (only needed to compute 'stress')
%
% OUTPUT:
%
%    none
%
% Copyright (C) 2010 Carlo de Falco, Rafael Vazquez
% Copyright (C) 2011, 2012, 2015 Rafael Vazquez
%
%    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/>.

function sp_to_vtk (u, space, geometry, npts, filename, fieldname, varargin)

  str1 = cat (2,'<?xml version="1.0"?> \n', ...
'<VTKFile type="Collection" version="0.1"> \n', ...
'<Collection> \n');

  str2 = cat (2, '<DataSet part="%d" file="%s.vts"/> \n');

  str3 = cat (2, ...
'</Collection>\n', ...
'</VTKFile> \n');

  if (length (filename) < 4 || ~strcmp (filename(end-3:end), '.pvd'))
    pvd_filename = cat (2, filename, '.pvd');
  else
    pvd_filename = filename;
    filename = filename (1:end-4);
  end

  fid = fopen (pvd_filename, 'w');
  if (fid < 0)
    error ('mp_sp_to_vtk: could not open file %s', pvd_filename);
  end

  fprintf (fid, str1);
  ind = union (find (filename == '/', 1, 'last'), find (filename == '\', 1, 'last')) + 1;
  if (isempty (ind)); ind = 1; end
  for iptc = 1:space.npatch
    filename_patch_without_path = cat (2, filename(ind:end), '_', num2str (iptc));
    filename_patch = cat (2, filename, '_', num2str (iptc));
    fprintf (fid, str2, iptc, filename_patch_without_path);
    if (isempty (space.dofs_ornt))
      sp_to_vtk (u(space.gnum{iptc}), space.sp_patch{iptc}, geometry(iptc), npts, ...
                           filename_patch, fieldname, varargin{:})
    else
      sp_to_vtk (u(space.gnum{iptc}) .* space.dofs_ornt{iptc}', space.sp_patch{iptc}, geometry(iptc), npts, ...
                           filename_patch, fieldname, varargin{:})
    end
  end
  fprintf (fid, str3);

  fclose (fid);

end

author: rafavzqz

repo: geopdes

save_path: github-repos/MATLAB/rafavzqz-geopdes


text:
function [output, layerOut] = FeatureTree2(visible, para, layer)
% remove target and cost function nodes if any
for i=1:length(layer)
    if strcmpi(layer{i}.name, 'target') || strcmpi(layer{i}.name, 'cross_entropy') || strcmpi(layer{i}.name, 'mse') || strcmpi(layer{i}.name, 'logistic')
        layer{i}.name = 'ignore';
    end
    if i > max(para.out_layer_idx)  % we don't compute layers that is after the last output layer
        layer{i}.name = 'ignore';
    end
end
para.NET.sentenceMinibatch = 1;
[minibatch] = MinibatchPackaging_tree4(visible, para);

output = {};
for utt_i = 1:minibatch.nBatch
    PrintProgress(utt_i, minibatch.nBatch, 100);
    batch_data = GetMinibatch2(minibatch, para, utt_i);

    % Use mode=3 to generate network output only
    if nargout>1
        [~, layerOut{utt_i}, output{utt_i}] = DNN_Cost10(layer, batch_data, para,3);
    else
        [~, ~, output{utt_i}] = DNN_Cost10(layer, batch_data, para,3);
    end
end
end

author: singaxiong

repo: SignalGraph


text:
function exact = p03_exact ( )

%*****************************************************************************80
%
%% P03_EXACT returns the exact integral for problem 3.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    17 January 2009
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Output, real EXACT, the value of the integral.
%
  exact = ( 16.0 / 3.0 ) * ( 2.0 - sqrt ( 2.0 ) );

  return
end

author: johannesgerer

repo: jburkardt-m


text:
function angles2d
%ANGLES2D Description of functions for manipulating angles.
%
%   Angles are normalized in an interval of width 2*PI. Most geom2d
%   functions return results in the [0 2*pi] interval, but it can be
%   convenient to consider the [-pi pi] interval as well. See the
%   normalizeAngle function to switch between conventions.
%
%   Angles are usually oriented. The default orientation is the CCW
%   (Counter-Clockwise) orientation.
%
%   See also
%   normalizeAngle, angleDiff, angleAbsDiff, angleSort
%   angle2Points, angle3Points, vectorAngle, lineAngle, edgeAngle
%

% ------
% Author: David Legland
% E-mail: david.legland@inrae.fr
% Created: 2010-03-31, using Matlab 7.4.0.287 (R2007a)
% Copyright 2010-2022 INRA - Cepia Software Platform

help('angles2d');

author: mattools

repo: matGeom

save_path: github-repos/MATLAB/mattools-matGeom


author: beckel

repo: nilm-eval

save_path: github-repos/MATLAB/beckel-nilm-eval


text:
% STD_LM_SETEEGFIELDS - set limo fields in eeg sets
%
% Usage:
%   >>   EEG = std_lm_seteegfields(STUDY,EEG,index,'datatype','Channels','format', 'cell')
%
% Inputs:
%      STUDY    - studyset structure containing some or all files in EEG
%      EEG      - EEG dataset
%      index    - index of dataset in STUDY.datasetinfo
%
% Optional inputs:
%       format    - 'matrix' (default) or 'cell'
%       datatype  - 'channels' or 'components'
%
% Outputs:
%      EEG    - updated EEG structure
%
% See also:
%
% Authors: Cyril Pernet, The university of Edinburgh
%          Arnaud Delorme, SCCN
%          Ramon Martinez-Cancino, SCCN

% Copyright (C) 2015  Ramon Martinez-Cancino,INC, SCCN
%
% This file is part of EEGLAB, see http://www.eeglab.org
% for the documentation and details.
%
% 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 COPYRIGHT HOLDERS AND CONTRIBUTORS "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 COPYRIGHT HOLDER OR CONTRIBUTORS 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 EEG = std_lm_seteegfields(STUDY,EEG,index,varargin)

%% Settings defaults
%  -----------------
try
    options = varargin;
    if ~isempty( varargin )
        for i = 1:2:numel(options)
            opt.(options{i}) = options{i+1};
        end
    else opt= []; end
catch
    error('limo_create_single_trials() error: calling convention {''key'', value, ... } error');
end

try, opt.format;           catch, opt.format        = 'cell';      end
try, opt.datatype;         catch, opt.datatype      = 'channels';  end

% Getting measureflags
% --------------------
measures   = {'erp','spec','ersp','timef','itc'};

% Getting prefix for channels/components
if strncmpi(opt.datatype,'chan',4)
    prefix = 'dat';
elseif strncmpi(opt.datatype,'comp',4)
    prefix = 'ica';
end

% Getting values
for i = 1:length(measures)
    opt.(measures{i}) = STUDY.etc.measureflags.([prefix measures{i}]);
end
% ---

fields = fieldnames(opt);
c = 1;
for i = 1: length(structfun(@numel,opt))
    if ~any ([isempty(eval(['opt.' fields{i}])),strcmp(fields(i),'format'), strcmp(fields(i),'datatype')])
    in_options{c}   = fields{i};
    in_options{c+1} = eval(['opt.' fields{i}]);
    c = c+2;
    end
end

%% Loading datase
%  --------------
path_tmp = rel2fullpath(STUDY.filepath,STUDY.datasetinfo(index).filepath);
name = fullfile(path_tmp, STUDY.datasetinfo(index).subject);
name2 = '';
if ~isempty(STUDY.datasetinfo(index).session) && length(STUDY.session) > 1
    name2 = name;
    name  = fullfile(path_tmp, [ STUDY.datasetinfo(index).subject sprintf('_ses-%2.2d',  STUDY.datasetinfo(index).session) ]);
end

%% Channels: update EEG.set file
%  -----------------------------
if strcmpi(opt.datatype,'channels')
    prefix = 'dat';
else
    prefix = 'ica';
end

% DATERP
if strcmp(opt.erp,'on')
    ext = [ prefix 'erp' ];
    EEG.etc.datafiles.(ext) = getfilename(name, name2, [ '.' ext ]);
    data = load('-mat',EEG.etc.datafiles.(ext));
    EEG.etc.timeerp = data.times;
end

% DATSPEC
if strcmp(opt.spec,'on')
    ext = [ prefix 'spec' ];
    EEG.etc.datafiles.(ext) = getfilename(name, name2, [ '.' ext ]);
    data = load('-mat',EEG.etc.datafiles.(ext));
    EEG.etc.freqspec = data.freqs;
end

% DAT TIMEF
if strcmp(opt.timef,'on')
    ext = [ prefix 'timef' ];
    EEG.etc.datafiles.(ext) = getfilename(name, name2, [ '.' ext ]);
    EEG.etc.datafiles.datersp  = EEG.etc.datafiles.(ext);
    data = load('-mat',EEG.etc.datafiles.(ext),'times','freqs');
    EEG.etc.timeersp = data.times;
    EEG.etc.freqersp = data.freqs;
end

% -------------------------------------------------------------------------
% -------------------------------------------------------------------------
function file_fullpath = rel2fullpath(studypath,filepath)
% Return full path if 'filepath' is a relative path. The output format will
% fit the one of 'filepath'. That means that if 'filepath' is a cell array,
% then the output will a cell array too, and the same if is a string.

nit = 1; if iscell(filepath), nit = length(filepath);end

for i = 1:nit
    if iscell(filepath),pathtmp = filepath{i}; else pathtmp = filepath; end
    if strfind(pathtmp(end),filesep), pathtmp = pathtmp(1:end-1); end % Getting rid of filesep at the end
    if strfind(pathtmp(1:2),['.' filesep])
        if iscell(filepath),
            file_fullpath{i} = fullfile(studypath,pathtmp(3:end));
        else
            file_fullpath = fullfile(studypath,pathtmp(3:end));
        end
    else
        if iscell(filepath),
            file_fullpath{i} = pathtmp;
        else
            file_fullpath = pathtmp;
        end
    end
end

function name = getfilename(name, name2, ext)

if ~exist([name ext],'file')
    if ~exist([name2 ext],'file')
        tmp = dir([name2 '*' ext ]);
        name = fullfile(tmp(1).folder,tmp(1).name);
        warning('couldn''t find a direct match between .set and .daterp, loading %s',name)
    else
        name = [name2 ext];
    end
else
    name = [name ext];
end

author: sccn

repo: eeglab

save_path: github-repos/MATLAB/sccn-eeglab


text:
function [ output_args ] = performEventFunctionAtDoubleHit( t)
%PERFORMEVENTFUNCTIONATDOUBLEHIT
% this function is used to detect the double touch in the z direction of the
% robot end effector.
global binaryValArray
global averageVal
threshhold=4;
[ f ] = getEEF_Force( t );
zforceVal=f{3}-3.3989;
zforceVal=-zforceVal;
if(zforceVal>0)
else
    zforceVal=0;
end

a=0.7;
averageVal=(averageVal*a+zforceVal*(1-a));

binaryVal=averageVal>threshhold;

binaryValArray(1)=[];
binaryValArray=[binaryValArray;binaryVal];

n=max(max(size(binaryValArray)));
sum=0;
    for i=1:(n-1)
        if((binaryValArray(i+1)-binaryValArray(i))==1)
            sum=sum+1;
        else
        end
    end
sum
    if(sum==2)
        EventFunctionAtDoubleHit();
        binaryValArray=binaryValArray*0;
    end
end

repo: KST-Kuka-Sunrise-Toolbox


text:
%% fnc_getSobolSetMatlab: give a set of sobol quasi-random
%
% Usage:
%   X = fnc_getSobolSetMatlab(dim, N)
%
% Inputs:
%    dim                number of variables, the MAX number of variables is 40
%    N                  number of samples
%
% Output:
%     X                matrix [N x dim] with the quasti-random samples
%
% ------------------------------------------------------------------------
%
%
% Author : Flavio Cannavo'
% e-mail: flavio(dot)cannavo(at)gmail(dot)com
% Release: 1.0
% Date   : 07-02-2011
%
% History:
% 1.0  07-02-2011  First release.
%
%
%

%%

function X = fnc_getSobolSetMatlab(dim, N)

p = sobolset(dim);
p = scramble(p,'MatousekAffineOwen');
X = net(p,N);

repo: mcbench-benchmarks


text:
function s=get_name(a)

s=[get_name(a.algorithm)];
eval_name

s=[s ' m=' num2str(a.m)];
if a.n~=30          s=[s ' n='  num2str(a.n)];      end;

author: canlab

repo: CanlabCore

save_path: github-repos/MATLAB/canlab-CanlabCore


text:
function triangle01_monomial_integral_test ( )

%*****************************************************************************80
%
%% TRIANGLE01_MONOMIAL_INTEGRAL_TEST estimates integrals over the unit triangle.
%
%  Location:
%
%    http://people.sc.fsu.edu/~jburkardt/m_src/triangle_integrals/triangle01_monomial_integral_test.m
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    18 April 2015
%
%  Author:
%
%    John Burkardt
%
  fprintf ( 1, '\n' );
  fprintf ( 1, 'TRIANGLE01_MONOMIAL_INTEGRAL_TEST\n' );
  fprintf ( 1, '  TRIANGLE01_MONOMIAL_INTEGRAL returns the integral Q of\n' );
  fprintf ( 1, '  a monomial X^I Y^J over the interior of the unit triangle.\n' );

  fprintf ( 1, '\n' );
  fprintf ( 1, '   I   J         Q(I,J)\n' );

  for d = 0 : 5
    fprintf ( 1, '\n' );
    for i = 0 : d
      j = d - i;
      q = triangle01_monomial_integral ( i, j );
      fprintf ( 1, '  %2d  %2d  %14.6g\n', i, j, q );
    end
  end

  return
end