text:
function [ center, eccent, parity ] = tree_arc_center ( nnode, inode, jnode )

%*****************************************************************************80
%
%% TREE_ARC_CENTER computes the center, eccentricity, and parity of a tree.
%
%  Discussion:
%
%    A tree is an undirected graph of N nodes, which uses N-1 edges,
%    and is connected.
%
%    A graph with N-1 edges is not guaranteed to be a tree, and so this
%    routine must first check that condition before proceeding.
%
%    The edge distance between two nodes I and J is the minimum number of
%    edges that must be traversed in a path from I and J.
%
%    The eccentricity of a node I is the maximum edge distance between
%    node I and the other nodes J in the graph.
%
%    The radius of a graph is the minimum eccentricity over all nodes
%    in the graph.
%
%    The diameter of a graph is the maximum eccentricity over all nodes
%    in the graph.
%
%    The center of a graph is the set of nodes whose eccentricity is
%    equal to the radius, that is, the set of nodes of minimum eccentricity.
%
%    For a tree, the center is either a single node, or a pair of
%    neighbor nodes.
%
%    The parity of the tree is 1 if the center is a single node, or 2 if
%    the center is 2 nodes.
%
%    The center of a tree can be found by removing all "leaves", that is,
%    nodes of degree 1.  This step is repeated until only 1 or 2 nodes
%    are left.
%
%    Thanks to Alexander Sax for pointing out that a previous version of the
%    code was failing when the tree had an odd parity, that is, a single
%    center node, 15 April 2013.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    28 June 2013
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, integer NNODE, the number of nodes.
%
%    Input, integer INODE(NNODE-1), JNODE(NNODE-1), the edges of
%    the tree.  Edge I connects nodes INODE(I) and JNODE(I).
%
%    Output, integer CENTER(2).  CENTER(1) is the index of the
%    first node in the center.  CENTER(2) is 0 if there is only one node
%    in the center, or else the index of the second node.
%
%    Output, integer ECCENT, the eccentricity of the nodes in
%    the center, and the radius of the the tree.
%
%    Output, integer PARITY, the parity of the tree, which is
%    normally 1 or 2.
%
  eccent = 0;
  center(1) = 0;
  center(2) = 0;
  parity = 0;

  if ( nnode <= 0 )
    fprintf ( 1, '\n' );
    fprintf ( 1, 'TREE_ARC_CENTER - Fatal error!\n' );
    fprintf ( 1, '  NNODE <= 0.\n' );
    error ( 'TREE_ARC_CENTER - Fatal error!' );
  elseif ( nnode == 1 )
    eccent = 0;
    center(1) = 1;
    center(2) = 0;
    parity = 1;
    return
  elseif ( nnode == 2 )
    eccent = 1;
    center(1) = 1;
    center(2) = 2;
    parity = 2;
    return
  end
%
%  Is this graph really a tree?
%
  nedge = nnode - 1;
  result = graph_arc_is_tree ( nedge, inode, jnode, nnode );

  if ( result == 0 )
    fprintf ( 1, '\n' );
    fprintf ( 1, 'TREE_ARC_CENTER - Fatal error!\n' );
    fprintf ( 1, '  This graph is NOT a tree.\n' );
    error ( 'TREE_ARC_CENTER - Fatal error!' );
  end
%
%  Compute the degrees.
%
  degree = graph_arc_degree ( nnode, nedge, inode, jnode );
%
%  Defoliate the tree.
%
  nnode2 = nnode;

  while ( 1 )

    eccent = eccent + 1;
%
%  Find and mark the leaves.
%
    nleaf = 0;

    for i = 1 : nnode

      if ( degree(i) == 1 )
        nleaf = nleaf + 1;
        list(nleaf) = i;
      end

    end
%
%  Delete the leaves.
%
    for ileaf = 1 : nleaf

      i = list(ileaf);

      iedge = 0;
      j = 0;

      while ( 1 )

        iedge = iedge + 1;

        if ( nedge < iedge )
          fprintf ( 1, '\n' );
          fprintf ( 1, 'TREE_ARC_CENTER - Fatal error!\n' );
          fprintf ( 1, '  Data or algorithm failure.\n' );
          error ( 'TREE_ARC_CENTER - Fatal error!' );
        end

        if ( inode(iedge) == i )
          j = jnode(iedge);
          inode(iedge) = - inode(iedge);
          jnode(iedge) = - jnode(iedge);
        elseif ( jnode(iedge) == i )
          j = inode(iedge);
          inode(iedge) = - inode(iedge);
          jnode(iedge) = - jnode(iedge);
        end

        if ( j ~= 0 )
          break
        end

      end

      degree(i) = -1;
      nnode2 = nnode2 - 1;
      degree(j) = degree(j) - 1;
%
%  If the other node has degree 0, we must have just finished
%  stripping all leaves from the tree, leaving a single node.
%  Don't kill it here.  It is our odd center.
%
%     if ( degree(j) == 0 )
%       nnode2 = nnode2 - 1;
%     end

    end
%
%  Find the remaining nodes.
%
    nnode2 = 0;

    for i = 1 : nnode

      if ( 0 <= degree(i) )
        nnode2 = nnode2 + 1;
        list(nnode2) = i;
      end

    end
%
%  If at least 3, more pruning is required.
%
    if ( nnode2 < 3 )
      break
    end

  end
%
%  If only one or two nodes left, we are done.
%
  parity = nnode2;

  center(1:nnode2) = list(1:nnode2);
  inode(1:nedge) = abs ( inode(1:nedge) );
  jnode(1:nedge) = abs ( jnode(1:nedge) );

  return
end


text:
%-----------------------------------------------------------------------
% Job configuration created by cfg_util (rev $Rev: 4899 $)
%-----------------------------------------------------------------------
matlabbatch{1}.menu_cfg{1}.menu_entry{1}.conf_entry.type = 'cfg_entry';
matlabbatch{1}.menu_cfg{1}.menu_entry{1}.conf_entry.name = 'Output Variable Name';
matlabbatch{1}.menu_cfg{1}.menu_entry{1}.conf_entry.tag = 'name';
matlabbatch{1}.menu_cfg{1}.menu_entry{1}.conf_entry.strtype = 's';
matlabbatch{1}.menu_cfg{1}.menu_entry{1}.conf_entry.extras = [];
matlabbatch{1}.menu_cfg{1}.menu_entry{1}.conf_entry.num = [1 Inf];
matlabbatch{1}.menu_cfg{1}.menu_entry{1}.conf_entry.check = @cfg_check_assignin;
matlabbatch{1}.menu_cfg{1}.menu_entry{1}.conf_entry.help = {'Enter a valid MATLAB variable name. The contents of the input to "Output Item" will be assigned to a variable of this name in MATLAB workspace.'};
matlabbatch{1}.menu_cfg{1}.menu_entry{1}.conf_entry.def = [];
matlabbatch{2}.menu_cfg{1}.menu_entry{1}.conf_entry.type = 'cfg_entry';
matlabbatch{2}.menu_cfg{1}.menu_entry{1}.conf_entry.name = 'Output Item';
matlabbatch{2}.menu_cfg{1}.menu_entry{1}.conf_entry.tag = 'output';
matlabbatch{2}.menu_cfg{1}.menu_entry{1}.conf_entry.strtype = 'e';
matlabbatch{2}.menu_cfg{1}.menu_entry{1}.conf_entry.extras = [];
matlabbatch{2}.menu_cfg{1}.menu_entry{1}.conf_entry.num = [];
matlabbatch{2}.menu_cfg{1}.menu_entry{1}.conf_entry.check = [];
matlabbatch{2}.menu_cfg{1}.menu_entry{1}.conf_entry.help = {'The contents that is passed to this input that will be assigned to the workspace variable whose name is given in "Output Variable Name".'};
matlabbatch{2}.menu_cfg{1}.menu_entry{1}.conf_entry.def = [];
matlabbatch{3}.menu_cfg{1}.menu_struct{1}.conf_exbranch.type = 'cfg_exbranch';
matlabbatch{3}.menu_cfg{1}.menu_struct{1}.conf_exbranch.name = 'Pass Output to Workspace';
matlabbatch{3}.menu_cfg{1}.menu_struct{1}.conf_exbranch.tag = 'cfg_assignin';
matlabbatch{3}.menu_cfg{1}.menu_struct{1}.conf_exbranch.val{1}(1) = cfg_dep;
matlabbatch{3}.menu_cfg{1}.menu_struct{1}.conf_exbranch.val{1}(1).tname = 'Val Item';
matlabbatch{3}.menu_cfg{1}.menu_struct{1}.conf_exbranch.val{1}(1).tgt_spec = {};
matlabbatch{3}.menu_cfg{1}.menu_struct{1}.conf_exbranch.val{1}(1).sname = 'Output Variable Name (cfg_entry)';
matlabbatch{3}.menu_cfg{1}.menu_struct{1}.conf_exbranch.val{1}(1).src_exbranch = substruct('.','val', '{}',{1}, '.','val', '{}',{1}, '.','val', '{}',{1});
matlabbatch{3}.menu_cfg{1}.menu_struct{1}.conf_exbranch.val{1}(1).src_output = substruct('()',{1});
matlabbatch{3}.menu_cfg{1}.menu_struct{1}.conf_exbranch.val{2}(1) = cfg_dep;
matlabbatch{3}.menu_cfg{1}.menu_struct{1}.conf_exbranch.val{2}(1).tname = 'Val Item';
matlabbatch{3}.menu_cfg{1}.menu_struct{1}.conf_exbranch.val{2}(1).tgt_spec = {};
matlabbatch{3}.menu_cfg{1}.menu_struct{1}.conf_exbranch.val{2}(1).sname = 'Output Item (cfg_entry)';
matlabbatch{3}.menu_cfg{1}.menu_struct{1}.conf_exbranch.val{2}(1).src_exbranch = substruct('.','val', '{}',{2}, '.','val', '{}',{1}, '.','val', '{}',{1});
matlabbatch{3}.menu_cfg{1}.menu_struct{1}.conf_exbranch.val{2}(1).src_output = substruct('()',{1});
matlabbatch{3}.menu_cfg{1}.menu_struct{1}.conf_exbranch.prog = @cfg_run_assignin;
matlabbatch{3}.menu_cfg{1}.menu_struct{1}.conf_exbranch.vout = [];
matlabbatch{3}.menu_cfg{1}.menu_struct{1}.conf_exbranch.check = [];
matlabbatch{3}.menu_cfg{1}.menu_struct{1}.conf_exbranch.help = {
                                                                'Assign a computation result to a workspace variable.'
                                                                'The value entered into "Output Item" will be assigned to a MATLAB workspace variable whose name is specified in "Output Variable Name". If this variable already exists, a new variable name will be generated. This can be useful to assess the results of computations with other MATLAB routines or for debugging.'
                                                                }';

author: spm


text:
function ind = VOChash_lookup_HH(hash,s)

hsize=numel(hash.key);
h=mod(str2double(s([5:7 9:end])),hsize)+1;
ind=hash.val{h}(strmatch(s,hash.key{h},'exact'));

author: aosokin


text:
% UTILS_OPTCREATE
% Create a structure of options given fieldName - value pairs in args as:
% args = {'field1', value1, 'field2', value,...};
% This function can be used as a template in different implementations.
% This can be done by:
% a) providing here an extra field {'optionsType', your_options_type}
% b) having an implementation of a function:
%  ['defaults_' your_options_type]
% which returns default values for this type of options
%
% COPYRIGHT: Andreas C. Damianou, 2012, 2013
% VARGPLVM
function options = util_optCreate(args1, args2)

options = [];

if nargin < 1
    return
end

if isstruct(args1)
    % In this case, the struct is initialised with another struct given in
    % arg1. Then, args2 (if present) can just alter specific fields or add.
    allfields = fieldnames(args1);
    for i=1:length(allfields)
        options.(allfields{i}) = args1.(allfields{i});
    end
    if nargin > 1
        for i = 1:2:length(args2)
            options.(args2{i}) = args2{i+1};
        end
    end
else
    for i = 1:2:length(args1)
        options.(args1{i}) = args1{i+1};
    end
end

%--- This can change depending on the implementation needs ----------

%-- Set with default values whatever is not already set
% Default values
if isfield(options, 'optionsType')
    fhandle = str2func(['defaults_' options.optionsType]);
    defaults = fhandle();
else
    defaults = {};
end
%----------------------------------------------------------------------

for i = 1:2:length(defaults)
    if ~isfield(options, defaults{i})
        options.(defaults{i}) = defaults{i+1};
    end
end

author: SheffieldML

repo: vargplvm


text:
%----UPDATE PLOT DATA/STYLE----%

function obj = updateData(obj, dataIndex)

try

    %-update plot based on TreatAs PlotOpts-%

    if ismember('pie3', lower(obj.PlotOptions.TreatAs))
        updatePie3(obj, dataIndex);
    elseif ismember('pcolor', lower(obj.PlotOptions.TreatAs))
        updatePColor(obj, dataIndex);
    elseif ismember('contour3', lower(obj.PlotOptions.TreatAs))
        updateContour3(obj, dataIndex);
    elseif ismember('ezpolar', lower(obj.PlotOptions.TreatAs))
        updateLineseries(obj, dataIndex);
    elseif ismember('polarhistogram', lower(obj.PlotOptions.TreatAs))
        updateHistogramPolar(obj, dataIndex);
    elseif ismember('coneplot', lower(obj.PlotOptions.TreatAs))
        updateConeplot(obj, dataIndex);
    elseif ismember('bar3', lower(obj.PlotOptions.TreatAs))
        updateBar3(obj, dataIndex);
    elseif ismember('bar3h', lower(obj.PlotOptions.TreatAs))
        updateBar3h(obj, dataIndex);
    elseif ismember('fmesh', lower(obj.PlotOptions.TreatAs))
        updateFmesh(obj, dataIndex);
    elseif ismember('surfc', lower(obj.PlotOptions.TreatAs))
        updateSurfc(obj, dataIndex);
    elseif ismember('meshc', lower(obj.PlotOptions.TreatAs))
        updateSurfc(obj, dataIndex);
    elseif ismember('surfl', lower(obj.PlotOptions.TreatAs))
        updateSurfl(obj, dataIndex);

    %-update plot based on plot call class-%

    else

        switch lower(obj.State.Plot(dataIndex).Class)

            %--SPIDER PLOT -> SPECIAL CASE--%
            case 'spider_plot_class'
                updateSpiderPlot(obj, dataIndex);

            %--GEOAXES -> SPECIAL CASE--%
            case 'geoaxes'
                UpdateGeoAxes(obj, dataIndex);

            %-EMULATE AXES -> SPECIAL CASE--%
            case 'nothing'
                updateOnlyAxes(obj, dataIndex);

            %--CORE PLOT OBJECTS--%
            case 'geobubble'
                updateGeobubble(obj, dataIndex);
            case 'scatterhistogram'
                updateScatterhistogram(obj, dataIndex);
            case 'wordcloud'
                updateWordcloud(obj, dataIndex);
            case 'heatmap'
                updateHeatmap(obj, dataIndex);
            case 'image'
                if ~obj.PlotOptions.Image3D
                    updateImage(obj, dataIndex);
                else
                    updateImage3D(obj, dataIndex);
                end
            case 'line'
                if obj.PlotlyDefaults.isGeoaxis
                    updateGeoPlot(obj, dataIndex);
                elseif ismember('polarplot', lower(obj.PlotOptions.TreatAs))
                    updatePolarplot(obj, dataIndex);
                elseif ismember('ternplot', lower(obj.PlotOptions.TreatAs))
                    updateTernaryPlot(obj, dataIndex);
                else
                    updateLineseries(obj, dataIndex);
                end
            case 'categoricalhistogram'
                updateCategoricalHistogram(obj, dataIndex);
            case 'histogram'
                updateHistogram(obj, dataIndex);
            case 'histogram2'
                updateHistogram2(obj, dataIndex);
            case 'patch'
                % check for histogram
                if isHistogram(obj,dataIndex)
                    updateHistogram(obj,dataIndex);
                elseif ismember('ternplotpro', lower(obj.PlotOptions.TreatAs))
                    updateTernaryPlotPro(obj, dataIndex);
                elseif ismember('ternpcolor', lower(obj.PlotOptions.TreatAs))
                    updateTernaryPlotPro(obj, dataIndex);
                elseif ismember('isosurface', lower(obj.PlotOptions.TreatAs))
                    updateIsosurface(obj, dataIndex);
                else
                    updatePatch(obj, dataIndex);
                end
            case 'rectangle'
                updateRectangle(obj,dataIndex);
            case 'surface'
                if ismember('surf', lower(obj.PlotOptions.TreatAs))
                    updateSurf(obj, dataIndex);
                elseif ismember('mesh', lower(obj.PlotOptions.TreatAs))
                    updateMesh(obj, dataIndex);
                elseif ismember('slice', lower(obj.PlotOptions.TreatAs))
                    updateSlice(obj, dataIndex);
                else
                    updateSurfaceplot(obj,dataIndex);
                end
            case {'functionsurface', 'parameterizedfunctionsurface'}
                updateFunctionSurface(obj,dataIndex);
            case 'implicitfunctionsurface'
                updateImplicitFunctionSurface(obj,dataIndex);

                %-GROUP PLOT OBJECTS-%
            case 'area'
                updateArea(obj, dataIndex);
            case 'areaseries'
                updateAreaseries(obj, dataIndex);
            case 'animatedline'
                updateAnimatedLine(obj, dataIndex);
            case 'bar'
                updateBar(obj, dataIndex);
            case 'barseries'
                updateBarseries(obj, dataIndex);
            case 'baseline'
                updateBaseline(obj, dataIndex);
            case {'contourgroup','contour'}
                if obj.PlotOptions.ContourProjection
                    updateContourProjection(obj,dataIndex);
                elseif ismember('terncontour', lower(obj.PlotOptions.TreatAs))
                    updateTernaryContour(obj, dataIndex);
                else
                    updateContourgroup(obj,dataIndex);
                end
            case 'functioncontour'
                updateFunctionContour(obj,dataIndex);
            case 'errorbar'
                updateErrorbar(obj,dataIndex);
            case 'errorbarseries'
                updateErrorbarseries(obj,dataIndex);
            case 'lineseries'
                updateLineseries(obj, dataIndex);
            case 'quiver'
                updateQuiver(obj, dataIndex);
            case 'quivergroup'
                updateQuivergroup(obj, dataIndex);
            case 'scatter'
                if ismember('scatterpolar', lower(obj.PlotOptions.TreatAs))
                    updateScatterPolar(obj, dataIndex);
                elseif obj.PlotlyDefaults.isGeoaxis
                    updateGeoScatter(obj, dataIndex);
                else
                    updateScatter(obj, dataIndex);
                end
            case 'scattergroup'
                updateScattergroup(obj, dataIndex);
            case 'stair'
                updateStair(obj, dataIndex);
            case 'stairseries'
                updateStairseries(obj, dataIndex);
            case 'stackedplot'
                updateStackedplot(obj, dataIndex);
            case 'stem'
                updateStem(obj, dataIndex);
            case 'stemseries'
                updateStemseries(obj, dataIndex);
            case 'surfaceplot'
                updateSurfaceplot(obj,dataIndex);
            case 'implicitfunctionline'
                updateLineseries(obj, dataIndex);

                %--Plotly supported MATLAB group plot objects--%
            case {'hggroup','group'}
                % check for boxplot
                if isBoxplot(obj, dataIndex)
                    updateBoxplot(obj, dataIndex);
                end
        end
    end

catch exception
    if obj.UserData.Verbose
        fprintf([exception.message '\nWe had trouble parsing the ' obj.State.Plot(dataIndex).Class ' object.\n',...
                 'This trace might not render properly.\n\n']);
    end
end

%------------------------AXIS/DATA CLEAN UP-------------------------------%

try
    %-AXIS INDEX-%
    axIndex = obj.getAxisIndex(obj.State.Plot(dataIndex).AssociatedAxis);

    %-CHECK FOR MULTIPLE AXES-%
    [xsource, ysource] = findSourceAxis(obj,axIndex);

    %-AXIS DATA-%
    eval(['xaxis = obj.layout.xaxis' num2str(xsource) ';']);
    eval(['yaxis = obj.layout.yaxis' num2str(ysource) ';']);

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

    % check for xaxis dates
    if strcmpi(xaxis.type, 'date')
        obj.data{dataIndex}.x =  convertDate(obj.data{dataIndex}.x);
    end

    % check for xaxis categories
    if strcmpi(xaxis.type, 'category') && ...
            ~strcmp(obj.data{dataIndex}.type,'box')
        obj.data{dataIndex}.x =  get(obj.State.Plot(dataIndex).AssociatedAxis,'XTickLabel');
        eval(['obj.layout.xaxis' num2str(xsource) '.autotick=true;']);
    end

    % check for yaxis dates
    if strcmpi(yaxis.type, 'date')
        obj.data{dataIndex}.y =  convertDate(obj.data{dataIndex}.y);
    end

    % check for yaxis categories
    if strcmpi(yaxis.type, 'category') && ...
            ~strcmp(obj.data{dataIndex}.type,'box')
        obj.data{dataIndex}.y =  get(obj.State.Plot(dataIndex).AssociatedAxis,'YTickLabel');
        eval(['obj.layout.yaxis' num2str(xsource) '.autotick=true;']);
    end
catch
    % TODO to the future
    % disp('catch at line 157 in updateData.m file')
end

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

try
    if obj.layout.isAnimation
        %- Play Button Options-%
        opts{1} = nan;
        opts{2}.frame.duration = obj.PlotOptions.FrameDuration;
        opts{2}.frame.redraw = true;
        opts{2}.fromcurrent = true;
        opts{2}.mode = 'immediate';
        opts{2}.transition.duration = obj.PlotOptions.FrameTransitionDuration;

        button{1}.label = '▶';
        button{1}.method = 'animate';
        button{1}.args = opts;

        opts{1} = {nan};
        opts{2}.transition.duration = 0;
        opts{2}.frame.duration = 0;

        button{2}.label = '◼';
        button{2}.method = 'animate';
        button{2}.args = opts;

        obj.layout.updatemenus{1}.type = 'buttons';
        obj.layout.updatemenus{1}.buttons = button;
        obj.layout.updatemenus{1}.pad.r = 70;
        obj.layout.updatemenus{1}.pad.t = 10;
        obj.layout.updatemenus{1}.direction = 'left';
        obj.layout.updatemenus{1}.showactive = true;
        obj.layout.updatemenus{1}.x = 0.01;
        obj.layout.updatemenus{1}.y = 0.01;
        obj.layout.updatemenus{1}.xanchor = 'left';
        obj.layout.updatemenus{1}.yanchor = 'top';

        obj.layout = rmfield(obj.layout,'isAnimation');
    end
catch
end

end

author: plotly


text:
function item = setval(item, val, dflag)

% function item = setval(item, val, dflag)
% set item.val{1} to val. Validity checks are performed through subsasgn.
% If val == {}, set item.val to {}.
% If dflag is true, and item.def is not empty, set the default setting for
% this item instead by calling feval(item.def{:}, val). If val == {}, use
% the string '<UNDEFINED>' as in a harvested tree. If dflag is true, but
% no item.def defined, set item.val{1} instead.
%
% This code is part of a batch job configuration system for MATLAB. See
%      help matlabbatch
% for a general overview.
%_______________________________________________________________________
% Copyright (C) 2007 Freiburg Brain Imaging

% Volkmar Glauche
% $Id: setval.m 5678 2013-10-11 14:58:04Z volkmar $

rev = '$Rev: 5678 $'; %#ok

if iscell(val) && isempty(val)
    if dflag
        if ~isempty(item.def)
            try
                feval(item.def, {'<UNDEFINED>'});
            catch
                cfg_message('matlabbatch:setval:defaults', ...
                            '%s: unable to set default value.', ...
                            subsasgn_checkstr(item, substruct('.','val')));
            end
        else
            item = subsasgn(item, substruct('.','val'), {});
        end
    else
        item = subsasgn(item, substruct('.','val'), {});
    end
else
    if dflag
        [sts, val1] = subsasgn_check(item, substruct('.','val'), {val});
        if sts
            if ~isempty(item.def)
                try
                    feval(item.def, val1);
                catch
                    cfg_message('matlabbatch:setval:defaults', ...
                                '%s: unable to set default value.', ...
                                subsasgn_checkstr(item, substruct('.','val')));
                end
            else
                item = subsasgn(item, substruct('.','val'), {val});
            end
        end
    else
        item = subsasgn(item, substruct('.','val'), {val});
    end
end

author: fieldtrip

repo: fieldtrip

save_path: github-repos/MATLAB/fieldtrip-fieldtrip


text:
function [source] = ft_source2full(source)

% FT_SOURCE2FULL recreates the grid locations outside the brain in the source
% reconstruction, so that the source volume again describes the full grid.
% This undoes the memory savings that can be achieved using FT_SOURCE2SPARSE
% and makes it possible again to plot  the source volume and save it to an
% external file.
%
% Use as
%   [source] = ft_source2full(source)
%
% See also FT_SOURCE2SPARSE

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

ft_defaults

if ~isfield(source, 'inside')  || ...
   ~isfield(source, 'outside') || ...
   ~isfield(source, 'dim')
  ft_error('one of the required fields is missing in the source structure');
end

if ~isfield(source, 'pos') && (~isfield(source, 'xgrid') || ~isfield(source, 'ygrid') || ...
                               ~isfield(source, 'zgrid'))
  ft_error('the input data needs at least a ''pos'' field, or ''x/y/zgrid''');
end

if isfield(source, 'xgrid'),
  xgrid = source.xgrid;
  ygrid = source.ygrid;
  zgrid = source.zgrid;
  sparsepos = source.pos;

  % recreate the positions of the dipole grid
  [X, Y, Z] = ndgrid(xgrid, ygrid, zgrid);
  pos = [X(:) Y(:) Z(:)];
else
  %FIXME this assumes that the voxel data are ordered as if in a regularly spaced 3D grid,
  %but with only the inside voxels present
  ft_warning('assuming the voxel data to be ordered as if in a regularly spaced 3D grid');
  xgrid = 1:source.dim(1);
  ygrid = 1:source.dim(2);
  zgrid = 1:source.dim(3);

  %establish a homogeneous transformation matrix from voxels to headspace based on the sparse positions
  sparsepos = source.pos;
  ok  = 0;
  cnt = 0;
  while ok==0,
    cnt  = cnt+1;
    dpos = sparsepos - sparsepos(cnt*ones(size(sparsepos,1),1),:);
    [srt, indx] = sort(sum(dpos.^2,2));
    srt    = dpos(indx,:);
    tmpsrt = abs(srt(2:7,:));
    csrt   = tmpsrt*tmpsrt';
    sel    = find(sum(csrt==0)>=2);
    if numel(sel)>=3,
      ok = 1;
    end
  end
  tmppos  = sparsepos(indx([1 sel(:)'+1]),:);
  tmpdpos = dpos(indx([1 sel(:)'+1]),:);

  % FIXME the following is a bit experimental and not fully tested yet it works in general case rotation
  M         = pinv(tmpdpos(2:4,:));

  % get rotation such that maxima are on diagonal and positive
  m(1) = find(M(1,:)==max(abs(M(1,:))));
  m(2) = find(M(2,:)==max(abs(M(2,:))));
  m(3) = find(M(3,:)==max(abs(M(3,:))));
  [srt, indx] = sort(m);
  M    = M(indx,:);
  M    = M*diag(sign(diag(M)));
  sparsepos = sparsepos*M;

  % translation
  T         = -min(sparsepos,[],1)+1;
  sparsepos = sparsepos + T(ones(size(sparsepos,1),1), :);

  % recreate the positions of the dipole grid
  [X, Y, Z] = ndgrid(xgrid, ygrid, zgrid);
  pos = [X(:) Y(:) Z(:)];
  pos = ft_warp_apply(inv([M T(:);0 0 0 1]), pos);
end

Nsparse = length(source.inside);
siz     = source.dim;
Nfull   = prod(siz);

% determine the size that each slice takes in memory
sx = 1;
sy = siz(1);
sz = siz(1) * siz(2);

if isfield(source, 'inside') && isfield(source, 'outside') && size(source.pos,1)==Nfull
  % it contains all source positions
  inside = source.inside;
  outside = source.outside;
else
  % it only contains the inside source positions, which are all inside the brain
  % reconstruct the original inside and outside grid locations
  inside = zeros(Nsparse,1);
  for i=1:Nsparse
    fx = find(xgrid==sparsepos(i,1));
    fy = find(ygrid==sparsepos(i,2));
    fz = find(zgrid==sparsepos(i,3));
      inside(i) = (fx-1)*sx + (fy-1)*sy + (fz-1)*sz + 1;
  end
  outside = setdiff([1:Nfull]', inside);
end

fprintf('total number of dipoles        : %d\n', length(inside)+length(outside));
fprintf('number of dipoles inside  brain: %d\n', length(inside));
fprintf('number of dipoles outside brain: %d\n', length(outside));

% determine whether the source is old or new style
fnames = fieldnames(source);
if any(~cellfun('isempty', strfind(fnames, 'dimord'))),
  stype = 'new';
else
  stype = 'old';
end

if strcmp(stype, 'old')
  % original code
  % first do the non-trial fields
  source.dim = [1 length(inside) 1]; %to fool parameterselection
  [param]    = parameterselection('all', source);
  trlparam   = strmatch('trial', param);
  sel        = setdiff(1:length(param), trlparam);
  ind=find(ismember(param,'inside')); % find the index of 'inside' field
  % because its position varies with isfield('plvspctrm') vs. 'cohspctrm'
  param      = param(sel(ind));

  for j = 1:length(param)
    dat = getsubfield(source, param{j});
    if islogical(dat)
      tmp         = false(1,Nfull);
      tmp(inside) = dat;
    elseif iscell(dat)
      tmp          = cell(1,Nfull);
      tmp(inside)  = dat;
      %tmp(outside) = nan;
    else
      tmp         = nan(1,Nfull);
      tmp(inside) = dat;
    end
    source = setsubfield(source, param{j}, tmp);
  end

  % then do the trial fields
  if     isfield(source, 'trial' )
    for j = 1:length(source.trial)
      tmpsource     = source.trial(j);
      tmpsource.dim = source.dim; % to fool parameterselection
      tmpparam      = parameterselection('all', tmpsource);
      for k = 1:length(tmpparam)
        dat = getsubfield(tmpsource, tmpparam{k});
        if islogical(dat)
          tmp         = false(1,Nfull);
          tmp(inside) = dat;
        elseif iscell(dat)
          tmp          = cell(1,Nfull);
          tmp(inside)  = dat;
          %tmp(outside) = nan;
        else
          tmp         = nan(1,Nfull);
          tmp(inside) = dat;
        end
        tmpsource = setsubfield(tmpsource, tmpparam{k}, tmp);
      end
      tmpsource       = rmfield(tmpsource, 'dim');
      source.trial(j) = tmpsource;
    end
  elseif isfield(source, 'trialA')
    for j = 1:length(source.trialA)
      tmpsource     = source.trialA(j);
      tmpsource.dim = source.dim; % to fool parameterselection
      tmpparam      = parameterselection('all', tmpsource);
      for k = 1:length(tmpparam)
        dat = getsubfield(tmpsource, tmpparam{k});
        if islogical(dat)
          tmp         = false(1,Nfull);
          tmp(inside) = dat;
        elseif iscell(dat)
          tmp          = cell(1,Nfull);
          tmp(inside)  = dat;
          %tmp(outside) = nan;
        else
          tmp         = nan(1,Nfull);
          tmp(inside) = dat;
        end
        tmpsource = setsubfield(tmpsource, tmpparam{k}, tmp);
      end
      tmpsource        = rmfield(tmpsource, 'dim');
      source.trialA(j) = tmpsource;
    end
  elseif isfield(source, 'trialB')
    for j = 1:length(source.trialB)
      tmpsource     = source.trialB(j);
      tmpsource.dim = source.dim; % to fool parameterselection
      tmpparam      = parameterselection('all', tmpsource);
      for k = 1:length(tmpparam)
        dat = getsubfield(tmpsource, tmpparam{k});
        if islogical(dat)
          tmp         = false(1,Nfull);
          tmp(inside) = dat;
        elseif iscell(dat)
          tmp          = cell(1,Nfull);
          tmp(inside)  = dat;
          %tmp(outside) = nan;
        else
          tmp         = nan(1,Nfull);
          tmp(inside) = dat;
        end
        tmpsource = setsubfield(tmpsource, tmpparam{k}, tmp);
      end
      tmpsource        = rmfield(tmpsource, 'dim');
      source.trialB(j) = tmpsource;
    end
  end

  % and finally do the coherence-like matrices (size Nvox X Nvox)
  fn = fieldnames(source);
  for i=1:length(fn)
    d = getfield(source, fn{i});
    m = size(d, 1);
    n = size(d, 2);
    if m==Nsparse && n==Nsparse
      tmp = nan(Nfull,Nfull);
      tmp(inside,inside) = d;
      source = setfield(source, fn{i}, tmp);
    end
  end

  % update the inside and outside definitions
  source.inside  = inside;
  source.outside = outside;
  source.pos     = pos;
  source.dim     = siz;
elseif strcmp(stype, 'new')
  % new style conversion
  fn = fieldnames(source);
  for i=1:numel(fn)
    if any(size(source.(fn{i}))==Nsparse)
      if iscell(source.(fn{i}))
        indx = find(size(source.(fn{i}))==Nsparse);
        if all(indx==1)
          tmp            = cell(Nfull,1);
          tmp(inside,1)  = source.(fn{i});
          source.(fn{i}) = tmp;
        elseif all(indx==2)
          tmp            = cell(1,Nfull);
          tmp(1,inside)  = source.(fn{i});
          source.(fn{i}) = tmp;
        else
          ft_warning('sparse to full conversion failed for field %s\n', fn{i});
        end
      else
        indx = find(size(source.(fn{i}))==Nsparse);
        if all(indx==1)
          tmpsiz = [size(source.(fn{i})) 1];
          tmp    = nan([Nfull tmpsiz(2:end)]);
          tmp(inside,:,:,:,:) = source.(fn{i});
        elseif all(indx==2)
          tmpsiz = [size(source.(fn{i})) 1];
          tmp    = nan([tmpsiz(1) Nfull tmpsiz(3:end)]);
          tmp(:,inside,:,:,:) = source.(fn{i});
        elseif all(indx==[1 2])
          % bivariate matrix
          tmpsiz = [size(source.(fn{i})) 1];
          tmp    = nan([Nfull Nfull tmpsiz(3:end)]);
          tmp(inside,inside,:,:,:) = source.(fn{i});
        else
          ft_warning('sparse to full conversion failed for field %s\n', fn{i});
        end
      end
      % nothing to do
    end
  end

  % update the inside and outside definitions and pos
  source.inside  = inside;
  source.outside = outside;
  source.pos     = pos;

end
cfg = [];
% add version information to the configuration
try
  % get the full name of the function
  cfg.version.name = mfilename('fullpath');
catch
  % required for compatibility with MATLAB versions prior to release 13 (6.5)
  [st, i] = dbstack;
  cfg.version.name = st(i);
end
cfg.version.id = '$Id$';
% remember the configuration details of the input data
try, cfg.previous = source.cfg; end
% remember the exact configuration details in the output
source.cfg = cfg;


text:
function locate_examples()
    % Changes current directory to the folder containing spx.
    root = spx.root_dir;
    lib = fileparts(root);
    examples = fullfile(lib, 'examples');
    cd (examples);
end

author: indigits

repo: sparse-plex


text:
function plot_cylinder_model(cylinder,Color,fig,nf,alp,Ind)

% ---------------------------------------------------------------------
% PLOT_CYLINDER_MODEL.M       Plots the given cylinder model
%
% Version 1.2.0
% Latest update     3 Aug 2021
%
% Copyright (C) 2013-2021 Pasi Raumonen
% ---------------------------------------------------------------------

% Plots the cylinder model.
% cylinder  Structure array containin the cylinder info
%               (radius, length, start, axis, BranchOrder)
% fig       Figure number
% nf        Number of facets in the cyliders (in the thickest cylinder,
%               scales down with radius to 4 which is the minimum)
% alp       Alpha value (1 = no trasparency, 0 = complete transparency)
% Color     If equals to "order", colors the cylinders based on branching
%               order, otherwise colors each branch with unique color
% Ind       Indexes of cylinders to be plotted from a subset of cylinders
%               (Optional, if not given then all cylinders are plotted)

% Changes from version 1.1.0 to 1.2.0, 3 Aug 2021:
% 1) Changed the surface plot ("patch") so that the edges are not plotted
%    with separate color, so the surface looks more smooth. Similarly added
%    shading. (These are added at the end of the file)
% 2) Added cylinder branch "Bran" and branch order "BOrd" vectors where the
%    coloring options are defined to prevent some errors

% Changes from version 1.0.0 to 1.1.0, 13 July 2020:
% 1) Added option for choosing the coloring based either on branch order or
%    unique color for each branch
% 2) Removed the possibility of the input "cylinder" being a matrix
% 3) Added default values for inputs

n = nargin;
if n < 5
  alp = 1;
  if n < 4
    nf = 20;
    if n < 3
      fig = 1;
      if n == 1
        Color = 'order';
      end
    end
  end
end

Rad = cylinder.radius;
Len = cylinder.length;
Sta = cylinder.start;
%Sta = Sta-Sta(1,:);
Axe = cylinder.axis;
if strcmp(Color,'order')
  BOrd = cylinder.BranchOrder;
  Bran = cylinder.branch;
end
if strcmp(Color,'branch')
  Bran = cylinder.branch;
  BOrd = cylinder.BranchOrder;
end

if nargin == 6
  Rad = Rad(Ind);
  Len = Len(Ind);
  Sta = Sta(Ind,:);
  Axe = Axe(Ind,:);
  BOrd = BOrd(Ind);
  if strcmp(Color,'branch')
    Bran = Bran(Ind);
  end
end

nc = size(Rad,1); % Number of cylinder

if strcmp(Color,'order')
  Color = 1;
  % Color the cylinders in branches based on the branch order
  col = [
    0.00  0.00  1.00
    0.00  0.50  0.00
    1.00  0.00  0.00
    0.00  0.75  0.75
    0.75  0.00  0.75
    0.75  0.75  0.00
    0.25  0.25  0.25
    0.75  0.25  0.25
    0.95  0.95  0.00
    0.25  0.25  0.75
    0.75  0.75  0.75
    0.00  1.00  0.00
    0.76  0.57  0.17
    0.54  0.63  0.22
    0.34  0.57  0.92
    1.00  0.10  0.60
    0.88  0.75  0.73
    0.10  0.49  0.47
    0.66  0.34  0.65
    0.99  0.41  0.23];
  col = repmat(col,[10,1]);
elseif strcmp(Color,'branch')
  Color = 0;
  % Color the cylinders in branches with an unique color of each branch
  N = double(max(Bran));
  col = rand(N,3);
  Par = cylinder.parent;
  for i = 2:nc
    if Par(i) > 0 && Bran(Par(i)) ~= Bran(i)
      C = col(Bran(Par(i)),:);
      c = col(Bran(i),:);
      while sum(abs(C-c)) < 0.2
        c = rand(1,3);
      end
      col(Bran(i),:) = c;
    end
  end
end

Cir = cell(nf,2);
for i = 4:nf
  B = [cos((1/i:1/i:1)*2*pi)' sin((1/i:1/i:1)*2*pi)' zeros(i,1)];
  T = [cos((1/i:1/i:1)*2*pi)' sin((1/i:1/i:1)*2*pi)' ones(i,1)];
  Cir{i,1} = [B; T];
  Cir{i,2} = [(1:1:i)' (i+1:1:2*i)' [(i+2:1:2*i)'; i+1] [(2:1:i)'; 1]];
end

Vert = zeros(2*nc*(nf+1),3);
Facets = zeros(nc*(nf+1),4);
fvd = zeros(nc*(nf+1),3);
t = 1;
f = 1;

% Scale, rotate and translate the standard cylinders
for i = 1:nc
  n = ceil(sqrt(Rad(i)/Rad(1))*nf);
  n = min(n,nf);
  n = max(n,4);
  C = Cir{n,1};
  % Scale
  C(:,1:2) = Rad(i)*C(:,1:2);
  C(n+1:end,3) = Len(i)*C(n+1:end,3);
  % Rotate
  ang = real(acos(Axe(i,3)));
  Axis = cross([0 0 1]',Axe(i,:)');
  Rot = rotation_matrix(Axis,ang);
  C = (Rot*C')';
  % Translate
  C = mat_vec_subtraction(C,-Sta(i,:));
  Vert(t:t+2*n-1,:) = C;
  Facets(f:f+n-1,:) = Cir{n,2}+t-1;
  if Color == 1
    fvd(f:f+n-1,:) = repmat(col(BOrd(i)+1,:),[n 1]);
  else
    fvd(f:f+n-1,:) = repmat(col(Bran(i),:),[n 1]);
  end
  t = t+2*n;
  f = f+n;
end
t = t-1;
f = f-1;
Vert = Vert(1:t,:);
Facets = Facets(1:f,:);
fvd = fvd(1:f,:);

figure(fig)
plot3(Vert(1,1),Vert(1,2),Vert(1,3))
patch('Vertices',Vert,'Faces',Facets,'FaceVertexCData',fvd,'FaceColor','flat')
alpha(alp)
axis equal
grid on
view(-37.5,30)

shading flat
lightangle(gca,-45,30)
lighting gouraud

author: InverseTampere

repo: TreeQSM


text:
% Fixed budget quantized kernel least-mean-square algorithm
%
% S. Zhao, B. Chen, P. Zhu, J. C. Principe, "Fixed budget quantized kernel
% least-mean-square algorithm", Signal Processing, Volume 93, Issue 9,
% September 2013, Pages 2759-2770,
% http://dx.doi.org/10.1016/j.sigpro.2013.02.012.
%
% Remark: significance calculation only implemented for Gaussian kernel.
%
% This file is part of the Kernel Adaptive Filtering Toolbox for Matlab.
% https://github.com/steven2358/kafbox/

classdef qklms_fb < kernel_adaptive_filter

    properties (GetAccess = 'public', SetAccess = 'private')
        eta = .9; % learning rate
        epsu = .1; % quantization threshold
        beta = .95; % forgetting factor for influence
        M = 500; % dictionary size (K in publication)
        kerneltype = 'gauss'; % kernel type
        kernelpar = 1; % kernel parameter
    end

    properties (GetAccess = 'public', SetAccess = 'private')
        dict = []; % codebook
        alpha = []; % expansion coefficients
        E = []; % significance vector (~ importance)
        lambda = []; % influence vector (~ how many samples are quantized
        % to each centre)
    end

    methods
        function kaf = qklms_fb(parameters) % constructor
            if (nargin > 0) % copy valid parameters
                for fn = fieldnames(parameters)'
                    if ismember(fn,fieldnames(kaf))
                        kaf.(fn{1}) = parameters.(fn{1});
                    end
                end
            end
        end

        function y_est = evaluate(kaf,x) % evaluate the algorithm
            if size(kaf.dict,1)>0
                k = kernel(kaf.dict,x,kaf.kerneltype,kaf.kernelpar);
                y_est = k'*kaf.alpha;
            else
                y_est = zeros(size(x,1),1);
            end
        end

        function train(kaf,x,y) % train the algorithm
            y_est = kaf.evaluate(x);
            err = y - y_est;

            m = size(kaf.dict,1);
            if m==0
                d2 = kaf.epsu^2 + 1; % force dictionary growth
            else
                [d2,j] = min(sum((kaf.dict - repmat(x,m,1)).^2,2));
            end

            c1 = pi*kaf.kernelpar^2/2; % pre-calculate

            if d2 <= kaf.epsu^2 % new basis under quantization threshold
                kaf.alpha(j) = kaf.alpha(j) + kaf.eta*err;

                % update significance for quantizing y to j-th centre (16)
                if m>1
                    inds = 1:m;
                    inds(j) = [];
                    kaf.E(inds) = kaf.beta*kaf.E(inds) + c1*abs(kaf.alpha(inds)).* ...
                        kernel(kaf.dict(inds,:),kaf.dict(j,:),kaf.kerneltype,kaf.kernelpar);
                end
                kaf.E(j) = abs(1 + kaf.eta*err/kaf.alpha(j)) * ...
                    kaf.beta*kaf.E(j) + ...
                    abs(kaf.alpha(j) + kaf.eta*err) * ...
                    c1*kernel(kaf.dict(j,:),kaf.dict(j,:),kaf.kerneltype,kaf.kernelpar);

                % update influence
                kaf.lambda = kaf.beta*kaf.lambda;
                kaf.lambda(j) = kaf.lambda(j) + 1;

            else % new basis not under quantization threshold
                if m < kaf.M % still room for extra centres
                    kaf.dict = [kaf.dict; x]; % add to codebook
                    kaf.alpha = [kaf.alpha; kaf.eta*err];

                    % update significance for addition of m+1-th centre (15)
                    kaf.E(m+1,1) = 0;
                    kaf.E = kaf.beta*kaf.E + abs(kaf.alpha(m+1)) * ...
                        kernel(kaf.dict,kaf.dict(m+1,:),kaf.kerneltype,kaf.kernelpar);

                    % update influence
                    kaf.lambda = kaf.beta*kaf.lambda;
                    % initial influence
                    kaf.lambda(m+1) = 1;

                else % no room for extra centres
                    [~,L] = min(kaf.E); % centre with lowest significance

                    % update significance for removing L-th centre (17)
                    kaf.E = kaf.E - kaf.lambda(L) * c1*abs(kaf.alpha) .* ...
                        kernel(kaf.dict,kaf.dict(L,:),kaf.kerneltype,kaf.kernelpar);

                    % prune
                    kaf.dict(L,:) = [];
                    kaf.alpha(L) = [];
                    kaf.E(L) = [];
                    kaf.lambda(L) = [];

                    % re-calculate error
                    y_est = kaf.evaluate(x);
                    err = y - y_est;

                    % grow
                    kaf.dict = [kaf.dict; x];
                    kaf.alpha = [kaf.alpha; kaf.eta*err];

                    % update significance for addition of m-th centre (15)
                    kaf.E(m,1) = 0;
                    kaf.E = kaf.beta*kaf.E + abs(kaf.alpha(m)) * ...
                        kernel(kaf.dict,kaf.dict(m,:),kaf.kerneltype,kaf.kernelpar);

                    % update influence
                    kaf.lambda = kaf.beta*kaf.lambda;
                    % initial influence
                    kaf.lambda(m) = 1;
                end
            end
        end
    end
end

repo: kafbox


text:
function fid = writefort15( f15dat, f15out, boudat )
%
if ( strcmp(strtrim(f15out), 'fort.15') )
    disp('Error: an output file name must not be fort.15') ;
    return ;
end

fid = fopen( f15out, 'w' ) ;

fprintf(fid, '%s\n', f15dat.rundes ) ; % RUNDES
fprintf(fid, '%s\n', f15dat.runid  ) ;  % RUNID

% NFOVER
L = length(f15dat.nfover);
dfmt = repmat('%g ',1, L);
fprintf( fid, [dfmt '        \t ! NFOVER \n'], f15dat.nfover ) ;

% NABOUT
fprintf(fid, '%d        \t ! NABOUT \n', f15dat.nabout ) ;

% NSCREEN
fprintf(fid, '%d        \t ! NSCREEN \n', f15dat.nscreen ) ;

% IHOT
fprintf(fid, '%d        \t ! IHOT \n', f15dat.ihot ) ;

% ICS
fprintf(fid, '%d        \t ! ICS \n', f15dat.ics ) ;

% IM
fprintf(fid, '%d        \t ! IM \n', f15dat.im ) ;

% IDEN
if f15dat.im == 20 || f15dat.im == 30
    fprintf( fid, '%d           \t ! IDEN \n', f15dat.iden ) ;
end

% NOLIBF
fprintf(fid, '%d        \t ! NOLIBF \n', f15dat.nolibf ) ;

% NOLIFA
fprintf(fid, '%d        \t ! NOLIFA \n', f15dat.nolifa ) ;

% NOLICA
fprintf(fid, '%d        \t ! NOLICA \n', f15dat.nolica ) ;

% NOLICAT
fprintf(fid, '%d        \t ! NOLICAT \n', f15dat.nolicat ) ;

% NWP
fprintf(fid, '%d        \t ! NWP \n', f15dat.nwp ) ;
if ( f15dat.nwp > 0 )
    for l = 1: f15dat.nwp
        fprintf(fid, '%s\n', f15dat.AttrName(l).name ) ;
    end
end

% NCOR
fprintf(fid, '%d   \t \t ! NCOR \n', f15dat.ncor ) ;

% NTIP
fprintf(fid, '%d   \t \t ! NTIP \n', f15dat.ntip ) ;

% NWS
fprintf(fid, '%d   \t \t ! NWS \n', f15dat.nws ) ;

% NRAMP
fprintf(fid, '%d   \t \t ! NRAMP \n', f15dat.nramp ) ;

% G
fprintf(fid, '%f   \t ! G  \n', f15dat.gravity ) ;

% TAU0
fprintf(fid, '%g   \t \t ! TAU0 \n', f15dat.tau0 ) ;

% Tau0FullDomainMin, Tau0FullDomainMax
if ( abs(f15dat.tau0 + 5.D0) < 1e-10 )
    fprintf(fid, '%f %f   \t ! Tau0FullDomainMin, Tau0FullDomainMax \n', f15dat.tau0minmax ) ;
end

% DTDP
fprintf( fid, '%g  \t \t ! DTDP \n', f15dat.dtdp ) ;

% STATIM
fprintf( fid, '%g  \t \t ! STATIM \n',  f15dat.statim ) ;

% REFTIM
fprintf( fid, '%g  \t \t ! REFTIM \n', f15dat.reftim ) ;

% WTIMINC
if f15dat.nws == 8
  fprintf( fid, '%d %d %d %d %d %g', f15dat.wtimnc ) ;
  fprintf( fid, '  \t ! YYYY MM DD HH24 StormNumber BLAdj \n' ) ;
elseif f15dat.nws >= 19
  fprintf( fid, '%d %d %d %d %d %g %d', f15dat.wtimnc ) ;
  fprintf( fid, '  \t ! YYYY MM DD HH24 StormNumber BLAdj geofactor \n' ) ;
elseif f15dat.nws > 0
  fprintf( fid, '%d  ', f15dat.wtimnc ) ;
  fprintf( fid, '  \t ! WTMINC \n' ) ;
end

% RNDY
fprintf( fid, '%g   \t \t ! RNDY \n', f15dat.rndy ) ;

% DRAMP
L = length(f15dat.dramp);
dfmt = repmat('%g ',1, L);
fprintf( fid, [dfmt ' \t \t ! DRAMP \n'],  f15dat.dramp ) ;

% A00, B00, C00
fprintf( fid, '%g %g %g  \t ! A00, B00, C00 \n', f15dat.a00b00c00 ) ;

% H0
len = length(f15dat.h0) ;
for k = 1: len
  fprintf( fid, '%g ', f15dat.h0(k) ) ;
end
fprintf( fid, '   \t ! H0, 2*dummy, VELMIN \n' ) ;

% SLAM0, SFEA0
fprintf( fid, '%f %f  \t \t ! SLAM0, SFEA0 \n', f15dat.slam ) ;

% CF
if ( f15dat.nolibf <= 2 )
    fprintf( fid, '%f ', f15dat.taucf ) ;
    fprintf( fid, '    \t ! CF \n' ) ;
end

% ESLM, ESLC
if ( f15dat.im <= 2 || f15dat.im == 10 || f15dat.im >= 111111 )
    fprintf( fid, '%f ', f15dat.elsm ) ;
    fprintf( fid, '    \t ! ELSM \n' ) ;
end

% CORI
fprintf( fid, '%f   \t ! CORI \n', f15dat.cori ) ;

% NTIF
fprintf( fid, '%d   \t \t ! NTIF \n', f15dat.ntif ) ;

% Tidal potential
for k = 1: f15dat.ntif
    fprintf( fid, '%s \n',  f15dat.tipotag(k).name ) ;
    fprintf( fid, '%f %16.9e %f %f %f', f15dat.tipotag(k).val ) ;
    fprintf( fid, '   \t !  TPK, AMIGT, ETRF, FFT, FACET \n' ) ;
end

% NBFR
fprintf( fid, '%d   \t \t ! NBFR \n', f15dat.nbfr ) ;
for k = 1: f15dat.nbfr
    fprintf( fid, '%s  \n', f15dat.bountag(k).name ) ;
    fprintf( fid, '%16.9e %f %f \n', f15dat.bountag(k).val ) ;
end

% Open boundary harmonic forcing
for k = 1: f15dat.nbfr
    fprintf(fid, '%s \n', f15dat.opealpha(k).name  ) ;
    fprintf(fid, '%16.9e %16.10g \n', f15dat.opealpha(k).val' ) ;
end

%  ANGINN
fprintf( fid, '%g   \t \t ! ANGINN \n', f15dat.anginn ) ;

% Land boundary
ibtype = [2 12 22 32 52] ;
sm = 0 ;

if ~isempty(boudat)
    for k = 1: length(ibtype)
        sm = sm + sum(~(boudat.ibtype - ibtype(k))) ;
    end
end
if ( sm > 0 )
    fprintf( fid, '%d  \t \t ! NFFR \n', f15dat.nffr ) ;

    nm = 0 ;
    for ib = 1: boudat.nbou
        ibty = boudat.ibtype(ib);

        switch ibty
            case {2,12,22,32,52}
                nm = nm + boudat.nvell(ib) ;
            otherwise
        end
    end

    for k = 1: f15dat.nffr
        fprintf( fid, '%s \n', f15dat.fbountag(k).name ) ;
        fprintf( fid, '%15.8e ', f15dat.fbounspec(k).val ) ;
        fprintf( fid, '\n') ;
    end

    for k = 1: f15dat.nffr
        fprintf( fid, '%s  \n',  f15dat.boualpha(k).name  ) ;

        % val = fscanf(fid, '%f ' ) ; % Must be revisit
        icnt = 0 ;
        for ib = 1: boudat.nbou
            ibty = boudat.ibtype(ib);

            switch ibty
                case {2,12,22,52}
                   for ir = 1: boudat.nvell(ib)
                      icnt = icnt + 1 ;

                      fprintf(fid, '%16.9e  ', f15dat.qnam(k).val(icnt,1:2) ) ;
                      fprintf(fid, ' \n' ) ;
                   end
                case 32
                   for ir = 1: boudat.nvell(ib)
                      icnt = icnt + 1 ;

                      fprintf(fid, '%15.8e ', f15dat.qnam(k).val(icnt,1:5) ) ;
                      fprintf(fid, ' \n' ) ;
                   end
                otherwise
            end
        end
    end
end

% NOUTE, TOUTSE, TOUTFE, NSPOOLE
fprintf( fid, '%d %g %g %d', f15dat.oute ) ;
fprintf( fid, '   \t ! NOUTE, TOUTSE, TOUTFE, NSPOOLE \n' ) ;

% NSTAE
fprintf( fid, '%d  \t \t ! NSTAE \n', f15dat.nstae ) ;

% STAE location
if ( f15dat.nstae > 0 )
    for k = 1: f15dat.nstae
        fprintf(fid, '%f %f %s \n', f15dat.elvstaloc(k,1:2), f15dat.elvstaname{k} ) ;
    end
end

% NOUTV, TOUTV, TOUTFV, NSPOOLV
fprintf( fid, '%d %g %g %d',  f15dat.outv ) ;
fprintf( fid, '  \t ! NOUTV, TOUTV, TOUTFV, NSPOOLV \n') ;

% NSTAV
fprintf( fid, '%d  \t \t ! NSTAV \n', f15dat.nstav ) ;

% STAV location
if ( f15dat.nstav > 0 )
    for k = 1: f15dat.nstav
        fprintf( fid, '%f %f %s \n', f15dat.velstaloc(k,1:2), f15dat.velstaname{k} ) ;
    end
end

% NOUTM, TOUTM, TOUTFM, NSPOOLM
if ( f15dat.nws ~= 0  )
    fprintf(fid, '%d %g %g %d', f15dat.outm ) ;
    fprintf(fid, '  \t ! NOUTM, TOUTM, TOUTFM, NSPOOLM \n') ;

    % NSTAM
    fprintf(fid, '%d   \t \t ! NSTAM \n', f15dat.nstam ) ;

    % STAM location
    if ( f15dat.nstam > 0 )
        for k = 1: f15dat.nstam
            fprintf( fid, '%f %f %s \n', f15dat.metstaloc(k,1:2), f15dat.metstaname{k} ) ;
        end
    end
end

% NOUTGE
fprintf( fid, '%d %g %g %d', f15dat.outge ) ;
fprintf( fid, '  \t ! NOUTGE, ... \n') ;

% NOUTGV
fprintf( fid, '%d %g %g %d', f15dat.outgv ) ;
fprintf( fid, '  \t ! NOUTGV, ... \n') ;

% NOUTGC
if ( f15dat.im == 10 )
    fprintf( fid, '%d %g %g %d', f15dat.outgc ) ;
    fprintf( fid, '  \t ! NOUTGC, ... \n' ) ;
end

% NOUTGM
if ( f15dat.nws ~= 0  )
    fprintf( fid, '%d %g %g %d', f15dat.outgm ) ;
    fprintf( fid, '  \t ! NOUTGM, ... \n' ) ;
end

% NFREQ
fprintf( fid, '%d  \t \t ! NFREQ \n', f15dat.nfreq ) ;
for k = 1: f15dat.nfreq
    fprintf( fid, '%s \n', f15dat.harfreq(k).name ) ;

    fprintf( fid, '%16.9e %f %f \n', f15dat.harfreq(k).val ) ;
end

% THAS, THAF, NHAINC, FMV
fprintf( fid, '%g %g %d %g  \t ! THAS, THAF, NHAINC, FMV \n', f15dat.outhar ) ;

% NHASE, NHASV, NHAGE, NHAGV
fprintf( fid, '%d ', f15dat.outhar_flag ) ;
fprintf( fid, '  \t ! NHASE, NHASV, NHAGE, NHAGV \n' ) ;

% NHSTAR, NHSINC
fprintf( fid, '%d %d  \t ! NHSTAR, NHSINC \n', f15dat.nhstar) ;

% ITITER, ISLDIA, CONVCR, ITMAX
fprintf( fid, '%d %d  %16.9e %d  \t ! ITITER, ISLDIA, CONVCR, ITMAX \n', f15dat.ititer ) ;

% Extra lines including NETCDF & namelist
for k = 1: f15dat.nextraline
    fprintf( fid, '%s\n', f15dat.extraline(k).msg ) ;
end

% If
if find(strcmp(fieldnames(f15dat),'controllist'),1)
for k = 1: length(f15dat.controllist)
    fprintf( fid, '! -- Begin %s Control Namelist -- \n', f15dat.controllist(k).type ) ;
    fprintf( fid, '&%sControl\n', f15dat.controllist(k).type ) ;
    for m = 1:length(f15dat.controllist(k).var)
        val = f15dat.controllist(k).var(m).val;
        if ~ischar(val); val = num2str(val); end
        fprintf( fid, '%s = %s,\n',f15dat.controllist(k).var(m).name,val) ;
    end
    fprintf( fid, '/\n') ;
    fprintf( fid, '! -- End %s Control Namelist -- \n', f15dat.controllist(k).type ) ;
end
end

fclose(fid) ;

end

author: CHLNDDEV


text:
function TooltipString=DoTellMeInfo(var)

Units = [];
Tips = '';
switch var

    %% Imaging Tab
    case 'BandWidth'
        Units ='Hz';
        Tips ='Full receiver bandwidth';
    case 'FOVFreq'
        Units ='m';
        Tips ='Field of view in the frequency encoding direction';
    case 'FOVPhase'
        Units ='m';
        Tips ='Field of view in the first phase encoding direction';
    case 'FlipAng'
        Units ='Degree';
        Tips ='Nominal flip angle of excitation pulse';
    case 'FreqDir'
        Tips ='Frequency encoding direction';
    case 'ResFreq'
        Tips ='Number of voxels in frequency encoding direction';
    case 'ResPhase'
        Tips ='Number of voxels in the first phase encoding direction';
    case 'ScanPlane'
        Tips ='The scanning plane';
    case 'SliceNum'
        Tips ='The number of encoding slice';
    case 'SliceThick'
        Units ='m';
        Tips ='The thickness of one slice';
    case 'TE'
        Units ='s';
        Tips ='The time of echo';
    case 'TEPerTR'
        Tips ='The number of echoes in multiple echo mode, using a number greater than one requires ''MultiEcho'' tab to be loaded';
    case 'TR'
        Units ='s';
        Tips ='The time of repetition';
        %% Advanced Tab
    case 'MasterTxCoil'
        Tips ='The master transmitting coil ID in multi RF transmitting mode';
    case 'MultiTransmit'
        Tips ='The flag for turning on and off multi RF transmitting mode, default mode is ''off'' for single RF transmitting';
    case 'NEX'
        Tips ='The number of excitation';
    case 'NoFreqAlias'
        Tips ='The flag for avoiding aliasing in frequency encoding direction, default ''on'' truncates object outside field of view in frequency encoding direction';
    case 'NoPhaseAlias'
        Tips ='The flag for avoiding aliasing in the first phase encoding direction, default ''on'' truncates object outside field of view in the first phase encoding direction';
    case 'NoSliceAlias'
        Tips ='The flag for avoiding aliasing in the second phase encoding (i.e. slice encoding) direction, default ''on'' truncates object outside field of view in slice encoding direction';
    case 'Shim'
        Tips ='Main static field shimming';
    case 'TEAnchor'
        Tips ='The flag for choosing TE time offset regarding the excitation RF pulse';
        %% Hardware Tab
    case 'B0'
        Units ='T';
        Tips ='Main static magnetic field strength';
    case 'B1Level'
        Units ='T';
        Tips ='A linear scale factor for B1. The input B1+ field with magnitude of this number produces nominal flip angle';
    case 'E1Level'
        Units ='T';
        Tips ='A linear scale factor for E1. When calculating spatial SAR, the input E1+ field is scaled by a factor of nominal RF amplitude divided by this number';
    case 'MaxGrad'
        Units ='T/m';
        Tips ='Maximum allowable gradient strength';
    case 'MaxSlewRate'
        Units ='T/m/s';
        Tips ='Maximum allowable gradient slew rate';
    case 'MinUpdRate'
        Units ='s';
        Tips ='Minimum update time on generating sequence waveform';
    case 'Model'
        Tips ='Model type';
    case 'NoiseLevel'
        Tips ='The level of adjustable noise, the higher the number, the more noise';
    case 'PulseType'
        Tips ='The type of generated sequence pulse';
    case 'SpinPerVoxel'
        Tips ='The number of spins in each voxel. Default one spin per voxel treats T2* equal to T2, use a number greater than one to simulate T2* effect based on T2Star input (linear simulation time cost)';
        %% Recon Tab
    case 'AutoRecon'
        Tips ='The flag for turning on and off automatic image reconstruction after MR signal acquisition';
    case 'ExternalEng'
        Tips ='User defined script for image reconstruction';
    case 'OutputType'
        Tips ='The type of output data including both simulated image and signal';
    case 'ReconEng'
        Tips ='The image reconstruction engine, choosing ''External'' uses external engine which requires ''ExternalEng'' to be provided';
    case 'ReconType'
        Tips ='The type of image reconstruction';
        %% CV Tab
    case 'CV1'
        Tips ='Controllable variable 1';
    case 'CV2'
        Tips ='Controllable variable 2';
    case 'CV3'
        Tips ='Controllable variable 3';
    case 'CV4'
        Tips ='Controllable variable 4';
    case 'CV5'
        Tips ='Controllable variable 5';
    case 'CV6'
        Tips ='Controllable variable 6';
    case 'CV7'
        Tips ='Controllable variable 7';
    case 'CV8'
        Tips ='Controllable variable 8';
    case 'CV9'
        Tips ='Controllable variable 9';
    case 'CV10'
        Tips ='Controllable variable 10';
    case 'CV11'
        Tips ='Controllable variable 11';
    case 'CV12'
        Tips ='Controllable variable 12';
    case 'CV13'
        Tips ='Controllable variable 13';
    case 'CV14'
        Tips ='Controllable variable 14';
        %% SpecialTech Tab
    case 'GRAPPA'
        Tips ='GRAPPA';
    case 'GM'
        Tips ='Generalized Multi-pool exchanging simulation';
    case 'FSE'
        Tips ='Fast Spin Echo';
    case 'EPI'
        Tips ='Echo Planar Imaging';
    case 'DummyPulse'
        Tips ='Dummy pulse';
    case 'CEST'
        Tips ='Chemical Exchange Saturation Transfer simulation';
    case 'PartialEcho'
        Tips ='Partial echo';
    case 'MultiEcho'
        Tips ='Multi echo';
    case 'MT'
        Tips ='Magnetization Transfer simulation';
    case 'ME'
        Tips ='Multiple pool spin Exchange simulation';
    case 'IRPrep'
        Tips ='Inversion recovery preparation';
    case 'T2Prep'
        Tips ='T2 decay preparation';
    case 'Gridding'
        Tips ='Non-Cartesian gridding';
    case 'Spiral'
        Tips ='Spiral imaging';
    case 'SENSE'
        Tips ='SENSE';
    case 'Radial'
        Tips ='Radial imaging';
    case 'RTRecon'
        Tips ='Real Time reconstruction';
    case 'ZeroFilling'
        Tips ='Zero filling k-space';
    case 'VIPR'
        Tips ='VIPR';
    case 'DP_Flag'
        Tips ='The flag for turning on and off dummy pulse';
    case 'DP_FlipAng'
        Units ='Degree';
        Tips ='The flip angle of excitation pulse for dummy pulse';
    case 'DP_Num'
        Tips ='The number of TRs for dummy pulse';
    case 'DP_TR'
        Units ='s';
        Tips ='The time of repetition for dummy pulse';
    case 'EPI_ESP'
        Units ='s';
        Tips ='The echo spacing for EPI';
    case 'EPI_ETL'
        Tips ='The echo train length for EPI';
    case 'EPI_EchoShifting'
        Tips ='The flag for turning on and off echo shifting';
    case 'EPI_ShotNum'
        Tips ='The number of EPI shots, multi shot EPI uses interleave mode';
    case 'FSE_ESP'
        Units ='s';
        Tips ='The echo spacing for FSE';
    case 'FSE_ETL'
        Tips ='The echo train length for FSE';
    case 'FSE_ShotNum'
        Tips ='The number of FSE shots, multi shot FSE uses interleave mode';
    case 'G_Deapodization'
        Tips ='The flag for turning on and off kernel deapodization (i.e. dividing reconstructed image with the iFFT of the gridding kernel)';
    case 'G_KernelSample'
        Tips ='The number of kernel sample point, the more sample points, the better kernel approximation';
    case 'G_KernelWidth'
        Tips ='The full width of kernel in the unit of gridding grid';
    case 'G_OverGrid'
        Tips ='The over gridding factor';
    case 'G_Truncation'
        Tips ='The flag for turning on and off image truncation for reconstructed image';
    case 'TI'
        Units ='s';
        Tips ='The time of inversion recovery';
    case 'MT_Flag'
        Tips ='The flag for turning on and off Magnetization Transfer simulation';
    case 'ME_Flag'
        Tips ='The flag for turning on and off Multiple pool spin Exchange simulation';
    case 'CEST_Flag'
        Tips ='The flag for turning on and off Chemical Exchange Saturation Transfer simulation';
    case 'GM_Flag'
        Tips ='The flag for turning on and off Generalized Multi-pool exchanging simulation';
    case 'RTR_Flag'
        Tips ='The flag for turning on and off real time reconstruction';
    case 'PlotK_Flag'
        Tips ='The flag for turning on and off real time k-space plotting';
    case 'DelayTime'
        Tips ='The delay time for refreshing graphics';
    case 'ME_TEs'
        Units ='s';
        Tips ='An array of multiple echo values';
    case 'R_AngPattern'
        Tips ='The pattern for sampling the angle in k-space';
    case 'R_AngRange'
        Tips ='The range of sampling angle';
    case 'R_SampPerSpoke'
        Tips ='The number of sampling points in each spoke';
    case 'R_SpokeNum'
        Tips ='The number of sampling spokes';
    case 'Sh_X'
        Tips ='The constant for X term';
    case 'Sh_Y'
        Tips ='The constant for Y term';
    case 'Sh_Z'
        Tips ='The constant for Z term';
    case 'Sh_ZX'
        Tips ='The constant for ZX term';
    case 'Sh_ZY'
        Tips ='The constant for ZY term';
    case 'Sh_Z2'
        Tips ='The constant for Z^2 term';
    case 'Sh_XYZ'
        Tips ='The constant for XYZ term';
    case 'Sh_X2_Y2'
        Tips ='The constant for (X^2)(Y^2) term';
    case 'S_ShotNum'
        Tips ='The number of spiral interleaves';
    case 'S_GradientEff'
        Tips ='A linear scale factor for adjusting maximum allowable gradient in spiral design';
    case 'S_F1'
        Tips ='A scale factor for varying FOV with k-space radius r, as FOV(r) = F0 + F1*r + F2*r*r (in variable density spiral design)';
    case 'S_F2'
        Tips ='A scale factor for varying FOV with k-space radius r, as FOV(r) = F0 + F1*r + F2*r*r (in variable density spiral design)';
    case 'tT2Prep'
        Units ='s';
        Tips ='The time of T2 decay preparation';
    case 'ZF_Kz'
        Tips ='The zero filling factor in Kz';
    case 'ZF_Ky'
        Tips ='The number of point in Ky after zero filling';
    case 'ZF_Kx'
        Tips ='The number of point in Kx after zero filling';
    case 'ChemShift'
        Units ='Hz/T';
        Tips ='The chemical shift of the spin';
    case 'Gyro'
        Units ='rad/s/T';
        Tips ='The gyromagnetic ratio of the spin';
    case 'Rho'
        Tips ='The spin density of the spin';
    case 'T1'
        Units ='s';
        Tips ='The longitudinal relaxation time';
    case 'T2'
        Units ='s';
        Tips ='The transverse relaxation time';
    case 'TypeNum'
        Tips ='The number of spin species';
    case 'ZCenter'
        Tips ='The index of the centeral spin in Z direction';
    case 'ZSpin'
        Tips ='The number of the spins in Z direction';
    case 'ZSpinGap'
        Units ='m';
        Tips ='The distance between adjacent spins in Z direction';
    case 'XCenter'
        Tips ='The index of the centeral spin in X direction';
    case 'XSpin'
        Tips ='The number of the spins in X direction';
    case 'XSpinGap'
        Units ='m';
        Tips ='The distance between adjacent spins in X direction';
    case 'Spat_Flag'
        Tips ='The flag to turn on and off 2D spatial RF analysis';
    case 'YCenter'
        Tips ='The index of the centeral spin in Y direction';
    case 'YSpin'
        Tips ='The number of the spins in Y direction';
    case 'YSpinGap'
        Units ='m';
        Tips ='The distance between adjacent spins in Y direction';
    case 'FreqRes'
        Tips ='The number of linear frequency sample points';
    case 'FreqUpLimit'
        Units ='Hz';
        Tips ='The upper limit of frequency range';
    case 'FreqDownLimit'
        Units ='Hz';
        Tips ='The lower limit of frequency range';
    case 'Freq_Flag'
        Tips ='The flag to turn on and off Spatial-Spectral RF analysis';
    case 'ConstantGrad'
        Units ='T/m';
        Tips ='The constant gradient applied when gradient tab is empty';
    case 'dB0'
        Units ='T';
        Tips ='The main static magnetic field offset';
        %% Pulse Waveform
    case 'tS'
        Units ='s';
        Tips ='The starting time point in TR section, any waveform timing in this pulse group is relative to this time point';
    case 'tE'
        Units ='s';
        Tips ='The ending time point in TR section, any waveform timing in this pulse group will be truncated after this time point';
    case 'TRStart'
        Tips ='The starting TR number';
    case 'TREnd'
        Tips ='The ending TR number';
    case 'Freq'
        Tips ='The occurrence frequency (e.g. 1 means occurring every TR section, 5 means occurring every 5 TR sections)';
    case 'Moments'
        Tips ='The flag for turning on and off the zeroth moment display for the gradient';
    case 'LineMarker'
        Tips ='The flag for turning on and off waveform line marker';
    case 'RenderMode'
        Tips ='The k-space rendering mode';
    case 'RenderPoint'
        Tips ='The flag for turning on and off k-space point rendering';
    case 'Tar'
        Tips ='The flag for turning on and off sequence deployment for Toppe';
    case 'PlayToppeMovie'
        Tips ='The flag for turning on and off Toppe movie playback';
    case 'NumTRSkip'
        Tips ='The number of TR sections to skip during Toppe movie playback';
    case 'PlayPulseqMovie'
        Tips ='The flag for turning on and off Pulseq movie playback';
    case 'ShowUnit'
        Tips ='The Pulseq display time unit';
    case 'ShowNum'
        Tips ='The number of time sections to show during Pulseq movie playback';
    case 'Apod'
        Tips ='Apodization methods for RF pulse';
    case 'FA'
        Units ='Degree';
        Tips ='Prescribed flip angle';
    case 'TBP'
        Tips ='The time bandwidth product of RF pulse';
    case 'dt'
        Units ='s';
        Tips ='The time interval of sample points';
    case 'rfPhase'
        Units ='rad';
        Tips ='RF pulse phase';
    case 'rfFreq'
        Units ='Hz';
        Tips ='RF pulse frequency offset';
    case 'tStart'
        Units ='s';
        Tips ='The starting time';
    case 'tEnd'
        Units ='s';
        Tips ='The ending time';
    case 'Switch'
        Tips ='The flag for turning on and off this pulse';
    case 'AnchorTE'
        Tips ='The flag for turning on and off TE reference, TE is calculated from this RF pulse if this flag is turned on';
    case 'Duplicates'
        Tips ='The number of the pulse duplicates, used for creating multiple pulses with the same shape';
    case 'DupSpacing'
        Units ='s';
        Tips ='The time spacing between pulse duplicates';
    case 'CoilID'
        Tips ='The ID of the coil element';
    case 'Notes'
        Tips ='The notes of this object';
    case 'PW'
        Tips ='The measure of the pulse width in Fermi RF pulse';
    case 'SLRPulseType'
        Tips ='The type of this SLR pulse, including ''st''(small tip angle pulse), ''ex''(excitation pulse), ''se''(spin-echo pulse), ''sat''(saturation pulse) and ''inv''(inversion pulse)';
    case 'FilterType'
        Tips ='The type of the applied filter design method, including ''ls''(least squares), ''min''(minimum phase), ''max''(maximum phase), ''pm''(Parks-McClellan equal ripple), and ''ms''(Hamming windowed sinc)';
    case 'PRipple'
        Tips ='The ripple factor at passband';
    case 'SRipple'
        Tips ='The ripple factor at stopband';
    case 'Adiab'
        Tips ='The adiabatic factor';
    case 'MaxB1'
        Units ='T';
        Tips ='The maximum B1 field';
    case 'MaxFreq'
        Units ='Hz';
        Tips ='The maximum RF frequency';
    case 'Lambda'
        Tips ='The lambda adiabatic factor';
    case 'Beta'
        Tips ='The beta adiabatic factor';
    case 'BIRFlag'
        Tips ='The type of BIR pulse, including ''BIR-1'', ''BIR-2'' and ''BIR-4''';
    case 'BIREFFlag'
        Tips ='The type of BIREF pulse, including ''BIREF-1'', ''BIREF-2a'' and ''BIREF-2b''';
    case 'rfGain'
        Tips ='The standard deviation of the normal distribution';
    case 'rfFile'
        Tips ='The path to the file that stores the RF pulse data, quoted using single quotes';
    case 't2Start'
        Units ='s';
        Tips ='The second gradient pulse starting time';
    case 't2End'
        Units ='s';
        Tips ='The second gradient pulse ending time';
    case 'tRamp'
        Units ='s';
        Tips ='The pulse ramp time from zero to plateau, assume symmetric ramp on both side';
    case 'GzAmp'
        Units ='T';
        Tips ='The amplitude of the Gz pulse';
    case 'Gz1Sign'
        Tips ='The polarity of the first gradient pulse, set 0 for nulling';
    case 'Gz2Sign'
        Tips ='The polarity of the second gradient pulse, set 0 for nulling';
    case 'Gz3Sign'
        Tips ='The polarity of the last gradient pulse, set 0 for nulling';
    case 'sRamp'
        Tips ='The sample points on the ramp, use the value of 2 for ignoring the area under the ramp, use values greater than 2 for counting the ramp area';
    case 'Area'
        Units ='1/m';
        Tips ='The area under this gradient pulse';
    case 'nCycles'
        Tips ='The number of cycles of phase across the pixel size';
    case 't1Start'
        Units ='s';
        Tips ='The first gradient pulse starting time';
    case 't1End'
        Units ='s';
        Tips ='The first gradient pulse ending time';
    case 'GzFile'
        Tips ='The path to the file that stores the Gz pulse data, quoted using single quotes';
    case 'GyFile'
        Tips ='The path to the file that stores the Gy pulse data, quoted using single quotes';
    case 'GxFile'
        Tips ='The path to the file that stores the Gx pulse data, quoted using single quotes';
    case 'ADCFile'
        Tips ='The path to the file that stores the ADC pulse data, quoted using single quotes';
    case 'GyAmp'
        Units ='T';
        Tips ='The amplitude of the Gy pulse';
    case 'Gy1Sign'
        Tips ='The polarity of the first gradient pulse, set 0 for nulling';
    case 'Gy2Sign'
        Tips ='The polarity of the second gradient pulse, set 0 for nulling';
    case 'Gy3Sign'
        Tips ='The polarity of the last gradient pulse, set 0 for nulling';
    case 't2Middle'
        Units ='s';
        Tips ='The second encoding gradient pulse middle time';
    case 't3Start'
        Units ='s';
        Tips ='The last gradient pulse starting time';
    case 'tMiddle'
        Units ='s';
        Tips ='The middle time of the pulse';
    case 'tOffset'
        Units ='s';
        Tips ='The time offset of the gradient pulse';
    case 'tGy1'
        Units ='s';
        Tips ='The duration of the first gradient pulse';
    case 'tGy2'
        Units ='s';
        Tips ='The duration of the second gradient pulse';
    case 'GxAmp'
        Units ='T';
        Tips ='The amplitude of the Gx pulse';
    case 'Gx1Sign'
        Tips ='The polarity of the first gradient pulse, set 0 for nulling';
    case 'Gx2Sign'
        Tips ='The polarity of the second gradient pulse, set 0 for nulling';
    case 'Gx3Sign'
        Tips ='The polarity of the last gradient pulse, set 0 for nulling';
    case 'sSample'
        Tips ='The number of linear sample points when ADC flag is 1';
    case 'Ext'
        Tips ='The Ext flag';
    case 'isVardens'
        Tips ='The flag for turning on and off variable density spiral design';
    case 'InOut'
        Tips ='Spiral in or out';
        %% Others
    case 'Name'
        Tips ='The name of the structure';
    case 'Type'
        Tips ='A description about the phantom type';
    case 'XDim'
        Tips ='The number of voxels in X direction';
    case 'YDim'
        Tips ='The number of voxels in Y direction';
    case 'ZDim'
        Tips ='The number of voxels in Z direction';
    case 'XDimRes'
        Units ='m';
        Tips ='The spatial resolution in X direction';
    case 'YDimRes'
        Units ='m';
        Tips ='The spatial resolution in Y direction';
    case 'ZDimRes'
        Units ='m';
        Tips ='The spatial resolution in Z direction';
    case 'Grid'
        Tips ='Turn on and off grid';
    case 'Box'
        Tips ='Turn on and off boundary box';
    case 'CameraTool'
        Tips ='Hide or show Matlab camera tool';
    case 'Color'
        Tips ='The display color';
    case 'Alpha'
        Tips ='The display transparency';
    case 'Radius'
        Units ='m';
        Tips ='The radius of the object';
    case 'CenterX'
        Units ='m';
        Tips ='The X coordinate of the object center';
    case 'CenterY'
        Units ='m';
        Tips ='The Y coordinate of the object center';
    case 'CenterZ'
        Units ='m';
        Tips ='The Z coordinate of the object center';
    case 'FaceNum'
        Tips ='The number of the faces for the object';
    case 'TypeIdx'
        Tips ='An index number of the spin species, used when the phantom has multiple spin species. The index must not exceed the ''TypeNum''';
    case 'TypeFlag'
        Tips ='A flag number for describing the type of the spin, 0 for free pool and 1 for bound pool';
    case 'LineShapeFlag'
        Tips ='A flag number for describing RF saturation line shape for bound proton pool, 0 for super-Lorentzian and 1 for Gaussian (:ToDo), the flag is ignored for free pool';
    case 'ECon'
        Units ='S/m';
        Tips ='A array with size of [1 3] for tissue electrical conductivity (optional)';
    case 'MassDen'
        Units ='kg/m^3';
        Tips ='The tissue mass density (optional)';
    case 'T2Star'
        Units ='s';
        Tips ='The T2* relaxation time';
    case 'K'
        Units ='1/s';
        Tips ='A array with the size of [1 TypeNum] for describing the exchange rate of the spin, ignored for regular phantom';
    case 'RadiusX'
        Units ='m';
        Tips ='The X semi-axis length of the ellipsoid';
    case 'RadiusY'
        Units ='m';
        Tips ='The Y semi-axis length of the ellipsoid';
    case 'RadiusZ'
        Units ='m';
        Tips ='The Z semi-axis length of the ellipsoid';
    case 'Length'
        Units ='m';
        Tips ='The length of the object';
    case 'Height'
        Units ='m';
        Tips ='The height of the pyramid';
    case 'Colormap'
        Tips ='The colormap for the field';
    case 'CLimDown'
        Tips ='The lower bound of color limits';
    case 'CLimUp'
        Tips ='The upper bound of color limits';
    case 'CoilDisplay'
        Tips ='The flag for turning on and off coil display';
    case 'CoilShow'
        Tips ='The flag for choosing active coil for field display';
    case 'Mode'
        Tips ='The B1 field display mode';
    case 'FieldType'
        Tips ='The flag to choose B1 field or E1 field, note E1 field only support ''Magnitude'' display mode';
    case 'Plane'
        Tips ='The flag for activating field slicing plane';
    case 'Azimuth'
        Units ='rad';
        Tips ='The azimuth angle of the plane';
    case 'Elevation'
        Units ='rad';
        Tips ='The elevation angle of the plane';
    case 'PosZ'
        Units ='m';
        Tips ='The Z position of object center';
    case 'PosY'
        Units ='m';
        Tips ='The Y position of object center';
    case 'PosX'
        Units ='m';
        Tips ='The X position of object center';
    case 'CurrentDir'
        Tips ='The current direction in the coil circle, 1 for clockwise, -1 for counterclockwise';
    case 'Scale'
        Tips ='The scale factor for the field amplitude';
    case 'Segment'
        Tips ='The number of line segments for approximating circle, MRiLab requires the same ''Segment'' for each coil circle';
    case 'Width'
        Units ='m';
        Tips ='The width of the object';
    case 'B1File'
        Tips ='The path to the file that stores the B1 field data, quoted using single quotes';
    case 'E1File'
        Tips ='The path to the file that stores the E1 field data, quoted using single quotes';
    case 'Interp'
        Tips ='The interpolation method';
    case 'GradZ'
        Tips ='The linear gradient in Z direction';
    case 'GradY'
        Tips ='The linear gradient in Y direction';
    case 'GradX'
        Tips ='The linear gradient in X direction';
    case 'DeltaZ'
        Units ='m';
        Tips ='The width of Gaussian function in Z direction';
    case 'DeltaY'
        Units ='m';
        Tips ='The width of Gaussian function in Y direction';
    case 'DeltaX'
        Units ='m';
        Tips ='The width of Gaussian function in X direction';
    case 'Equation'
        Tips ='A field described with a symbolic equation';
    case 'MagFile'
        Tips ='The path to the file that stores the dB0 field data, quoted using single quotes';
    case 'GradLine'
        Tips ='The gradient sequence line';
    case 'DispMode'
        Tips ='The display mode';
    case 'GradZEqu'
        Tips ='A symbolic equation for gradient field vector in Z direction';
    case 'GradYEqu'
        Tips ='A symbolic equation for gradient field vector in Y direction';
    case 'GradXEqu'
        Tips ='A symbolic equation for gradient field vector in X direction';
    case 'GradFile'
        Tips ='The path to the file that stores the gradient field data, quoted using single quotes';
    case 'Object'
        Tips ='The object model, currently only supports ''Sphere''';
    case 'ViewPoint'
        Tips ='A default view point';
    case 'ZoomOut'
        Tips ='A factor of view zoom out';
    case 'Sample'
        Tips ='The sample steps between two adjacent positions during movement';
    case 'Repeat'
        Tips ='The repeat time of playback';
    case 'Direction'
        Tips ='A vector describing translation direction in 3D space';
    case 'Displacement'
        Units ='m';
        Tips ='An equation of translation displacement pattern with respect to time';
    case 'Axis'
        Tips ='A vector describing rotation axis in 3D space';
    case 'Angle'
        Units ='rad';
        Tips ='An equation of rotation angle with respect to time';
    case 'LocZ'
        Tips ='The Z location of the selected voxel';
    case 'LocY'
        Tips ='The Y location of the selected voxel';
    case 'LocX'
        Tips ='The X location of the selected voxel';
    case 'WindowSize'
        Tips ='The window width of the spin evolution plot';
    case 'ISOHighlight'
        Tips ='The flag for turning on and off isocenter mark';
    case 'Axes'
        Tips ='The flag for turning on and off axes label';
    case 'N_Gram'
        Units ='g';
        Tips ='The number to specify averaged N-gram SAR, set to 0 indicating unaveraged spatial SAR';
    case 'N_Second'
        Units ='s';
        Tips ='The nominal time window for SAR calculation';
end

if ~isempty(Units)
    TooltipString = [var '(' Units '): ' Tips];
else
    TooltipString = [var ': ' Tips];
end

end

author: leoliuf

repo: MRiLab

save_path: github-repos/MATLAB/leoliuf-MRiLab


text:
function desktopPath = generateDesktopPath(dirnameSrc)
if ~exist('dirnameSrc','var')
    dirnameSrc = filesepStandard(pwd);
end

if ispc()
    if ~exist([dirnameSrc, 'desktopPath.txt'],'file')
        system(sprintf('call "%sgenerateDesktopPath.bat"', dirnameSrc));
    end

    if exist([dirnameSrc, 'desktopPath.txt'],'file')
        fid = fopen([dirnameSrc, 'desktopPath.txt'],'rt');
        line = fgetl(fid);
        line(line=='"')='';
        desktopPath = strtrim(line);
        fclose(fid);
    else
        desktopPath = sprintf('%%userprofile%%');
    end
else
    desktopPath = fullpath('~/Desktop');
end
try
    fprintf('Desktop path:  %s\n', desktopPath);
    desktopPath = filesepStandard(desktopPath);
    desktopPath = desktopPath(1:end-1);
catch
    fprintf('ERROR: Desktop path not found\n');
    desktopPath = '';
end

author: BUNPC


text:
classdef Resampler < Algorithm

    properties (Access = public)
        samplingInterval = 1;
    end

    methods (Access = public)

        function obj = Resampler(samplingInterval)
            if nargin > 0
                obj.samplingInterval = samplingInterval;
            end
            obj.name = 'resampler';
            obj.inputPort = DataType.kSignal;
            obj.outputPort = DataType.kSignal;
        end

        function dataFiltered = compute(obj,data)
            dataFiltered = data(1:obj.samplingInterval:end,:);
        end

        function str = toString(obj)
            str = sprintf('%s_%d',obj.name,obj.samplingInterval);
        end

        function editableProperties = getEditableProperties(obj)
            editableProperties = Property('samplingInterval',obj.samplingInterval,1,100);
        end
    end
end

author: avenix

repo: WDK

save_path: github-repos/MATLAB/avenix-WDK


text:
function CPD = maximize_params(CPD, temp)
% MAXIMIZE_PARAMS Set the params of a tabular node to their ML/MAP values.
% CPD = maximize_params(CPD, temp)

if ~adjustable_CPD(CPD), return; end

%assert(approxeq(sum(CPD.counts(:)), CPD.nsamples)); % false!
switch CPD.prior_type
 case 'none',
  counts = reshape(CPD.counts, size(CPD.CPT));
  CPD.CPT = mk_stochastic(counts);
 case 'dirichlet',
  counts = reshape(CPD.counts, size(CPD.CPT));
  CPD.CPT = mk_stochastic(counts + CPD.dirichlet);

 % case 'entropic',
%   % For an HMM,
%   % CPT(i,j) = pr(X(t)=j | X(t-1)=i) = transprob(i,j)
%   % counts(i,j) = E #(X(t-1)=i, X(t)=j) = exp_num_trans(i,j)
%   Z = 1-temp;
%   fam_sz = CPD.sizes;
%   psz = prod(fam_sz(1:end-1));
%   ssz = fam_sz(end);
%   counts = reshape(CPD.counts, psz, ssz);
%   CPT = zeros(psz, ssz);
%   for i=CPD.entropic_pcases(:)'
%     [CPT(i,:), logpost] = entropic_map_estimate(counts(i,:), Z);
%   end
%   non_entropic_pcases = mysetdiff(1:psz, CPD.entropic_pcases);
%   for i=non_entropic_pcases(:)'
%     CPT(i,:) = mk_stochastic(counts(i,:));
%   end
%   %for i=1:psz
%   %  [CPT(i,:), logpost] = entropic_map(counts(i,:), Z);
%   %end
%   if CPD.trim & (temp < 2) % at high temps, we would trim everything!
%     % grad(j) = d log lik / d theta(i ->j)
%     % CPT(i,j) = 0 => counts(i,j) = 0
%     % so we can safely replace 0s by 1s in the denominator
%     denom = CPT(i,:) + (CPT(i,:)==0);
%     grad = counts(i,:) ./ denom;
%     trim = find(CPT(i,:) <= exp(-(1/Z)*grad)); % eqn 32
%     if ~isempty(trim)
%       CPT(i,trim) = 0;
%       if all(CPD.trimmed_trans(i,trim)==0) % trimming for 1st time
% 	disp(['trimming CPT(' num2str(i) ',' num2str(trim) ')'])
%       end
%       CPD.trimmed_trans(i,trim) = 1;
%     end
%   end
%   CPD.CPT = myreshape(CPT, CPD.sizes);
end

author: bayesnet

repo: bnt

save_path: github-repos/MATLAB/bayesnet-bnt


text:
%{
 * Copyright (C) 2020-2030, The Regents of The University of Michigan.
 * All rights reserved.
 * This software was developed in the Biped Lab (https://www.biped.solutions/)
 * under the direction of Jessy Grizzle, grizzle@umich.edu. This software may
 * be available under alternative licensing terms; contact the address above.
 * 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 OWNER 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.
 * The views and conclusions contained in the software and documentation are those
 * of the authors and should not be interpreted as representing official policies,
 * either expressed or implied, of the Regents of The University of Michigan.
 *
 * AUTHOR: Bruce JK Huang (bjhuang[at]umich.edu)
 * WEBSITE: https://www.brucerobot.com/
%}

function opt = optimizeCost(opt, X, target_size, box_width)

    R_init = rotx(opt.rpy_init(1)) * roty(opt.rpy_init(2)) * rotz(opt.rpy_init(3));
    opt.H_init(1:3, 1:3) = R_init;

    if opt.UseCentroid
        centroid = mean(X, 2);
        opt.H_init(1:3, 4) = -centroid(1:3);
    else
        opt.H_init(1:3, 4) = opt.T_init;
    end

    switch opt.method
        case 'Customize'
            opt.metric = "PointToAxis";
            opt.unit = "L1-inspired";
            tic;
            optimizeCustomizedCost(opt);
            opt.computation_time = toc;
        case 'Constraint Customize'
            opt.metric = "PointToAxis";
            opt.unit = "L1-inspired";
            tic;
            opt = optimizeConstraintCustomizedCost(opt, X, target_size, box_width);
            opt.computation_time = toc;
        case 'The Wall'
            opt.metric = "PointToAxis";
            opt.unit = "L1-inspired";
            tic;
            opt = optimizeBoundaries(opt, opts, X);
            opt.computation_time = toc;
        case 'Coherent Point Drift'
            opt.metric = "--";
            opt.unit = "RMSE";
            tic;
            optimizeCPD(opt);
            opt.computation_time = toc;
        case 'Iterative Closest Point (point)'
            opt.metric = "PointToPoint";
            opt.unit = "RMSE";
            tic;
            optimizeICPPoint(opt);
            opt.computation_time = toc;
        case 'Iterative Closest Point (plane)'
            opt.metric = "PointToPlane";
            opt.unit = "RMSE";
            tic;
            optimizeICPPlane(opt);
            opt.computation_time = toc;
        case 'Normal-distributions Transform'
            opt.metric = "--";
            opt.unit = "RMSE";
            tic;
            optimizeNDT(opt);
            opt.computation_time = toc;
        case 'GICP-SE3'
            tic;
            optimizeGICP_SE3(opt)
            opt.computation_time = toc;
        case 'GICP-SE3 (plane)'
            tic;
            optimizeGICPPlane_SE3(opt)
            opt.computation_time = toc;
        case 'GICP-SE3-costimized'
            tic;
            optimizeCustomizedGICP_SE3(opt);
            opt.computation_time = toc;
        case 'Two Hollow Strips'
            tic;
            optimizeHollowStrips(opt);
            opt.computation_time = toc;
        case '3D IoU'
        case 'Project'
            tic;
            OptimizeUsingRobustNormalVector(opt);
            opt.computation_time = toc;
        otherwise
            disp('no such optimization method')
            return
    end
end

author: UMich-BipedLab


text:
function wts = cliqf ( nt, t, kind, alpha, beta, a, b, lo )

%*****************************************************************************80
%
%% CLIQF computes a classical quadrature formula, with optional printing.
%
%  Discussion:
%
%    This routine computes all the weights of an interpolatory
%    quadrature formula with
%    1. only simple knots and
%    2. a classical weight function with any valid A and B, and
%    3. optionally prints the knots and weights and a check of the moments.
%
%    To evaluate this quadrature formula for a given function F,
%    call routine EIQFS.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    05 January 2010
%
%  Author:
%
%    Original FORTRAN77 version by Sylvan Elhay, Jaroslav Kautsky.
%    MATLAB version by John Burkardt.
%
%  Reference:
%
%    Sylvan Elhay, Jaroslav Kautsky,
%    Algorithm 655: IQPACK, FORTRAN Subroutines for the Weights of
%    Interpolatory Quadrature,
%    ACM Transactions on Mathematical Software,
%    Volume 13, Number 4, December 1987, pages 399-415.
%
%  Parameters:
%
%    Input, integer NT, the number of knots.
%
%    Input, real T(NT), the knots.
%
%    Input, integer KIND, the rule.
%    1, Legendre,             (a,b)       1.0
%    2, Chebyshev Type 1,     (a,b)       ((b-x)*(x-a))^(-0.5)
%    3, Gegenbauer,           (a,b)       ((b-x)*(x-a))^alpha
%    4, Jacobi,               (a,b)       (b-x)^alpha*(x-a)^beta
%    5, Generalized Laguerre, (a,+oo)     (x-a)^alpha*exp(-b*(x-a))
%    6, Generalized Hermite,  (-oo,+oo)   |x-a|^alpha*exp(-b*(x-a)^2)
%    7, Exponential,          (a,b)       |x-(a+b)/2.0|^alpha
%    8, Rational,             (a,+oo)     (x-a)^alpha*(x+b)^beta
%    9, Chebyshev Type 2,     (a,b)       ((b-x)*(x-a))^(+0.5)
%
%    Input, real ALPHA, the value of Alpha, if needed.
%
%    Input, real BETA, the value of Beta, if needed.
%
%    Input, real A, B, the interval endpoints.
%
%    Input, integer LO, indicates what is to be done.
%    > 0, compute and print weights and moments check.
%    = 0, compute weights.
%    < 0, compute and print weights.
%
%    Output, real WTS(NT), the weights.
%
  key = 1;
  mlt = zeros(nt,1);
  mlt(1:nt) = 1;
  ndx = zeros(nt,1);

  wts = ciqf ( nt, t, mlt, nt, ndx, key, kind, alpha, beta, a, b, lo );

  return
end

author: johannesgerer

repo: jburkardt-m


text:
function test_approx_test02 ( )

%*****************************************************************************80
%
%% TEST_APPROX_TEST02 shows how P00_STORY can be called.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    12 August 2011
%
%  Author:
%
%    John Burkardt
%
  fprintf ( 1, '\n' );
  fprintf ( 1, 'TEST_APPROX_TEST02\n' );
  fprintf ( 1, '  P00_STORY prints the problem "story".\n' );

  prob_num = p00_prob_num ( );

  for prob = 1 : prob_num

    fprintf ( 1, '\n' );
    fprintf ( 1, '  Problem %d\n', prob );

    p00_story ( prob );

  end

  return
end


text:
function varargout=beta_w(varargin)

switch class(varargin{1})

    case 'double'
        varargout{1} = beta(varargin{2},varargin{1});

    case 'char'
        operator = CreateBasicOperator('convex','decreasing','positive','callback');
        operator.domain = [0 inf];
        varargout{1} = [];
        varargout{2} = operator;
        varargout{3} = varargin{3};

    otherwise
        error(['SDPVAR/' upper(mfilename) ' called with weird argument']);
end

author: yalmip

repo: YALMIP

save_path: github-repos/MATLAB/yalmip-YALMIP


text:
% Given a reference multidimensional signal y and a series of penalty terms P(x,lambda,d,p), solves the generalized Total Variation
% proximity operator
%
%        min_x 0.5 ||x-y||^2 + sum_i P(x,lambda_i,d_i,p_i) .
%
% where P(x,lambda_i,d_i,p_i) = lambda_i * sum_j TV(x(d_i)_j,p_i) with x(d)_j every possible 1-dimensional slice of x
% following the dimension d_i, TV(x,p) the TV-Lp prox operator for x.
%
% Inputs:
%   - y: input signal.
%   - lambdas: vector of lambda penalties of each penalty term.
%   - ds: vector of dimensions of application of each penalty term.
%   - norms: vector of norms of each penalty term (1, 2 or inf are accepted).
%   - [threads]: number of threads to use (default: 1)
%
% Outputs:
%   - x: solution of the proximity problem.
%   - info: statistical info of the run algorithm:
%       info.iters: number of major iterations run.
%       info.stop: value of the stopping criterion.
%
% Examples:
%
%   - Filter 2D signal using TV-L1 norm:
%        prox_TVgen(X,[lambda lambda],[1 2],[1 1])
%
%   - Filter 2D signal using TV-L2 norm:
%        prox_TVgen(X,[lambda lambda],[1 2],[2 2])
%
%   - Filter 2D signal using TV-L1 norm for the rows, TV-L2 for the columns, and different penalties:
%        prox_TVgen(X,[lambdaRows lambdaCols],[1 2],[1 2])
%
%   - Filter 1D signal using both TV-L1 and TV-L2 norms:
%        prox_TVgen(X,[lambda1 lambda2],[1 1],[1 2])
%
%   - Filter 3D signal using TV-L1 norm:
%        prox_TVgen(X,[lambda lambda lambda],[1 2 3],[1 1 1])
%
%   - Filter 3D signal using TV-L2 norm, not penalizing over the second dimension:
%        prox_TVgen(X,[lambda lambda],[1 3],[2 2])
%
%   - Filter 2D signal using both TV-L1 and TV-L2 norms:
%        prox_TVgen(X,[lambda1 lambda1 lambda2 lambda2],[1 2 1 2],[1 1 2 2])
%
%   - ... and so on, any combination of norms and dimensions is possible.
function [x,info] = prox_TVgen(y,lambdas,ds,norms,threads)
    % Check inputs
    if nargin < 5, threads = 1; end;
    if length(lambdas) ~= length(ds) || length(ds) ~= length(norms)
        fprintf(1,'ERROR (prox_TVgen): arguments defining penalties differ in length.\n');
        x = y;
        return;
    end;
    if sum(norms == 1 | norms == 2) ~= length(norms)
        fprintf(1,'ERROR (prox_TVgen): unacceptable norms requested. Available norms: 1, 2.\n');
        x = y;
        return;
    end;

    % Invoke C solver
    [x,in] = solveTVgen_PDykstrac(y,lambdas,ds,norms,threads);
    info.iters = in(1);
    info.stop = in(2);
end

author: beckel

repo: nilm-eval

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


text:
function skullStruct = readLeksellSkullFile(filename)
%"readLeksellSkullFile"
%   Uses the docodeLeksellData function to decode the data from the Leksell
%   files and places the data into meaningful structures. This
%   information is not currently used in CERR.
%
%KRK 05/30/07
%
%Usage:
%   function skullStruct = readLeksellSkullFile(filename)
%
% Copyright 2010, Joseph O. Deasy, on behalf of the CERR development team.
%
% This file is part of The Computational Environment for Radiotherapy Research (CERR).
%
% CERR development has been led by:  Aditya Apte, Divya Khullar, James Alaly, and Joseph O. Deasy.
%
% CERR has been financially supported by the US National Institutes of Health under multiple grants.
%
% CERR is distributed under the terms of the Lesser GNU Public License.
%
%     This version of CERR 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.
%
% CERR 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 CERR.  If not, see <http://www.gnu.org/licenses/>.

fid = fopen(filename, 'r', 'b');

data = decodeLeksellData(fid);

fclose(fid);

% If data only has one cell, it's empty since there is always one cell that
%   contains no meaningful data
if length(data) < 2
    skullStruct = [];
    return;
end

% Get meaningful data
skullData = data{1};

% Get the row size of the skull data
skullStruct.colLen = skullData{1};
% Get the column size of the skull data
skullStruct.rowLen = skullData{2};
% Get the spherical coordinate angles (theta = colAngles, phi =
%   rowAngles).  Ex: rowAngles(2) corresponds to the entire 2nd row
%   in the skullStruct.radii matrix
skullStruct.rowAngles = skullData{3};
skullStruct.colAngles = skullData{4};
% Get the sphereical coordinate radii (rho)
skullStruct.radii = reshape(skullData{5}, [skullData{2}, skullData{1}])';
% skullData{6} is unknown, it is the same size as skullData{5}, but is
%   filled with small floating point numbers.
skullStruct.mystery_value_1 = reshape(skullData{6}, [skullData{2}, skullData{1}])';

author: cerr

repo: CERR

save_path: github-repos/MATLAB/cerr-CERR