text:
%kScrollbar 'Scrollbar Object Menuform'
% This MatLab function was automatically generated by a converter (KhorosToMatLab) from the Khoros Scrollbar.pane file
%
% Parameters:
%
% Example: kScrollbar( {})
%
% Khoros helpfile follows below:
%
% .begin tagged
%
% .item "The Scrollbar GUI Object"
% The Scrollbar object allows users to reposition data that is too large to
% fit in a viewing window.  The Scrollbar consists of a long, thin area called the
% scrollbar trough, within which is a smaller, movable button
% called the slider.  Scrollbars may be horizontal or vertical.
% Horizontal scrollbars have the minimum value at the top and
% the maximum value at the bottom; vertical scrollbars have the
% minimum value at the left and the maximum value at the right.
%
% With a horizontal scrollbar, clicking in the scrollbar trough
% to the right of the slider will cause the scrollbar to
% increment in value; clicking in the scrollbar trough to the
% left of the slider will cause the scrollbar to decrement in
% value.  Similarly with a vertical scrollbar, clicking above
% the slider will decrement the value while clicking below the
% slider will increment the value.  The slider itself can be
% dragged back and forth (with a horizontal scrollbar) or up and
% down (with a vertical scrollbar) to change the value.
%
% .item "Minimum Value"
% The minimum world coordinate value that is to be represented by the scrollbar.
%
% .item "Maximum Value"
% The maximum world coordinate value that is to be represented by the scrollbar.
%
% .item "Default Value"
% The default world coordinate value that is to be represented by the scrollbar.
%
% .item "Increment Value"
% The value by which the scrollbar slider will be incremented
% or decremented when the user clicks in the scrollbar trough.
%
% .item "Orientation"
% The scrollbar may be vertical or horizontal.
% .end tagged

function varargout = kScrollbar(varargin)
Inputs={};
if nargin ==0
  arglist={'',''};
elseif nargin ==1
  arglist=varargin{1};
else error('Usage: [out1,..] = kScrollbar(arglist).');
end
if size(arglist,2)~=2
  error('arglist must be of form {''ParameterTag1'',value1;''ParameterTag2'',value2}')
 end
narglist={};
maxval={};
minval={};
istoggle=[];
was_set=istoggle * 0;
paramtype={};
% identify the input arrays and assign them to the arguments as stated by the user
if ~iscell(Inputs)
Inputs = {Inputs};
end
NumReqOutputs=0; nextinput=1; nextoutput=1;
  for ii=1:size(arglist,1)
  wasmatched=0;
  for jj=1:size(narglist,1)
   if strcmp(arglist{ii,1},narglist{jj,1})  % a given argument was matched to the possible arguments
     wasmatched = 1;
     was_set(jj) = 1;
     if strcmp(narglist{jj,2}, '__input')
      if (nextinput > length(Inputs))
        error(['Input ' narglist{jj,1} ' has no corresponding input!']);
      end
      narglist{jj,2} = 'OK_in';
      nextinput = nextinput + 1;
     elseif strcmp(narglist{jj,2}, '__output')
      if (nextoutput > nargout)
        error(['Output nr. ' narglist{jj,1} ' is not present in the assignment list of outputs !']);
      end
      if (isempty(arglist{ii,2}))
        narglist{jj,2} = 'OK_out';
      else
        narglist{jj,2} = arglist{ii,2};
      end

      nextoutput = nextoutput + 1;
      if (minval{jj} == 0)
         NumReqOutputs = NumReqOutputs - 1;
      end
     elseif isstr(arglist{ii,2})
      narglist{jj,2} = arglist{ii,2};
     else
        if strcmp(paramtype{jj}, 'Integer') & (round(arglist{ii,2}) ~= arglist{ii,2})
            error(['Argument ' arglist{ii,1} ' is of integer type but non-integer number ' arglist{ii,2} ' was supplied']);
        end
        if (minval{jj} ~= 0 | maxval{jj} ~= 0)
          if (minval{jj} == 1 & maxval{jj} == 1 & arglist{ii,2} < 0)
            error(['Argument ' arglist{ii,1} ' must be bigger or equal to zero!']);
          elseif (minval{jj} == -1 & maxval{jj} == -1 & arglist{ii,2} > 0)
            error(['Argument ' arglist{ii,1} ' must be smaller or equal to zero!']);
          elseif (minval{jj} == 2 & maxval{jj} == 2 & arglist{ii,2} <= 0)
            error(['Argument ' arglist{ii,1} ' must be bigger than zero!']);
          elseif (minval{jj} == -2 & maxval{jj} == -2 & arglist{ii,2} >= 0)
            error(['Argument ' arglist{ii,1} ' must be smaller than zero!']);
          elseif (minval{jj} ~= maxval{jj} & arglist{ii,2} < minval{jj})
            error(['Argument ' arglist{ii,1} ' must be bigger than ' num2str(minval{jj})]);
          elseif (minval{jj} ~= maxval{jj} & arglist{ii,2} > maxval{jj})
            error(['Argument ' arglist{ii,1} ' must be smaller than ' num2str(maxval{jj})]);
          end
        end
     end
     if ~strcmp(narglist{jj,2},'OK_out') &  ~strcmp(narglist{jj,2},'OK_in')
       narglist{jj,2} = arglist{ii,2};
     end
   end
   end
   if (wasmatched == 0 & ~strcmp(arglist{ii,1},''))
        error(['Argument ' arglist{ii,1} ' is not a valid argument for this function']);
   end
end
% match the remaining inputs/outputs to the unused arguments and test for missing required inputs
 for jj=1:size(narglist,1)
     if  strcmp(paramtype{jj}, 'Toggle')
        if (narglist{jj,2} ==0)
          narglist{jj,1} = '';
        end;
        narglist{jj,2} = '';
     end;
     if  ~strcmp(narglist{jj,2},'__input') && ~strcmp(narglist{jj,2},'__output') && istoggle(jj) && ~ was_set(jj)
          narglist{jj,1} = '';
          narglist{jj,2} = '';
     end;
     if strcmp(narglist{jj,2}, '__input')
      if (minval{jj} == 0)  % meaning this input is required
        if (nextinput > size(Inputs))
           error(['Required input ' narglist{jj,1} ' has no corresponding input in the list!']);
        else
          narglist{jj,2} = 'OK_in';
          nextinput = nextinput + 1;
        end
      else  % this is an optional input
        if (nextinput <= length(Inputs))
          narglist{jj,2} = 'OK_in';
          nextinput = nextinput + 1;
        else
          narglist{jj,1} = '';
          narglist{jj,2} = '';
        end;
      end;
     else
     if strcmp(narglist{jj,2}, '__output')
      if (minval{jj} == 0) % this is a required output
        if (nextoutput > nargout & nargout > 1)
           error(['Required output ' narglist{jj,1} ' is not stated in the assignment list!']);
        else
          narglist{jj,2} = 'OK_out';
          nextoutput = nextoutput + 1;
          NumReqOutputs = NumReqOutputs-1;
        end
      else % this is an optional output
        if (nargout - nextoutput >= NumReqOutputs)
          narglist{jj,2} = 'OK_out';
          nextoutput = nextoutput + 1;
        else
          narglist{jj,1} = '';
          narglist{jj,2} = '';
        end;
      end
     end
  end
end
if nargout
   varargout = cell(1,nargout);
else
  varargout = cell(0);
end
global KhorosRoot
if exist('KhorosRoot') && ~isempty(KhorosRoot)
w=['"' KhorosRoot];
else
if ispc
  w='"C:\Program Files\dip\khorosBin\';
else
[s,w] = system('which cantata');
w=['"' w(1:end-8)];
end
end
callKhoros([w 'Scrollbar" '],Inputs,narglist);


text:
%% lighten underlay volume - for white background
function Z = lighten_underlay_edges(Z, maxpercent)

Zi = Z(:);

mindat = min(Zi);
maxdat = max(Zi);

Z(Z == mindat) = maxdat;

% higher ending values of k will be more 'softening'; lower =
% more dark edges

myprctiles = prctile(Zi, [0:1:maxpercent]); % speed up!!

for k = 1:maxpercent % for each of the darkest/lowest-value percentiles

    w = 1 - k/100; % weights, for wtd average of orig and inverted, 1 is all inverted
    % (bilinear weighting)

    %isdark = Zi < prctile(Zi, k) & Zi >= prctile(Zi, k - 1);
    isdark = Zi < myprctiles(k + 1) & Zi >= myprctiles(k);

    Z(isdark) = maxdat;

    %Z(isdark) = w .* (maxdat) + (1 - w) .* Z(isdark);

end % the weight value loop

end

author: canlab

repo: CanlabCore

save_path: github-repos/MATLAB/canlab-CanlabCore


text:
% [succeeded,aviname] = make_ctrax_result_movie('param',value,...)
% Options:
% 'moviename': name of raw movie to annotate
% 'trxname': name of mat file containing trajectories, or struct of
% trajectories
% 'aviname': name of movie to output to
% 'colors': colors to plot each fly
% 'nzoomr': number of rows of zoom fly boxes
% 'nzoomc': number of columns of zoom fly boxes
% 'boxradius': radius of the zoom fly box in pixels
% 'taillength': number of frames of trajectory to plot behind each fly
% 'zoomflies': flies to zoom in on
% 'fps': frames per second of output movie
% 'maxnframes': number of frames to output
% 'firstframe': first frame to output
% 'compression': compressor to use when outputting (pc only, will be set to
% 'none' for linux).
% 'figpos': position of figure
% if any parameters are not given, the user will be prompted for these
function [succeeded,aviname,figpos,height,width] = make_jaaba_result_movie_plotallflies(expdir,varargin)

succeeded = false;
defaults.boxradius = 1.5;
defaults.taillength = 30;
defaults.fps = 20;
defaults.zoomflies = [];
defaults.nzoomr = 5;
defaults.nzoomc = 3;
defaults.compression = 'None';
allowedcompressions = {'Indeo3', 'Indeo5', 'Cinepak', 'MSVC', 'RLE', 'None','Uncompressed AVI','Motion JPEG AVI'};
useVideoWriter = exist('VideoWriter','file');
mencoderoptions = '';
mencoder_maxnframes = inf;

% hand-chosen colors
% behavior2color = {'WingExtension',[.4392,.1882,.6275]
%   'Stop',[0,0,.5156]
%   'Copulation',[0,0,.875]
%   'AttemptCopulation',[0,.2344,1]
%   'Backup',[0,.9688,1]
%   'Crabwalk',[.3281,1,.6719]
%   'PivotHead',[.6875,1,.3125]
%   'Walk',[0,.75,0]
%   'PivotTail',[.9,.9,0]
%   'Jump',[1,.5938,0]
%   'Righting',[1,.2188,0]
%   'Chase',[.8594,0,0]
%   'Touch',[.5,0,0]
%   'WingGrooming',[.6353,.1804,.4510]};

% behavior2color = {'WingGrooming',[1,.4,.755]
%   'Stop',[0,0,.5156]
%   'WingExtension',[0,0,.875]
%   'AttemptCopulation',[0,.2344,1]
%   'Copulation',[0,.9688,1]
%   'Righting',[.3281,1,.6719]
%   'Backup',[0,.75,0]
%   'Crabwalk',[.9,.9,0]
%   'PivotHead',[1,.75,0]
%   'PivotTail',[1,.5938,0]
%   'Touch',[1,.2188,0]
%   'Chase',[.8594,0,0]
%   'Walk',[.5,0,0]
%   'Jump',[.6353,.1804,.4510]};

%  'Crabwalk',[0,.6094,1]

[moviefilestr,trxfilestr,aviname,colors,behavior_cm,...
  zoomflies,nzoomr,nzoomc,boxradius,...
  taillength,fps,maxnframes,firstframes,compression,figpos,movietitle,...
  useVideoWriter,mencoderoptions,mencoder_maxnframes,...
  avifileTempDataFile,titletext,dynamicflyselection,...
  doshowsex,doplotwings,...
  classifierparamsfiles,...
  behaviors,...
  scorefns,...
  score_norm_prctile,...
  behaviorweight,...
  printnone,...
  DEBUG,...
  animaltype,...
  behavior2color] = ...
  myparse(varargin,'moviefilestr','movie.ufmf','trxfilestr','registered_trx.mat',...
  'aviname','','colors',[],'behavior_cm',@jet,...
  'zoomflies',[],'nzoomr',defaults.nzoomr,'nzoomc',defaults.nzoomc,...
  'boxradius',nan,'taillength',defaults.taillength,'fps',defaults.fps,...
  'maxnframes',inf,'firstframes',1,'compression',defaults.compression,...
  'figpos',[50,50,1500,1000],'movietitle','','useVideoWriter',useVideoWriter,...
  'mencoderoptions',mencoderoptions,'mencoder_maxnframes',mencoder_maxnframes,...
  'avifileTempDataFile','',...
  'titletext',true,...
  'dynamicflyselection',true,...
  'doshowsex',true,...
  'doplotwings',true,...
  'classifierparamsfiles',{},...
  'behaviors',{},...
  'scorefns',{},...
  'score_norm_prctile',[],...
  'behaviorweight',[],...
  'printnone',false,...
  'debug',false,...
  'animaltype','fly',...
  'behavior2color',{});

if ~ischar(compression),
  compression = '';
end
if ~isempty(compression) && ~any(strcmpi(compression,allowedcompressions)),
  fprintf('Unknown compressor %s\n',compression);
  compression = '';
end

if isempty(behaviors) || isempty(scorefns),
  classifierparams = ReadClassifierParamsFile(classifierparamsfiles,'isrelativepath',true);
  for j = 1:numel(classifierparams),
    if iscell(classifierparams(j).behaviors.names),
      behaviors{end+1} = sprintf('%s_',classifierparams(j).behaviors.names{:}); %#ok<AGROW>
      behaviors{end} = behaviors{end}(1:end-1);
    else
      behaviors{end+1} = classifierparams(j).behaviors.names; %#ok<AGROW>
    end
    behaviors{end} = regexprep(behaviors{end},'(_|^)([a-z])','${upper($2)}');
    scorefns{end+1} = classifierparams(j).file.scorefilename; %#ok<AGROW>
  end
end

nbehaviors = numel(behaviors);

if isempty(behavior2color),
  behavior2color = [behaviors(:),mat2cell(jet(nbehaviors)*.9,ones(nbehaviors,1),3)];
end

idxremove = ~ismember(lower(behaviors),lower(behavior2color(:,1)));
if any(idxremove),
  fprintf('Ignoring the following behaviors:\n');
  fprintf('  %s\n',behaviors{idxremove});
  behaviors(idxremove) = [];
  scorefns(idxremove) = [];
end

idxremove = ~ismember(lower(behavior2color(:,1)),lower(behaviors));
if any(idxremove),
  fprintf('Ignoring the following behavior colors:\n');
  fprintf('  %s\n',behavior2color{idxremove,1});
  behavior2color(idxremove,:) = [];
end

nbehaviors = numel(behaviors);

%% open files

trxname = fullfile(expdir,trxfilestr);
[trx,trxname,loadsucceeded,timestamps] = load_tracks(trxname);
if ~loadsucceeded,
  error('Could not load trx from file %s',trxname);
end

if ~isempty(moviefilestr)
  moviename = fullfile(expdir,moviefilestr);
  [readframe,nframes,fid] = get_readframe_fcn(moviename);
  im = readframe(1);
  [nr,nc,ncolors] = size(im);
  showgrid = false;
else
  maxx = max([trx.x]+[trx.a]*2);
  maxy = max([trx.y]+[trx.a]*2);
  nr = ceil(maxy);
  nc = ceil(maxx);
  ncolors = 1;
  nframes = max([trx.endframe]);
  readframe = @(x)(uint8(256*ones(nr,nc)));
  fid = 0;
  showgrid = true;
  % grid width
  gridwidth = nanmean([trx.a])*5;

  % create new grid
  handles.guidata.bkgdgrid = nan(2,1);
  xgrid = gridwidth/2:gridwidth:nc;
  xgrid1 = [xgrid;xgrid;nan(1,numel(xgrid))];
  xgrid2 = [zeros(1,numel(xgrid));nr+ones(1,numel(xgrid));nan(1,numel(xgrid))];
  ygrid = gridwidth/2:gridwidth:nr;
  ygrid2 = [ygrid;ygrid;nan(1,numel(ygrid))];
  ygrid1 = [zeros(1,numel(ygrid));nc+ones(1,numel(ygrid));nan(1,numel(ygrid))];
%   for i = 1:numel(handles.guidata.axes_previews),
%     holdstate = ishold(handles.guidata.axes_previews(i));
%     hold(handles.guidata.axes_previews(i),'on');
%     handles.guidata.bkgdgrid(1,i) = plot(handles.guidata.axes_previews(i),xgrid1(:),xgrid2(:),'--','Color',[.7,.7,.7],'LineWidth',.5,'HitTest','off');
%     handles.guidata.bkgdgrid(2,i) = plot(handles.guidata.axes_previews(i),ygrid1(:),ygrid2(:),'--','Color',[.7,.7,.7],'LineWidth',.5,'HitTest','off');
%     if ~holdstate,
%       hold(handles.guidata.axes_previews(i),'off');
%     end
%
%   end

end

nids = length(trx);
trxendframes = [trx.endframe];
trxfirstframes = [trx.firstframe];
T0 = min(trxfirstframes);
T1 = max(trxendframes);
endframes = min(T1,firstframes+maxnframes-1);

scores = cell(1,nids);
score_norms = nan(1,nbehaviors);
for j = 1:nids,
  scores{j} = nan(nbehaviors,trx(j).nframes);
end
isscore = false(1,nbehaviors);
for i = 1:nbehaviors,
  scorefilename = fullfile(expdir,scorefns{i});
  if ~exist(scorefilename,'file'),
    warning('Scores file %s does not exist, skipping',scorefilename);
    continue;
  end
  isscore(i) = true;
  tmp = load(scorefilename);
  score_norms(i) = tmp.allScores.scoreNorm;
  for j = 1:nids,
    scores{j}(i,:) = tmp.allScores.scores{j}(trxfirstframes(j):trxendframes(j));
  end
end
for j = 1:nids,
  scores{j} = scores{j}(isscore,:);
end
behaviors = behaviors(isscore);
scorefns = scorefns(isscore);
nbehaviors = numel(behaviors);

%% output avi name
haveaviname = false;
if ischar(aviname) && ~isempty(aviname)
  [~,ext] = splitext(aviname);
  if strcmpi(ext,'.avi'),
    haveaviname = true;
  end
end
if ~haveaviname,
  fprintf('Choose avi file to output annotated version of %s\n',moviename);
  [base,~] = splitext(movienameonly);
  aviname = [moviepath,'ctraxresults_',base,'.avi'];
  helpmsg = {};
  helpmsg{1} = 'Choose avi file to write annotated movie to.';
  helpmsg{2} = sprintf('Raw movie input: %s',moviename);
  helpmsg{3} = sprintf('Trx file name input: %s',trxname);
  [avinameonly,avipath] = uiputfilehelp('*.mat',sprintf('Choose output avi for %s',movienameonly),aviname,'helpmsg',helpmsg);
  if ~ischar(avinameonly),
    return;
  end
  aviname = [avipath,avinameonly];
end

%% more parameter parsing

doshowsex = doshowsex && isfield(trx,'sex') && ~any(cellfun(@isempty,{trx.sex}));
if doshowsex,
  sexes = {};
  for i = 1:numel(trx),
    sexes = union(sexes,unique(trx(i).sex));
  end
  sexes = upper(sexes);
  sexes = sort(sexes);
  sexmarkers = {'none','*','x','o','+','d','s','p','h'};
  if numel(sexes) > numel(sexmarkers),
    sexmarkers = repmat(sexmarkers,[1,ceil(numel(sexes)/numel(sexmarkers))]);
  end
  sexmarkers = sexmarkers(1:numel(sexes));
end

doplotwings = doplotwings && all(isfield(trx,{'xwingl','ywingl','xwingr','ywingr'}));

if isnan(boxradius),
  meana = 4*nanmean([trx.a]);
  boxradius = ceil(defaults.boxradius*meana);
end

%% which flies to zoom in on

if ~isempty(zoomflies),
  nzoom = numel(zoomflies);
  if max(zoomflies) > nids || min(zoomflies) < 1,
    error('Illegal values for zoomflies');
  end
  nzoomc = ceil(nzoom/nzoomr);
else
  nzoom = nzoomr*nzoomc;
end

% choose some random flies to zoom in on
nframesoverlap = zeros(1,numel(trx));
for i = 1:numel(firstframes),
  nframesoverlap = nframesoverlap + ...
      max(0,min(endframes(i),trxendframes)-max(firstframes(i),trxfirstframes) + 1);
end
fliesmaybeplot = find(nframesoverlap > 0);

if isempty(zoomflies),
  if numel(fliesmaybeplot) < nzoom,
    zoomflies = [fliesmaybeplot,nan(1,nzoom-length(fliesmaybeplot))];
    fprintf('Not enough flies to plot\n');
  else

    if dynamicflyselection,

      % choose flies to plot in each frame
      fliesplotperframe = nan(nzoom,nframes);

      % flies currently chosen
      fliesplotcurr = nan(nzoom,1);
      allfliesplotted = [];

      % loop through all frames

      for i = 1:numel(firstframes),
        for f = firstframes(i):endframes(i),

          % flies that are currently alive
          isalive = trxfirstframes <= f & trxendframes >= f;
          fliesalive = find(isalive);

          % remove newly dead flies
          openzoomboxes = isnan(fliesplotcurr) | ~ismember(fliesplotcurr,fliesalive);
          fliesplotcurr(openzoomboxes) = nan;

          % how many flies do we need to choose
          nflieschoose = nnz(openzoomboxes);

          % flies we can choose to add
          fliesleft = setdiff(fliesalive,fliesplotcurr);

          nfliesshort = max(0,nflieschoose-numel(fliesleft));
          if nfliesshort > 0,
            newfliesplot = fliesleft;
          else
            [~,order] = sort(-[trx(fliesleft).endframe]);
            newfliesplot = fliesleft(order(1:nflieschoose));
          end
          allfliesplotted = [allfliesplotted,newfliesplot]; %#ok<AGROW>
          fliesplotcurr(openzoomboxes) = [newfliesplot,nan(1,nfliesshort)];
          fliesplotperframe(:,f) = fliesplotcurr;

        end
      end

      zoomflies = fliesplotperframe(:,firstframes(1));

    else
      fliesmaybeplot = fliesmaybeplot(randperm(length(fliesmaybeplot)));
      [~,flieswithmostframes] = sort(-nframesoverlap(fliesmaybeplot));
      zoomflies = sort(fliesmaybeplot(flieswithmostframes(1:nzoom)));
    end
  end
elseif nzoom > length(zoomflies),
  zoomflies = [zoomflies,nan(1,nzoom-length(zoomflies))];
end
if nzoomr>0
  zoomflies = reshape(zoomflies,[nzoomr,nzoomc]);
end
rowszoom = floor(nr/nzoomr);

%% choose colors for each behavior

nadd = nbehaviors - size(colors,1);
colors = cat(1,colors,nan(nadd,3));
[ism,idx] = ismember(lower(behaviors),lower(behavior2color(:,1)));
ism = ism & any(isnan(colors),2)';
if any(ism),
  colors(ism,:) = cat(1,behavior2color{idx(ism),2});
end

ism = any(isnan(colors),2)';
if any(ism),
  tmp = behavior_cm(256);
  d = inf(256,1);
  for i = find(~ism),
    d = min(d,sum(abs(bsxfun(@minus,tmp,colors(i,:))),2));
  end
  for i = find(ism),
    [~,j] = max(d);
    colors(i,:) = tmp(j,:);
    d = min(d,sum(abs(bsxfun(@minus,tmp,colors(i,:))),2));
  end
end

%% get the score normalizations

if ~isempty(score_norm_prctile),
  %   score_norms = nan(nbehaviors,1);
  for i = 1:nbehaviors,
    tmp = [];
    for j = 1:nids,
      tmp = [tmp,scores{j}(i,:)]; %#ok<AGROW>
    end
    score_norms(i) = prctile(abs(tmp),score_norm_prctile);
  end
end

for fly = 1:nids,
  for behaviori = 1:nbehaviors,
    scores{fly}(behaviori,:) = min(1,max(-1,scores{fly}(behaviori,:)/score_norms(behaviori)));
  end
end

%% choose weights for each behavior: upweight if less common

if isempty(behaviorweight),
  behaviorweight = nan(nbehaviors,1);
  for behaviori = 1:nbehaviors,
    tmp = [];
    for j = 1:nids,
      tmp = [tmp,scores{j}(behaviori,:)]; %#ok<AGROW>
    end
    behaviorweight(behaviori) = nnz(tmp<=0);
  end
  behaviorweight = behaviorweight / sum(behaviorweight);
end
if numel(behaviorweight) == 1,
  behaviorweight = ones(nbehaviors,1);
end
behaviorweight = behaviorweight(:);

%% initialize plots

figure(1);
clf;
hold on;
hax = gca;
set(hax,'position',[0,0,1,1]);
axis off;
isdisplay = ispc || ~strcmpi(get(1,'XDisplay'),'nodisplay');

% corners of zoom boxes in plotted image coords
x0 = nc+(0:nzoomc-1)*rowszoom+1;
y0 = (0:nzoomr-1)*rowszoom+1;
x1 = x0 + rowszoom - 1;
y1 = y0 + rowszoom - 1;

% relative frame offset
nframesoff = trxfirstframes - 1;

% pre-allocate
himzoom = zeros(nzoomr,nzoomc);
htail = zeros(2*taillength+1,nids);
htri = zeros(1,nids);
hwing = zeros(1,nids);
hsexmarker = zeros(1,nids);
hsexmarker_zoom = zeros(1,nids);
scalefactor = rowszoom / (2*boxradius+1);
hzoom = zeros(nzoomr,nzoomc);
hzoomwing = zeros(nzoomr,nzoomc);
htextzoom = zeros(nzoomr,nzoomc);
hzoombox = zeros(nzoomr,nzoomc);
if strcmpi(animaltype,'larvacontour')
  hextra = zeros(nids,2);
  hextrazoom = zeros(nzoomr,nzoomc,2);
end

% initialize plots
him = image([1,nc],[1,nr],zeros([nr,nc,3]));
axis image;
axis tight;

if nzoom >0
  axis([.5,x1(end)+.5,.5,y1(end)+.5]);
else
  axis([.5,nc+.5,.5,nr+.5]);
end
axis off;

if titletext && isdisplay,
  htext = text(.5,.5,'frame info','Parent',hax,'BackgroundColor','k','Color','g','VerticalAlignment','bottom','interpreter','none');
end

for i = 1:nzoomr,
  for j = 1:nzoomc,
    himzoom(i,j) = image([x0(j),x1(j)],[y0(i),y1(i)],repmat(uint8(123),[boxradius*2+1,boxradius*2+1,3]));
  end
end

for fly = 1:nids,
  htail(:,fly) = plot(nan(2,2*taillength+1),nan(2,2*taillength+1),'k-');
  if doshowsex,
    hsexmarker(fly) = plot(nan,nan,'k.');
  end
  htri(fly) = plot(nan,nan,'k-');
  if doplotwings,
    hwing(fly) = plot(nan,nan,'k.-');
  end
  if strcmpi(animaltype,'larvacontour')
    hextra(fly,1) = plot(nan,nan,'Marker','.','LineStyle','-');
    hextra(fly,2) = plot(nan,nan,'Marker','o','LineStyle','None');
  end

end

for i = 1:nzoomr,
  for j = 1:nzoomc,
    hzoom(i,j) = plot(nan,nan,'k-');
    hzoomwing(i,j) = plot(nan,nan,'k.-');
    xtmp = x0(j) + [0,1]*(x1(j)-x0(j)) + [1,-1];
    ytmp = y0(i) + [0,1]*(y1(i)-y0(i)) + [1,-1];
    hzoombox(i,j) = plot(xtmp([1,2,2,1,1]),ytmp([1,1,2,2,1]),'k-','LineWidth',2,...
      'Visible','off');
    htextzoom(i,j) = text((x0(j)+x1(j))/2,.95*y0(i)+.05*y1(i),'',...
      'color','k','horizontalalignment','center',...
      'verticalalignment','bottom','fontweight','bold',...
      'FontSize',18);
    if doshowsex,
      hsexmarker_zoom(i,j) = plot(nan,nan,'k.');
    end

    if strcmpi(animaltype,'larvacontour')
      hextrazoom(i,j,1) = plot(nan,nan,'Marker','.','LineStyle','-');
      hextrazoom(i,j,2) = plot(nan,nan,'Marker','o','LineStyle','None');
    end

  end
end

mencoder_nframes = 0;

tailx = nan(2*taillength+1,nids);
taily = nan(2*taillength+1,nids);
tailcolor = zeros(2*taillength+1,3,nids);
tailisbehavior = false(2*taillength+1,nids);

%%

for segi = 1:numel(firstframes),
  firstframe = firstframes(segi);
  endframe = endframes(segi);

  frame = firstframe - 1;
  idx = frame - nframesoff;
  tailend = 1;
  for fly = 1:nids,
    tailx(:,fly) = padgrab(trx(fly).x,nan,1,1,idx(fly)-taillength,idx(fly)+taillength);
    taily(:,fly) = padgrab(trx(fly).y,nan,1,1,idx(fly)-taillength,idx(fly)+taillength);
    tmpscores = padgrab(scores{fly},0,1,nbehaviors,idx(fly)-taillength,idx(fly)+taillength);
    tailisbehavior(:,fly) = any(tmpscores > 0,1);
    w = max(0,bsxfun(@times,tmpscores,behaviorweight));
    w = permute(w,[1,3,2]);
    z = sum(w,1);
    z(z==0) = 1;
    w = bsxfun(@rdivide,w,z);
    tailcolortmp = sum(bsxfun(@times,w,colors),1);
    tailcolortmp = max(0,min(1,tailcolortmp));
    tailcolor(:,:,fly) = permute(tailcolortmp,[3,2,1]);
    for i = 1:taillength*2,
      set(htail(i,fly),'XData',tailx([i,i+1],fly),'YData',taily([i,i+1],fly),'Color',tailcolor(i,:,fly));
      if tailisbehavior(i,fly),
        set(htail(i,fly),'LineWidth',2);
      else
        set(htail(i,fly),'LineWidth',1);
      end
    end
  end
  %tailcolor = max(0,min(1,tailcolor));

  for frame = firstframe:endframe,
    if mod(frame - firstframe,20) == 0,
      fprintf('frame %d\n',frame);
    end

    % relative frame
    idx = frame - nframesoff;

    isalive = frame >= trxfirstframes & ...
      frame <= trxendframes;

    % draw the unzoomed image
    im = uint8(readframe(frame));
    if ncolors == 1,
      im = repmat(im,[1,1,3]);
    end
    set(him,'cdata',im);

    % draw frame number text box
    % text doesn't show up in no display mode
    if titletext && isdisplay,
      framestr = sprintf('Frame %d, t = %.2f s',frame,timestamps(frame)-timestamps(1));
      if ~isempty(movietitle),
        framestr = {framestr,movietitle}; %#ok<AGROW>
      end
      set(htext,'String',framestr);
    end

    % draw the zoomed image
    if dynamicflyselection && exist('fliesplotperframe','var'),
      zoomflies = reshape(fliesplotperframe(:,frame),[nzoomr,nzoomc]);
    end
    for i = 1:nzoomr,
      for j = 1:nzoomc,
        fly = zoomflies(i,j);

        % fly not visible?
        if isnan(fly) || ~isalive(fly),
          set(himzoom(i,j),'cdata',repmat(uint8(123),[boxradius*2+1,boxradius*2+1,3]));
          continue;
        end

        % grab a box around (x,y)
        x = round(trx(fly).x(idx(fly)));
        y = round(trx(fly).y(idx(fly)));
        boxradx1 = min(boxradius,x-1);
        boxradx2 = min(boxradius,size(im,2)-x);
        boxrady1 = min(boxradius,y-1);
        boxrady2 = min(boxradius,size(im,1)-y);
        box = uint8(zeros(2*boxradius+1));
        box(boxradius+1-boxrady1:boxradius+1+boxrady2,...
          boxradius+1-boxradx1:boxradius+1+boxradx2) = ...
          im(y-boxrady1:y+boxrady2,x-boxradx1:x+boxradx2);
        set(himzoom(i,j),'cdata',repmat(box,[1,1,3]));

      end
    end

    idxfuture = frame + taillength - nframesoff;
    colorfuture = zeros(nids,3);
    isalivefuture = frame + taillength >= trxfirstframes & ...
      frame + taillength <= trxendframes;

    tailx = cat(1,tailx(2:end,:),nan(1,nids));
    taily = cat(1,taily(2:end,:),nan(1,nids));
    tailisbehavior = cat(1,tailisbehavior(2:end,:),false(1,nids));
    for fly = 1:nids,
      if ~isalivefuture(fly), continue; end
      w = max(0,scores{fly}(:,idxfuture(fly))).*behaviorweight;
      tailisbehavior(end,fly) = any(scores{fly}(:,idxfuture(fly))>0,1);
      z = sum(w);
      if z == 0,
        colorfuture(fly,:) = [0,0,0];
      else
        w = w /z;
        %colortmp(fly,:) = sum(bsxfun(@times,w,colors),1)/z*max(scores{fly}(:,idx(fly)));
        colorfuture(fly,:) = sum(bsxfun(@times,w,colors),1);
      end
      tailx(end,fly) = trx(fly).x(idxfuture(fly));
      taily(end,fly) = trx(fly).y(idxfuture(fly));
    end
    colorfuture = max(0,min(1,colorfuture));

    tailcolor = cat(1,tailcolor(2:end,:,:),permute(colorfuture,[3,2,1]));
    colortmp = permute(tailcolor(taillength+1,:,:),[3,2,1]);
    isbehavior = tailisbehavior(taillength+1,:);

    % sanity check
    if DEBUG,
      tmpcolortmp = zeros(nids,3);
      tmptailx = nan(2*taillength+1,nids);
      tmptaily = nan(2*taillength+1,nids);
      tmpisbehavior = false(1,nids);
      for fly = 1:nids,
        if ~isalive(fly), continue; end
        w = max(0,scores{fly}(:,idx(fly))).*behaviorweight;
        z = sum(w);
        if z == 0,
          tmpcolortmp(fly,:) = [0,0,0];
        else
          w = w /z;
          %colortmp(fly,:) = sum(bsxfun(@times,w,colors),1)/z*max(scores{fly}(:,idx(fly)));
          tmpcolortmp(fly,:) = sum(bsxfun(@times,w,colors),1);
        end
        tmptailx(:,fly) = padgrab(trx(fly).x,nan,1,1,idx(fly)-taillength,idx(fly)+taillength);
        tmptaily(:,fly) = padgrab(trx(fly).y,nan,1,1,idx(fly)-taillength,idx(fly)+taillength);

        tmpisbehavior(fly) = any(scores{fly}(:,idx(fly)) > 0);

      end
      if ~all(isnan(tailx(:)) & isnan(tmptailx(:)) | (tailx(:) == tmptailx(:)))
        error('tailx does not pass sanity check');
      end
      if ~all(isnan(taily(:)) & isnan(tmptaily(:)) | (taily(:) == tmptaily(:)))
        error('taily does not pass sanity check');
      end
      if ~all(colortmp(:) == tmpcolortmp(:)),
        error('colortmp does not pass sanity check');
      end
      if ~all(isbehavior(:) == tmpisbehavior(:)),
        error('isbehavior does not pass sanity check');
      end
    end

    % plot the zoomed out position
    for fly = 1:nids,
      if isalive(fly),
        set(htail(tailend,fly),'XData',tailx(end-1:end,fly),...
          'YData',taily(end-1:end,fly),...
          'Color',tailcolor(end-1,:,fly));
        if tailisbehavior(end-1,fly),
          set(htail(tailend,fly),'LineWidth',2);
        else
          set(htail(tailend,fly),'LineWidth',1);
        end
%         i0 = max(1,idx(fly)-taillength);
%         set(htail(fly),'xdata',trx(fly).x(i0:idx(fly)),...
%           'ydata',trx(fly).y(i0:idx(fly)),'color',colortmp(fly,:));
        if doshowsex,
          sexi = find(strcmpi(trx(fly).sex{idx(fly)},sexes),1);
          sexmarker = sexmarkers{sexi};
          sexcolorcurr = [1,1,1];
          set(hsexmarker(fly),'xdata',trx(fly).x(idx(fly)),...
            'ydata',trx(fly).y(idx(fly)),...
            'color',sexcolorcurr,...
            'marker',sexmarker,...
            'markerfacecolor',sexcolorcurr);
        end
        if strcmpi(animaltype,'fly'),
          updatefly(htri(fly),trx(fly),idx(fly));
        elseif strcmpi(animaltype,'larvacontour'),
          pos.xcontour = trx(fly).xcontour{idx(fly) };
          pos.ycontour = trx(fly).ycontour{idx(fly) };
          pos.xspine = trx(fly).xspine(:,idx(fly) );
          pos.yspine = trx(fly).yspine(:,idx(fly));

          updatelarvacontour(htri(fly),hextra(fly,:),pos);

        end
        set(htri(fly),'Color',colortmp(fly,:));
        if strcmpi(animaltype,'larvacontour'),
          set(hextra(fly,:),'Color',colortmp(fly,:));
        end
        if isbehavior(fly),
          set(htri(fly),'Linewidth',2);
        else
          set(htri(fly),'Linewidth',1);
        end
        if doplotwings,
          xwing = [trx(fly).xwingl(idx(fly)),trx(fly).x(idx(fly)),trx(fly).xwingr(idx(fly))];
          ywing = [trx(fly).ywingl(idx(fly)),trx(fly).y(idx(fly)),trx(fly).ywingr(idx(fly))];
          set(hwing(fly),'XData',xwing,'YData',ywing,'Color',colortmp(fly,:));
        end
      else
        set(htail(:,fly),'xdata',nan,'ydata',nan);
        set(htri(fly),'xdata',[],'ydata',[]);
        if doshowsex,
          set(hsexmarker(fly),'xdata',[],'ydata',[]);
        end
        if doplotwings,
          set(hwing(fly),'XData',[],'YData',[]);
        end
        if strcmpi(animaltype,'larvacontour'),
          set(hextra(fly,1),'XData',[],'YData',[]);
          set(hextra(fly,2),'XData',[],'YData',[]);
        end
      end
    end
    tailend = tailend + 1;
    if tailend >= 2*taillength+1,
      tailend = 1;
    end

    % plot the zoomed views
    for i = 1:nzoomr,
      for j = 1:nzoomc,
        fly = zoomflies(i,j);
        if ~isnan(fly) && isalive(fly),
          x = trx(fly).x(idx(fly));
          y = trx(fly).y(idx(fly));
          x = boxradius + (x - round(x))+.5;
          y = boxradius + (y - round(y))+.5;
          x = x * scalefactor;
          y = y * scalefactor;
          x = x + x0(j) - 1;
          y = y + y0(i) - 1;
          a = trx(fly).a(idx(fly))*scalefactor;
          b = trx(fly).b(idx(fly))*scalefactor;
          theta = trx(fly).theta(idx(fly));
          if strcmpi(animaltype,'larvacontour')
            pos.xcontour = trx(fly).xcontour{idx(fly)};
            pos.ycontour = trx(fly).ycontour{idx(fly)};
            pos.xspine = trx(fly).xspine(:,idx(fly));
            pos.yspine = trx(fly).yspine(:,idx(fly));
            cx = trx(fly).x(idx(fly));
            cy = trx(fly).y(idx(fly));
            pos.xcontour = pos.xcontour + boxradius - round(cx) + 0.5;
            pos.ycontour = pos.ycontour + boxradius - round(cy) + 0.5;
            pos.xspine = pos.xspine + boxradius - round(cx) + 0.5;
            pos.yspine = pos.yspine + boxradius - round(cy) + 0.5;

            pos.xcontour = pos.xcontour*scalefactor + x0(j) -1;
            pos.ycontour = pos.ycontour*scalefactor + y0(i) -1;
            pos.xspine = pos.xspine*scalefactor + x0(j) -1;
            pos.yspine = pos.yspine*scalefactor + y0(i) -1;

          end

          labelis = find(scores{fly}(:,idx(fly)) > 0);
          tmpscores = scores{fly}(labelis,idx(fly));
          s = cell(1,numel(labelis));
          for iii = 1:numel(labelis),
            s{iii} = sprintf('\\color[rgb]{%f %f %f}%s',tmpscores(iii)*colors(labelis(iii),:),behaviors{labelis(iii)});
          end

          % try to fit on two lines, with at most two behaviors per-line
          if numel(labelis) > 2,
            nperline = 2;
            s0 = s;
            s = cell(1,ceil(numel(labelis)/nperline));
            for linei = 1:numel(s),
              s{linei} = sprintf('%s, ',s0{2*(linei-1)+1:min(numel(labelis),2*linei)});
              s{linei} = s{linei}(1:end-2);
            end
          end

          if printnone && isempty(s),
            s = {'None'};
          end

%           if doshowsex,
%             s = sprintf('%d, %s',fly,trx(fly).sex{idx(fly)});
%           else
%             s = sprintf('%d',fly);
%           end
          if doplotwings,
            xwingl = trx(fly).xwingl(idx(fly)) - round(trx(fly).x(idx(fly))) + boxradius + .5;
            ywingl = trx(fly).ywingl(idx(fly)) - round(trx(fly).y(idx(fly))) + boxradius + .5;
            xwingl = xwingl * scalefactor;
            ywingl = ywingl * scalefactor;
            xwingl = xwingl + x0(j) - 1;
            ywingl = ywingl + y0(i) - 1;
            xwingr = trx(fly).xwingr(idx(fly)) - round(trx(fly).x(idx(fly))) + boxradius + .5;
            ywingr = trx(fly).ywingr(idx(fly)) - round(trx(fly).y(idx(fly))) + boxradius + .5;
            xwingr = xwingr * scalefactor;
            ywingr = ywingr * scalefactor;
            xwingr = xwingr + x0(j) - 1;
            ywingr = ywingr + y0(i) - 1;
            xwing = [xwingl,x,xwingr];
            ywing = [ywingl,y,ywingr];
          end

          if strcmpi(animaltype,'fly'),
            updatefly(hzoom(i,j),x,y,theta,a,b);
          elseif strcmpi(animaltype,'larvacontour'),
            updatelarvacontour(hzoom(i,j),hextrazoom(i,j,:),pos);
          end
          if doplotwings,
            set(hzoomwing(i,j),'XData',xwing,'YData',ywing,'Color',colortmp(fly,:));
          end
          set(htextzoom(i,j),'string',s,'color',colortmp(fly,:));
          if isbehavior(fly),
            set(hzoombox(i,j),'Visible','on','color',colortmp(fly,:));
            set(hzoom(i,j),'Linewidth',2);
          else
            set(hzoombox(i,j),'Visible','off');
            set(hzoom(i,j),'Linewidth',1);
          end
          if doshowsex,
            sexcolorcurr = [1,1,1];
            sexi = find(strcmpi(trx(fly).sex{idx(fly)},sexes),1);
            sexmarker = sexmarkers{sexi};
            set(hsexmarker_zoom(i,j),'Visible','on','XData',x,'YData',y,'color',sexcolorcurr,'MarkerFaceColor',sexcolorcurr,'Marker',sexmarker);
          end

%           if doshowsex,
%             sexi = find(strcmpi(trx(fly).sex{idx(fly)},sexes),1);
%             sexmarker = sexmarkers{sexi};
%             set(htextzoom(i,j),'xdata',x,...
%               'ydata',y,...
%               'color',colortmp(fly,:),...
%               'marker',sexmarker,...
%               'markerfacecolor',colortmp(fly,:));
%           end
          set(hzoom(i,j),'color',colortmp(fly,:));
        else
          set(hzoom(i,j),'xdata',[],'ydata',[]);
          if doplotwings,
            set(hzoomwing(i,j),'XData',[],'YData',[]);
          end
          set(htextzoom(i,j),'string','');
          if doshowsex,
            set(hsexmarker_zoom(i,j),'Visible','off');
          end
%           if doshowsex,
%             set(htextzoom(i,j),'xdata',[],'ydata',[]);
%           end
        end
      end
    end

    if frame == firstframes(1),
      if ~isempty(figpos),
        set(1,'Position',figpos);
      else
        input('Resize figure 1 to the desired size, hit enter when done.');
        figpos = get(1,'Position');
      end
      if ~DEBUG,
        set(1,'visible','off');
        if useVideoWriter,
          if strcmpi(compression,'None') || strcmpi(compression,'Uncompressed AVI'),
            profile = 'Uncompressed AVI';
          else
            profile = 'Motion JPEG AVI';
          end
          aviobj = VideoWriter(aviname,profile); %#ok<TNMLP>
          set(aviobj,'FrameRate',fps);
          if ~strcmpi(profile,'Uncompressed AVI'),
            set(aviobj,'Quality',100);
          end
          open(aviobj);
        else
          if isempty(avifileTempDataFile),
            aviobj = avifile(aviname,'fps',fps,'quality',100,'compression',compression);  %#ok<REMFF1>
          else
            aviobj = myavifile(aviname,'fps',fps,'quality',100,'compression',compression,...
              'TempDataFile',avifileTempDataFile);
            fprintf('Temporary data file for avi writing: %s\n',aviobj.TempDataFile);
          end
        end
      end
    end

    if ~DEBUG,
    if frame == firstframes(1),
      fr = getframe_invisible(hax);
      [height,width,~] = size(fr);
      gfdata = getframe_initialize(hax);
      [fr,height,width] = getframe_invisible_nocheck(gfdata,[height,width],false,false);
%       height = ceil(height/4)*4;
%       width = ceil(width/4)*4;
%       fr = getframe_invisible(hax,[height,width]);
    else
      fr = getframe_invisible_nocheck(gfdata,[height,width],false);
    end
    if useVideoWriter,
      writeVideo(aviobj,fr);
    else
      aviobj = addframe(aviobj,fr);
    end
    set(1,'Position',figpos);
    else
      drawnow;
%       if any(isbehavior(:)),
%         input('');
%       end
    end
    %drawnow;

  end

end

%% clean up

fprintf('Finishing AVI...\n');
if ~DEBUG,
if useVideoWriter,
  close(aviobj);
else
  aviobj = close(aviobj); %#ok<NASGU>
end
fprintf('Cleanup...\n');
getframe_cleanup(gfdata);
end
if fid > 0,
  fclose(fid);
end

succeeded = true;

author: kristinbranson

repo: JAABA


text:
%% Intialise Parameters
function init_params()
global El_lim_low El_lim_up LB UB El_range PopulationSize P Generations iteration;
global El vel NumVer pbest_loc gbest_loc pbest_val EV gbest_val_record;
global FitVal gbest_val ev av d n p t delta_t k_pb_i k_pb_f k_gb_i k_gb_f w
El_lim_low(1,1:NumVer)=LB;
El_lim_up(1,1:NumVer)=UB;
El_range(1,1:NumVer)= El_lim_up - El_lim_low;
P=PopulationSize;% number of particles
iteration=Generations;% number of iterations
% Intialise Pbest of all agents & Gbest
El =zeros(1,NumVer,P);
vel =zeros(1,NumVer,P);
pbest_loc=zeros(1,NumVer,P);
gbest_loc=zeros(1,NumVer);
pbest_val=zeros(P);
EV=P*iteration;
gbest_val_record(1:EV) = 0;
FitVal = 0.0;
gbest_val = 0.0;
ev = 0;
av = 0;
d=0; n=0; p=0; t=0;
delta_t=1;
% Set parameters of the optimization problem
k_pb_i = 2.5;
k_pb_f = 0.5;
k_gb_i = 0.5;
k_gb_f = 2.5;
w = 0.9;
end

author: Sable

repo: mcbench-benchmarks


text:
Q = [1 1];                      % State weights
R = 0.5;                        % du weights
Ru = 0.5;                       % u weights
B = [0; 1];
C = eye(Nvar);
D = 0;
LB = -20*ones(N,1);        % Lower bound of control input
UB = 20*ones(N,1);         % Upper bound of control input

author: eurika-kaiser

repo: SINDY-MPC


text:
%script test recursive
%%%%% quantitative differences: because of the way the level is updated
%%%%% (should be analytically equivalent but not numerically). Small
%%%%% deiations for 2-ToM (but accumulates) larger for 3-ToM
x=zeros(2,101); xtest=x;
y1=VBA_random ('Bernoulli', 0.7, 100, 1);
inF.player=1;
inF.lev=0;
Par=0; %theta
for i=1:100
    x(:,i+1)=RecToMfunction(x(:,i),Par,y1(i),inF);
    xtest(:,i+1)=evolution0bisND(xtest(:,i),Par,y1(i),inF);
end
%%
N=1000
xtest=zeros(8,N+1); x=zeros(9,N+1);
y1=VBA_random('Bernoulli', 0.7, N, 1);
y2=VBA_random('Bernoulli', 0.5, N, 1);
HS=cat(3,[1,0;0,1],[0,1;1,0]);
inF.player=1;
inF.lev=1;
inF.game=HS;
inF.fobs=@ObsRecGen;
inF.indParev=1;
inF.indParobs=1;
inF.dummyPar=[1;0];
in={HS,1};
Par=0; %theta
for i=1:N
    [x(:,i+1),indlev]=RecToMfunction(x(:,i),Par,[y2(i);y1(i)],inF);
    xtest(:,i+1)=evolution1bisND(xtest(:,i),Par,[y2(i);y1(i)],in);
end

%%%%compare
indcorr1=[1:2,6,8,7,9,3:4];
xcomp=x(indcorr1,:);
find(abs(xcomp-xtest)>1e-15)
%%
N=50;

y1=VBA_random('Bernoulli',0.7,100,1);
y2=VBA_random('Bernoulli',0.5,100,1);
%%
xtest=zeros(23,N+1); x=zeros(26,N+1);
HS=cat(3,[1,0;0,1],[0,1;1,0]);
inF.player=1;
inF.lev=2;
inF.game=HS;
inF.fobs=@ObsRecGen;
inF.indParev=1;
inF.indParobs=1;
inF.dummyPar=[1;0];
in={HS,1};
Par=0; %theta
for i=1:N
    [x(:,i+1),indlev]=RecToMfunction(x(:,i),Par,[y2(i);y1(i)],inF);
    xtest(:,i+1)=evolution2bisND(xtest(:,i),Par,[y2(i);y1(i)],in);
end
indcorr2=[2:5,10+indcorr1,20,21,7,9,8,10,23,25,24,26,1];
xcomp=x(indcorr2,:);
xcomp(end,:)=-xcomp(end,:);
%%%%compare
find(abs(xcomp-xtest)>1e-5)
%%
N=50;
xtest=zeros(53,N+1); x=zeros(60,N+1);
y1=VBA_random('Bernoulli',0.7,100,1);
y2=VBA_random('Bernoulli',0.5,100,1);
HS=cat(3,[1,0;0,1],[0,1;1,0]);
inF.player=1;
inF.lev=3;
inF.game=HS;
inF.fobs=@ObsRecGen;
inF.indParev=1;
inF.indParobs=1;
inF.dummyPar=[1;0];
in={HS,1};
Par=0; %theta
for i=1:3
    [x(:,i+1),indlev]=RecToMfunction(x(:,i),Par,[y2(i);y1(i)],inF);
    xtest(:,i+1)=evolution3bisND(xtest(:,i),Par,[y2(i);y1(i)],in);
end
 indcorr3=[3:6,11+indcorr1,21,22,27+indcorr2,54,55,8,10,9,11,24,26,25,27,57,59,58,60,1,2];
 xcomp=x(indcorr3,:);
 xcomp(37,:)=-xcomp(37,:); %estimation level 2

%%%%compare
find(abs(xcomp-xtest)>1e-4)

%%
N=10;
inF.lev=7;
%xtest=zeros(53,N+1);
x=zeros(sizeXrec(inF.lev,2),N+1);
y1=VBA_random('Bernoulli',0.7,100,1);
y2=VBA_random('Bernoulli',0.5,100,1);

HS=cat(3,[1,0;0,1],[0,1;1,0]);
inF.player=1;

inF.game=HS;
inF.fobs=@ObsRecGen;
inF.indParev=1;
inF.indParobs=1;
inF.dummyPar=[1;0];
in={HS,1};
Par=0; %theta
for i=1:N
    [x(:,i+1),indlev]=RecToMfunction(x(:,i),Par,[y2(i);y1(i)],inF);
   % xtest(:,i+1)=evolution3bisND(xtest(:,i),Par,[y2(i);y1(i)],in);
end
%  indcorr3=[3:6,11+indcorr1,21,22,27+indcorr2,54,55,8,10,9,11,24,26,25,27,57,59,58,60,1,2];
%  xcomp=x(indcorr3,:);
%  xcomp(37,:)=-xcomp(37,:); %estimation level 2
%
% %%%%compare
% find(abs(xcomp-xtest)>1e-10)
%%
N=100;
 x=zeros(264,N+1);
y1=VBA_random('Bernoulli',0.7,100,1);
y2=VBA_random('Bernoulli',0.5,100,1);
HS=cat(3,[1,0;0,1],[0,1;1,0]);
inF.player=1;
inF.lev=5;
inF.game=HS;
inF.fobs=@ObsRecGen;
inF.indParev=1;
inF.indParobs=1;
inF.dummyPar=[1;0];
in={HS,1};
Par=0; %theta
for i=1:N
    [x(:,i+1),indlev]=RecToMfunction(x(:,i),Par,[y2(i);y1(i)],inF);
end

author: MBB-team

repo: VBA-toolbox


text:
function M = makeSincMovie

%   Copyright 2008-2009 The MathWorks, Inc.

% Make suitable X and Y matricies to plot the sinc function over
[X, Y] = meshgrid(linspace(-3*pi, 3*pi, 50), linspace(-3*pi, 3*pi, 50));
% From X & Y make the matrix R of radius
R = sqrt(X.^2 + Y.^2);
% Calculate the sinc function from R
F = sin(R)./R;
% Make an amplitude vector that represents the various amplitudes we will
% movie the function over
A = linspace(-1, 1, 30);
% Iterate over each amplitude
for j = 1:numel(A)
    % Plot the surface with colouring defined by the function with no
    % amplitude modification
    surf(X, Y, A(j)*F, F);
    % Always use the same axis otherwise the movie will look funny
    axis([-10 10 -10 10 -1 1]);
    % Then turn off the axis - don't want tosee them in the movie. Also
    % remove the lines around each face of the surface
    axis off
    shading flat
    % Finally get this frame and the symmetric one so that when we show the
    % whole movie it ends at the same point it began
    M(j) = getframe; %#ok<AGROW>
    M(2*numel(A) - j) = M(j); %#ok<AGROW>
end


text:
function obj = setParent(obj, parent)
    % set the parent link frame
    %
    % Parameters:
    % parent: the parent link frame object @type char

    assert(ischar(parent),...
        'The parent link name must be a character vector.');
    obj.Parent = parent;

end

author: ayonga

repo: frost-dev

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


text:
%sim_dDiv.m
%Jamie Near, 2020.
%
% USAGE:
% d_out = sim_dDiv(d_in,factor)
%
% DESCRIPTION:
% Divide a density matrix by a scalar.  This function is necessary becuase
% the density matrix is a cell array, so cannot be scaled using simple
% division.
%
% INPUTS:
% d_in      = input density matrix to be divided.
% factor    = scalar factor to divide by.
%
% OUTPUTS:
% d_out     = output, result of d_in/factor.

function d_out = sim_dDiv(d_in,factor)

d_out=cell(size(d_in));
for m=1:length(d_in)
    d_out{m}=d_in{m}./factor;
end

author: CIC-methods

repo: FID-A

save_path: github-repos/MATLAB/CIC-methods-FID-A


text:
% Author: Guosheng Lin (guosheng.lin@gmail.com)

% perform segmentation prediction on user provided images.
% specify the location of your images, e.g., using the following
% folder which contains serveral example images:
% ds_config.img_data_dir='../datasets/custom_data';

function demo_refinenet_test_example_eyth()
 folders = {'vid4','vid6','vid9'};

 for i = 1:length(folders)
    rng('shuffle');

    addpath('./my_utils');
    dir_matConvNet='../libs/matconvnet/matlab';
    run(fullfile(dir_matConvNet, 'vl_setupnn.m'));

    run_config=[];
    ds_config=[];

    run_config.use_gpu=true;
    % run_config.use_gpu=false;
    run_config.gpu_idx=1;

    % result dir:
    result_name=[char(folders(i))]

    result_dir=fullfile('../cache_data', 'eyth_on_eyth_val/val_1', result_name);

    % the folder that contains testing images:
    ds_config.img_data_dir=strcat('../datasets/eyth/JPEGImages/test/',char(folders(i)));

    % using a trained model which is trained on VOC 2012
    % run_config.trained_model_path='../model_trained/refinenet_res101_voc2012.mat';
    % ds_config.class_info=gen_class_info_voc();

    % using the object parsing model
    run_config.trained_model_path='../model_trained/refinenet_res101_eyth.mat';
    ds_config.class_info=gen_class_info_ego();

    % for voc trained model, control the size of input images
    run_config.input_img_short_edge_min=450;
    run_config.input_img_short_edge_max=600;

    % set the input image scales, useful for multi-scale evaluation
    % e.g. using multiple scale settings (1.0 0.8 0.6) and average the resulting score maps.
    run_config.input_img_scale=1.0;

    run_config.gen_net_opts_fn=@gen_net_opts_model_type1;

    run_config.run_evaonly=true;
    ds_config.use_custom_data=true;
    ds_config.use_dummy_gt=true;
    run_config.use_dummy_gt=ds_config.use_dummy_gt;

    ds_config.ds_name='tmp_data';
    run_config.root_cache_dir=result_dir;
    mkdir_notexist(run_config.root_cache_dir);

    run_config.model_name=result_name;

    diary_dir=run_config.root_cache_dir;
    mkdir_notexist(diary_dir);
    diary(fullfile(diary_dir, 'output.txt'));
    diary on

    run_dir_name=fileparts(mfilename('fullpath'));
    [~, run_dir_name]=fileparts(run_dir_name);
    run_config.run_dir_name=run_dir_name;
    run_config.run_file_name=mfilename();

    ds_info=gen_dataset_info(ds_config);
    my_diary_flush();

    train_opts=run_config.gen_net_opts_fn(run_config, ds_info.class_info);

    imdb=my_gen_imdb(train_opts, ds_info);

    data_norm_info=[];
    data_norm_info.image_mean=128;

    imdb.ref.data_norm_info=data_norm_info;

    if run_config.use_gpu
        gpu_num=gpuDeviceCount;
        if gpu_num>=1
            gpuDevice(run_config.gpu_idx);
        else
            error('no gpu found!');
        end
    end

    [net_config, net_exp_info]=prepare_running_model(train_opts);

    my_net_tool(train_opts, imdb, net_config, net_exp_info);

    fprintf('\n\n--------------------------------------------------\n\n');
    disp('results are saved in:');
    disp(run_config.root_cache_dir);

    my_diary_flush();
    diary off

 end
end

author: aurooj

repo: Hand-Segmentation-in-the-Wild


text:
function p = high_card_fun ( m, n )

%*****************************************************************************80
%
%% HIGH_CARD_FUN estimates the value of breaking the deck at location K.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    23 May 2012
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, integer M, the number of cards in the deck.
%
%    Input, integer N, the number of trials.
%
%    Output, real P(M), the estimated probability of picking the correct
%    high card by discarding so many cards and taking the next card that
%    is higher.
%
  fprintf ( 1, '\n' );
  fprintf ( 1, 'HIGH_CARD_FUN\n' );
  fprintf ( 1, '  Using N=%d cards and T=%d trials,\n', m, n );
  fprintf ( 1, '  estimate the chances of correctly picking the highest card\n' );
  fprintf ( 1, '  by looking at cards 0 through K-1, and then taking the first\n' );
  fprintf ( 1, '  subsequent card that is bigger.\n' );

  p = zeros ( m, 1 );

  parfor i = 1 : m

    p(i) = high_card_simulation ( m, n, i - 1 );

  end

  return
end
function p = high_card_simulation ( deck_size, trial_num, skip_num )

%*****************************************************************************80
%
%% HIGH_CARD_SIMULATION simulates a game of choosing the highest card in a deck.
%
%  Discussion:
%
%    You are given a deck of DECK_SIZE cards.
%
%    Your goal is to select the high card.  For convenience, we can assume
%    the cards are a permutation of the integers from 1 to DECK_SIZE, but in
%    fact the user mustn't see such values or else it's obvious which is the
%    largest card.
%
%    However, your choice is made under the following rules:  You may turn over
%    one card at a time.  When a card is turned over, you may declare that to be
%    your choice, or else turn over another card.  If you have not chosen a card
%    by the end, then your choice is the final card.
%
%    If you have no idea what to do, and simply decide in advance to pick
%    a card "at random", that is, for example, you decide to pick the 15th card
%    before having seen any cards, then your probability of winning is 1/DECK_SIZE.
%
%    The question is, can you do better than that?
%
%    Your strategy is as follows: always look at the first SKIP_NUM cards without
%    choosing them.  Then choose the very next card you encounter that is larger
%    than the cards you skipped.
%
%    Using this program, you can easily see that skipping 5 cards is much better
%    than picking one at random, skipping 10 is even better, and so on.  Of course,
%    you can't skip too many cards, and in fact, the results seem to be best for
%    somewhere around 30 to 35 cards skipped.  For problems like this, the
%    optimal value is somewhere around 1 / e, where E is the base of the natural
%    logarithm system.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    22 March 2009
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, integer DECK_SIZE, the number of cards in the deck.
%    2 <= DECK_SIZE.  Default value is 52;
%
%    Input, integer TRIAL_NUM, the number of times we will simulate this process.
%    Default value is 100.
%
%    Input, integer SKIP_NUM, the number of initial cards you plan to examine
%    but will NOT select.  If SKIP_NUM is 0, you don't look at any cards first.
%    0 <= SKIP_NUM < DECK_SIZE.  Default value is DECK_SIZE/3.
%
%    Output, real P, the estimated probability that your strategy of skipping
%    SKIP_NUM cards and then selecting the next card that is bigger, will
%    result in choosing the highest card.
%
  if ( nargin < 3 )
    trial_num = 100;
  end

  if ( nargin < 2 )
    skip_num = deck_size / 3;
  end

  if ( nargin < 1 )
    deck_size = 52;
  end
%
%  Make sure we got integers.
%
  deck_size = floor ( deck_size );
  skip_num = floor ( skip_num );
  trial_num = floor ( trial_num );
%
%  Check values.
%
  if ( deck_size < 2 )
    fprintf ( 1, '\n' );
    fprintf ( 1, 'HIGH_CARD_SIMULATION - Fatal error!\n' );
    fprintf ( 1, '  DECK_SIZE must be at least 2.\n' );
    fprintf ( 1, '  Your value was %d\n', deck_size );
    error ( 'HIGH_CARD_SIMULATION - Fatal error!' );
  end

  if ( skip_num < 0 )
    skip_num = 0;
  end

  if ( deck_size <= skip_num )
    fprintf ( 1, '\n' );
    fprintf ( 1, 'HIGH_CARD_SIMULATION - Fatal error!\n' );
    fprintf ( 1, '  SKIP_NUM must be less than DECK_SIZE.\n' );
    fprintf ( 1, '  Your values were DECK_SIZE = %d, SKIP_NUM = %d\n', deck_size, skip_num );
    error ( 'HIGH_CARD_SIMULATION - Fatal error!' );
  end

  if ( trial_num < 1 )
    fprintf ( 1, '\n' );
    fprintf ( 1, 'HIGH_CARD_SIMULATION - Fatal error!\n' );
    fprintf ( 1, '  TRIAL_NUM must be at least 1.\n' );
    fprintf ( 1, '  Your value was %d\n', trial_num );
    error ( 'HIGH_CARD_SIMULATION - Fatal error!' );
  end

  correct = 0;

  for trial = 1 : trial_num

    cards = permutation_random ( deck_size );

    if ( 1 <= skip_num )
      skip_max = max ( cards(1:skip_num) );
    else
      skip_max = -Inf;
    end

    true_max = max ( cards(1:deck_size) );
%
%  In case you don't encounter a card larger than SKIP_MAX,
%  we'll assume you pick the last card in the deck, even though
%  you know it's a loser.
%
    choice = cards(deck_size);
%
%  Turn over the remaining cards in the deck, but stop
%  immediately when you find one bigger than SKIP_MAX.
%
    for card = skip_num + 1 : deck_size
      if ( skip_max < cards(card) )
        choice = cards(card);
        break;
      end
    end
%
%  Record successful choices.
%
    if ( choice == true_max )
      correct = correct + 1;
    end

  end
%
%  Estimate the probability.
%
  p = correct / trial_num;

  return
end
function p = permutation_random ( n )

%*****************************************************************************80
%
%% PERMUTATION_RANDOM returns a random permutation.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    22 March 2009
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, integer N, the number of objects to permute.
%
%    Output, integer P(N), a permutation of the integers from 1 to N.
%
  p = ( 1 : n );

  for i = 1 : n - 1

    k = i4_random ( i, n );

    p1 = p(i);
    p(i) = p(k);
    p(k) = p1;

  end

  return
end
function value = i4_random ( lo, hi )

%*****************************************************************************80
%
%% I4_RANDOM returns a random integer between LO and HI.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    22 March 2009
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, integer LO, HI, the limits.
%
%    Output, integer VALUE, the random integer.
%
  r = ( hi + 1 - lo ) * rand ( );
  value = lo + floor ( r );

  return
end

author: johannesgerer

repo: jburkardt-m


text:
function test_null_left ( )

%*****************************************************************************80
%
%% TEST_NULL_LEFT tests left null vectors.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    15 March 2015
%
%  Author:
%
%    John Burkardt
%
  fprintf ( 1, '\n' );
  fprintf ( 1, 'TEST_NULL_LEFT\n' );
  fprintf ( 1, '  A = a test matrix of order M by N;\n' );
  fprintf ( 1, '  x = an M vector, candidate for a left null vector.\n' );
  fprintf ( 1, '\n' );
  fprintf ( 1, '  ||A|| = Frobenius norm of A.\n' );
  fprintf ( 1, '  ||x|| = L2 norm of x.\n' );
  fprintf ( 1, '  ||A''*x||/||x|| = L2 norm of A''*x over L2 norm of x.\n' );
  fprintf ( 1, '\n' );
  fprintf ( 1, '  Title                    M     N      ' );
  fprintf ( 1, '||A||            ||x||        ||A''*x||/||x||\n' );
  fprintf ( 1, '\n' );
%
%  A123
%
  title = 'A123';
  m = 3;
  n = 3;
  a = a123 ( );
  x = a123_null_left ( );
  error_l2 = r8mat_is_null_left ( m, n, a, x );
  norm_a_frobenius = r8mat_norm_fro ( m, n, a );
  norm_x_l2 = r8vec_norm_l2 ( m, x );
  fprintf ( 1, '  %-20s  %4d  %4d  %14g  %14g  %10.2g\n', ...
    title, m, n, norm_a_frobenius, norm_x_l2, error_l2 );
%
%  CHEBY_DIFF1
%
  title = 'CHEBY_DIFF1';
  m = 5;
  n = m;
  a = cheby_diff1 ( n );
  x = cheby_diff1_null_left ( m, n );
  error_l2 = r8mat_is_null_left ( m, n, a, x );
  norm_a_frobenius = r8mat_norm_fro ( m, n, a );
  norm_x_l2 = r8vec_norm_l2 ( m, x );
  fprintf ( 1, '  %-20s  %4d  %4d  %14g  %14g  %10.2g\n', ...
    title, m, n, norm_a_frobenius, norm_x_l2, error_l2 );
%
%  CREATION
%
  title = 'CREATION';
  m = 5;
  n = m;
  a = creation ( m, n );
  x = creation_null_left ( m, n );
  error_l2 = r8mat_is_null_left ( m, n, a, x );
  norm_a_frobenius = r8mat_norm_fro ( m, n, a );
  norm_x_l2 = r8vec_norm_l2 ( m, x );
  fprintf ( 1, '  %-20s  %4d  %4d  %14g  %14g  %10.2g\n', ...
    title, m, n, norm_a_frobenius, norm_x_l2, error_l2 );
%
%  DIF1
%  Only has null vectors for M odd.
%
  title = 'DIF1';
  m = 5;
  n = m;
  a = dif1 ( m, n );
  x = dif1_null_left ( m, n );
  error_l2 = r8mat_is_null_left ( m, n, a, x );
  norm_a_frobenius = r8mat_norm_fro ( m, n, a );
  norm_x_l2 = r8vec_norm_l2 ( m, x );
  fprintf ( 1, '  %-20s  %4d  %4d  %14g  %14g  %10.2g\n', ...
    title, m, n, norm_a_frobenius, norm_x_l2, error_l2 );
%
%  DIF1CYCLIC
%
  title = 'DIF1CYCLIC';
  m = 5;
  n = m;
  a = dif1cyclic ( n );
  x = dif1cyclic_null_left ( m, n );
  error_l2 = r8mat_is_null_left ( m, n, a, x );
  norm_a_frobenius = r8mat_norm_fro ( m, n, a );
  norm_x_l2 = r8vec_norm_l2 ( m, x );
  fprintf ( 1, '  %-20s  %4d  %4d  %14g  %14g  %10.2g\n', ...
    title, m, n, norm_a_frobenius, norm_x_l2, error_l2 );
%
%  DIF2CYCLIC
%
  title = 'DIF2CYCLIC';
  m = 5;
  n = m;
  a = dif2cyclic ( n );
  x = dif2cyclic_null_left ( m, n );
  error_l2 = r8mat_is_null_left ( m, n, a, x );
  norm_a_frobenius = r8mat_norm_fro ( m, n, a );
  norm_x_l2 = r8vec_norm_l2 ( m, x );
  fprintf ( 1, '  %-20s  %4d  %4d  %14g  %14g  %10.2g\n', ...
    title, m, n, norm_a_frobenius, norm_x_l2, error_l2 );
%
%  EBERLEIN
%
  title = 'EBERLEIN';
  m = 5;
  n = 5;
  r8_lo = -5.0;
  r8_hi = +5.0;
  seed = 123456789;
  [ alpha, seed ] = r8_uniform_ab ( r8_lo, r8_hi, seed );
  a = eberlein ( alpha, n );
  x = eberlein_null_left ( m, n );
  error_l2 = r8mat_is_null_left ( m, n, a, x );
  norm_a_frobenius = r8mat_norm_fro ( m, n, a );
  norm_x_l2 = r8vec_norm_l2 ( m, x );
  fprintf ( 1, '  %-20s  %4d  %4d  %14g  %14g  %10.2g\n', ...
    title, m, n, norm_a_frobenius, norm_x_l2, error_l2 );
%
%  FIBONACCI1
%
  title = 'FIBONACCI1';
  m = 5;
  n = m;
  r8_lo = -5.0;
  r8_hi = +5.0;
  seed = 123456789;
  [ f1, seed ] = r8_uniform_ab ( r8_lo, r8_hi, seed );
  [ f2, seed ] = r8_uniform_ab ( r8_lo, r8_hi, seed );
  a = fibonacci1 ( n, f1, f2 );
  x = fibonacci1_null_left ( m, n, f1, f2 );
  error_l2 = r8mat_is_null_left ( m, n, a, x );
  norm_a_frobenius = r8mat_norm_fro ( m, n, a );
  norm_x_l2 = r8vec_norm_l2 ( m, x );
  fprintf ( 1, '  %-20s  %4d  %4d  %14g  %14g  %10.2g\n', ...
    title, m, n, norm_a_frobenius, norm_x_l2, error_l2 );
%
%  LAUCHLI
%
  title = 'LAUCHLI';
  m = 6;
  n = m - 1;
  r8_lo = -5.0;
  r8_hi = +5.0;
  seed = 123456789;
 [ alpha, seed ] = r8_uniform_ab ( r8_lo, r8_hi, seed );
  a = lauchli ( alpha, m, n );
  x = lauchli_null_left ( alpha, m, n );
  error_l2 = r8mat_is_null_left ( m, n, a, x );
  norm_a_frobenius = r8mat_norm_fro ( m, n, a );
  norm_x_l2 = r8vec_norm_l2 ( m, x );
  fprintf ( 1, '  %-20s  %4d  %4d  %14g  %14g  %10.2g\n', ...
    title, m, n, norm_a_frobenius, norm_x_l2, error_l2 );
%
%  LINE_ADJ
%
  title = 'LINE_ADJ';
  m = 7;
  n = m;
  a = line_adj ( n );
  x = line_adj_null_left ( m, n );
  error_l2 = r8mat_is_null_left ( m, n, a, x );
  norm_a_frobenius = r8mat_norm_fro ( m, n, a );
  norm_x_l2 = r8vec_norm_l2 ( m, x );
  fprintf ( 1, '  %-20s  %4d  %4d  %14g  %14g  %10.2g\n', ...
    title, m, n, norm_a_frobenius, norm_x_l2, error_l2 );
%
%  MOLER2
%
  title = 'MOLER2';
  m = 5;
  n = 5;
  a = moler2 ( );
  x = moler2_null_left ( );
  error_l2 = r8mat_is_null_left ( m, n, a, x );
  norm_a_frobenius = r8mat_norm_fro ( m, n, a );
  norm_x_l2 = r8vec_norm_l2 ( m, x );
  fprintf ( 1, '  %-20s  %4d  %4d  %14g  %14g  %10.2g\n', ...
    title, m, n, norm_a_frobenius, norm_x_l2, error_l2 );
%
%  ONE
%
  title = 'ONE';
  m = 5;
  n = 5;
  a = one ( m, n );
  x = one_null_left ( m, n );
  error_l2 = r8mat_is_null_left ( m, n, a, x );
  norm_a_frobenius = r8mat_norm_fro ( m, n, a );
  norm_x_l2 = r8vec_norm_l2 ( m, x );
  fprintf ( 1, '  %-20s  %4d  %4d  %14g  %14g  %10.2g\n', ...
    title, m, n, norm_a_frobenius, norm_x_l2, error_l2 );
%
%  RING_ADJ
%  N must be a multiple of 4 for there to be a null vector.
%
  title = 'RING_ADJ';
  m = 12;
  n = 12;
  a = ring_adj ( n );
  x = ring_adj_null_left ( m, n );
  error_l2 = r8mat_is_null_left ( m, n, a, x );
  norm_a_frobenius = r8mat_norm_fro ( m, n, a );
  norm_x_l2 = r8vec_norm_l2 ( m, x );
  fprintf ( 1, '  %-20s  %4d  %4d  %14g  %14g  %10.2g\n', ...
    title, m, n, norm_a_frobenius, norm_x_l2, error_l2 );
%
%  ROSSER1
%
  title = 'ROSSER1';
  m = 8;
  n = 8;
  a = rosser1 ( );
  x = rosser1_null_left ( );
  error_l2 = r8mat_is_null_left ( m, n, a, x );
  norm_a_frobenius = r8mat_norm_fro ( m, n, a );
  norm_x_l2 = r8vec_norm_l2 ( m, x );
  fprintf ( 1, '  %-20s  %4d  %4d  %14g  %14g  %10.2g\n', ...
    title, m, n, norm_a_frobenius, norm_x_l2, error_l2 );
%
%  ZERO
%
  title = 'ZERO';
  m = 5;
  n = 5;
  a = zero ( m, n );
  x = zero_null_left ( m, n );
  error_l2 = r8mat_is_null_left ( m, n, a, x );
  norm_a_frobenius = r8mat_norm_fro ( m, n, a );
  norm_x_l2 = r8vec_norm_l2 ( m, x );
  fprintf ( 1, '  %-20s  %4d  %4d  %14g  %14g  %10.2g\n', ...
    title, m, n, norm_a_frobenius, norm_x_l2, error_l2 );

  return
end


text:
function [n, freqs] = plx_adchan_freqs(filename)
% plx_adchan_freq(filename): read the per-channel frequencies for analog channels from a .plx or .pl2 file
%
% [n, freqs] = plx_adchan_freq(filename)
%
% INPUT:
%   filename - if empty string, will use File Open dialog
%
% OUTPUT:
%   freqs - array of frequencies
%   n - number of channels

if nargin ~= 1
    error 'expected 1 input argument';
end

[ filename, isPl2 ] = internalPL2ResolveFilenamePlx( filename );
if isPl2 == 1
    pl2 = PL2GetFileIndex(filename);
    n = numel(pl2.AnalogChannels);
    if n > 0
        freqs = zeros(n,1);
        for i=1:n
            freqs(i,1) = pl2.AnalogChannels{i}.SamplesPerSecond;
        end
    end
    return;
end

[n,freqs] = mexPlex(12, filename);

end

author: buzsakilab

repo: buzcode

save_path: github-repos/MATLAB/buzsakilab-buzcode


text:
function ExampleAnalog()
% function OutputAnalogVector = ExampleAnalog()

% Prompt for the correct DLL
disp(' ');  % Blank line
DLLName = input('DLL Name: ', 's');

% Load the appropriate DLL
[nsresult] = ns_SetLibrary(DLLName);
if (nsresult ~= 0)
    disp('DLL was not found!');
    return
end

% Find out the data file from user
disp(' ');  % Blank line
filename = input('Data file: ', 's');

% Load data file and display some info about the file
% Open data file
[nsresult, hfile] = ns_OpenFile(filename);
if (nsresult ~= 0)
    disp('Data file did not open!');
    return
end
clear filename;

% Get file information
[nsresult, FileInfo] = ns_GetFileInfo(hfile);
% Gives you EntityCount, TimeStampResolution and TimeSpan
if (nsresult ~= 0)
    disp('Data file information did not load!');
    return
end

% Define some variables needed for firing rates
stepsize = 0.02; % seconds
stepsize1 = 0.1; % seconds
if FileInfo.TimeSpan > 150  % Limit the timespan shown in the graphs to 150 seconds
    totaltime = 150;
else
    totaltime = FileInfo.TimeSpan; % seconds
end
time = 0 : stepsize : totaltime;   % Initialize time axis for gaussian plot

% Build catalogue of entities
[nsresult, EntityInfo] = ns_GetEntityInfo(hfile, [1 : 1 : FileInfo.EntityCount]);

NeuralList = find([EntityInfo.EntityType] == 4);    % List of EntityIDs needed to retrieve the information and data
SegmentList = find([EntityInfo.EntityType] == 3);
AnalogList = find([EntityInfo.EntityType] == 2);
EventList = find([EntityInfo.EntityType] == 1);

% How many of a particular entity do we have
cNeural = length(NeuralList);
cSegment = length(SegmentList);
cAnalog = length(AnalogList);
cEvent = length(EventList);

clear FileInfo;

if (cNeural == 0)
    disp('No neural events available!');
end

if (cSegment == 0)
    disp('No segment entities available!');
end

if (cAnalog == 0)
    disp('No analog entities available!');
end

if (cEvent == 0)
    disp('No event entities available!');
end

% prompt to get the number of points
max_count = input('How many data points (x-axis) do you want max? ');

% Have user pick a channels or channels for further analysis
% Show the user how many channels are available
disp(' ');
disp(['There are ' num2str(cAnalog) ' analog channels.']);
disp(' ');

yes = input('Show first all analog channels? y/n ', 's');
if (yes == 'y')
    for i = 1 : 1 : length(AnalogList)
        chan = AnalogList(i);
        count = min(max_count, EntityInfo(chan).ItemCount);
        % Get the fist data points of the waveform and show it
        [nsresult, ContinuousCount, wave] = ns_GetAnalogData(hfile, AnalogList(i), 1, count);
        figure;
        plot(wave);
    end
    return;
end

disp(['First analog entity: ' num2str(AnalogList(1)) ]);
channel = input('Which data channels would you like to display? (e.g. 1 or [1 2 3])');
if (EntityInfo(channel).EntityType == 2)  % Have to check that the selected channel actually exists
else
    disp('Channel is not of type analog');
    return
end
clear cNeural cSegment;

% Throw away entity infos we don't need to save memory
EntityInfo = rmfield(EntityInfo, 'EntityType');

%
% Load the waveform data and do the analysis
%

count = min(max_count, EntityInfo(channel).ItemCount);

% Get the fist data points of the waveform and show it
[nsresult, ContinuousCount, wave] = ns_GetAnalogData(hfile, channel, 1, count);
plot(wave);
return;

author: fieldtrip

repo: fieldtrip

save_path: github-repos/MATLAB/fieldtrip-fieldtrip


text:
% ------------------------------------------------------------------------
%  Copyright (C)
%  Universitat Politecnica de Catalunya BarcelonaTech (UPC) - Spain
%  University of California Berkeley (UCB) - USA
%
%  Jordi Pont-Tuset <jordi.pont@upc.edu>
%  Pablo Arbelaez <arbelaez@berkeley.edu>
%  June 2014
% ------------------------------------------------------------------------
% This file is part of the MCG package presented in:
%    Arbelaez P, Pont-Tuset J, Barron J, Marques F, Malik J,
%    "Multiscale Combinatorial Grouping,"
%    Computer Vision and Pattern Recognition (CVPR) 2014.
% Please consider citing the paper if you use this code.
% ------------------------------------------------------------------------
function [ bdry ]  = seg2bdry_wt(seg, fmt)
if nargin<2, fmt = 'doubleSize'; end;

if ~strcmp(fmt,'imageSize') && ~strcmp(fmt,'doubleSize'),
    error('possible values for fmt are: imageSize and doubleSize');
end

[tx, ty, nch] = size(seg);

if nch ~=1,
    error('seg must be a scalar image');
end

bdry = zeros(2*tx+1, 2*ty+1);

edgels_v = abs( seg(1:end-1, :) - seg(2:end, :) );
edgels_v(end+1, :) = 0;
edgels_h = abs( seg(:, 1:end-1) - seg(:, 2:end) );
edgels_h(:, end+1) = 0;

bdry(3:2:end, 2:2:end) = edgels_v;
bdry(2:2:end, 3:2:end) = edgels_h;
bdry(3:2:end-1, 3:2:end-1)= max ( max(edgels_h(1:end-1, 1:end-1), edgels_h(2:end, 1:end-1)), max(edgels_v(1:end-1,1:end-1), edgels_v(1:end-1,2:end)) );

bdry(1, :) = bdry(2, :);
bdry(:, 1) = bdry(:, 2);
bdry(end, :) = bdry(end-1, :);
bdry(:, end) = bdry(:, end-1);

if strcmp(fmt,'imageSize'),
    bdry = bdry(3:2:end, 3:2:end);
end

author: s-gupta

repo: rcnn-depth

save_path: github-repos/MATLAB/s-gupta-rcnn-depth


text:
 function G = Gnearest(nx, ny, nb, na, mask)
%|function G = Gnearest(nx, ny, nb, na, mask)
%|
%| Generate a geometric system matrix for tomographic projection
%| based on simple nearest-neighbor "interpolation."
%| This system model is completely inadequate for reconstructing
%| real tomographic data, but is marginally useful for simple simulations.
%|
%| The image size is nx * ny.
%| The sinogram size is nb * na (n_radial_bins X n_angles).
%| Returned G is [nb * na, nx * ny] - this size is need for .wtf saving.
%|
%| Jeff Fessler

if nargin < 1, nx = 64; end
if nargin < 2, ny = nx; end
if nargin < 3, nb = nx; end
if nargin < 4, na = floor(nb * pi/2); end
if nargin < 5, mask = true(nx,ny); end

warning 'this function is obsolete.  use Gtomo2_strip() instead'

% pixel centers
[x,y] = ndgrid([0:nx-1] - (nx-1)/2, [0:ny-1] - (ny-1)/2);
x = x(mask(:));
y = y(mask(:));
np = length(x);		% sum(mask(:)) - total # of support pixels

angle = [0:na-1]'/na * pi;
tau = cos(angle) * x' + sin(angle) * y';	% [na,np] projected pixel center
ii = round(tau + (nb-1)/2);			% counting from 1 (matlab)
good = ii(:) >= 1 & ii(:) <= nb;
if any(~good), warning 'FOV too small', end

ii = ii + [0:na-1]' * nb * ones(1,np);

%np = sum(mask(:));
%nc = np;	jj = 1:np;		% compact G
nc = nx * ny;	jj = find(mask(:))';	% all-column G
jj = jj(ones(1,na),:);

G = sparse(ii(good), jj(good), ones(sum(good),1), nb*na, nc);

if 0
%	subplot(121), im(embed(sum(G)', mask))		% for compact
	subplot(121), im(reshape(sum(G), nx, ny))	% for all-column
	subplot(122), im(reshape(sum(G'), nb, na))
end

author: JeffFessler

repo: mirt

save_path: github-repos/MATLAB/JeffFessler-mirt


text:
function should_evaluate = exp_beginfun(setting, varargin)
% Begin the definition of a function with lazy evaluation behavior.
% Evaluate = exp_beginfun(AttributeSetting, Attributes...)
%
% General
% =======
% If a function uses exp_beginfun/exp_endfun, MATLAB code that uses these functions does not "call"
% them in an imperative manner, but instead builds symbolic expressions that can be stored,
% transmitted, or manipulated before/if they are being evaluated. Evaluation is generally done via
% exp_eval() (or exp_eval_optimized()). Almost all uses of exp_beginfun by functions look like this:
%
%   function [res1,res2,...] = my_function(arg1,arg2, ...)
%   if ~exp_beginfun('somesetting') return; end
%
%   ... <actual implementation> ...
%
%   exp_endfun;
%
% ... where 'somesetting' specifies the behavioral setting for the function. Generally, all
% functions that operate on data sets and modify the signal content in a way that needs to be
% replicated online have the setting 'filter', and all functions that have no equivalent during
% online processing (e.g., marker transcription) have the setting 'editing'. Functions that do
% processing that cannot be reproduced online have the setting 'offline'. See filters/flt_* and
% dataset_ops/set_* for some examples. A setting is a shortcut name for an assignment of attributes
% associated with the function declaration (attributes can also be directly specified to
% exp_beginfun).
%
%
% Filter and Editing Functions
% ============================
% Functions of type 'filter' should be written just like regular EEGLAB signal processing functions
% (i.e., they receive a data set and some options, and produce a transformed version of the data
% set), but functions that are supposed to be usable online should be "chunkable", i.e. they should
% produce the same sequence of output samples when they are called for consecutive segments (chunks)
% of the data set. Thus, the function should dump its internal state in an additional output (the
% 'state' variable) and be able to resume from this state when called on the next chunk. Chunks may
% well be just a single sample long. During online processing, state will be passed in automatically.
%
% In some cases, this strict behavior would be prohibitively inefficient - in this case, the
% function may refrain to non-adaptive behavior if the chunks are too short. In particular, in
% BCILAB, it may be assumed that each filter is first called on a complete calibration data set in a
% single call, and is then either called on a complete test data set in one call (the next 'chunk'),
% or piecemeal on small blocks of an incoming online stream. ICA, for example, adapts itself only on
% the calibration data set (i.e. the first chunk) and uses that decomposition for all remaining ones
% (i.e. during online use).
%
% Functions of type 'editing' should be written just like regular EEGLAB data set editing functions
% (e.g. pop_chanedit), BUT they must not change the data or srate fields.
%
% Rationale/Background
% ====================
% For the purpose of data processing function declarations, the exp_beginfun/exp_endfun contract
% handles the book-keeping necessary to be able to reproduce the filter (and entire filter chain)
% 1:1 online, and offers additional functionality, such as transparent caching of results (for later
% near-instantaneuous retrieval) and the annotation of all data sets with an executable expression
% which, when invoked, reproduces the data set 1:1, as well as checks to guarantee that no invalid
% "tampering" is done by the user on data sets (such as running a forward-backward filter in EEGLAB
% or manually scaling the data).
%
% In its general form, the functions exp_beginfun/exp_endfun are part of the machinery that enables
% symbolic computation in MATLAB (in the style of Mathematica), and the
% 'filter'/'editing'/'symbolic'/'annotate'/'default' argument is a short-cut which selects default
% attributes of the function being defined (which can be overridden). Some of these attributes are
% customizable plugin "slots" (plugins to be applied to each argument or all inputs or all outputs
% of the function being defined). 'Filter'/'editing' imply such hooks as default attributes, which
% in turn implement the tamper-proofing checks and some post-processing on the emitted data sets.
% The post-processing adds fields (subfields of the .tracking field) to the (otherwise arbitrary)
% data set structs that are produced by the function being defined, most importantly the
% .tracking.expression, which is an executable (recursive) expression representation which mirrors
% each data set in symbolic form, and .tracking.online_expression, in the case of filters, which
% defines the arguments to be used for the respective filter during online processing. These
% operations are usually completely transparent to the user (and locked away in the .tracking
% field), and make it exceptionally simple to author components for the toolbox.
%
% The symbolic-computation functionality implemented by the expressions/exp_* functions is largely
% non-overlapping with what the MATLAB Symbolic toolbox offers: the Symbolic toolbox exclusively
% implements mathematics (e.g. derivatives, integrals, simplifications), while the expression system
% of BCILAB exclusively implements expression manipulation, such as substitution, pattern matching,
% hashing, memoization, and reflection. These facilities allow to implement very high-level behavior
% in functions such as bci_train, such as generic programming (e.g. utl_crossval, utl_searchmodel,
% utl_gridsearch being non-intrusively adapted to EEG sets, STUDY sets, etc.), optimizations (e.g.
% in-memory caching and on-disk caching of intermediate results, computational shortcuts, and more
% in the future), with relatively little coding effort.
%
% While, by default, a function is invoked online with the same parameters that were used for
% offline processing (aside from a 'state' parameter for stateful functions), it is is possible to
% specify that certain parameters should have different values during online processing. This is
% enabled by the attribute 'append_online' (list of arguments to be appended in the online case),
% which can be specified either to exp_beginfun or exp_endfun. Alternatively, the entire expression
% can be subsituted by setting the 'set_online' attribute (expression to use during online
% processing).
%
% In:
%  AttributeSetting: baseline attribute assignment for the function; must be specified. Below are
%                    typical scenarios:
%                    * defining a function that should not evaluate immediately, but only after
%                      evaluation by an exp_eval() clause:
%                      'default' : delayed execution, input fingerprint checks
%
%                    * defining a function that processes data sets and should have
%                      delayed-evaluation semantics:
%                      'filter'  : functions that do non-trivial processing of data sets, and which
%                                  should support online replication of some offline processing;
%                                  attributes: delayed (1), add_impure (1), set_online
%                                  ('reproduce'), also performs EEGLAB-aware input data set checks &
%                                  transforms note: when an offline processing is replicated online,
%                                  the function is called successively on short EEGLAB data sets,
%                                  with inputs identical to the input parameters as they were
%                                  snapshot after execution of the function's body (by exp_endfun)
%                      'editing' : functions that do processing of data sets which has no equivalent
%                                  in online processing, i.e. which pass the input data through
%                                  unmodified; attributes: delayed (1), add_impure (1), set_online
%                                  ('passthrough'), also performs EEGLAB-aware input data set checks
%                      'offline' : functions that do processing of data sets which cannot be
%                                  reproduced online; these generate an error, when used within
%                                  online BCIss; attributes: delayed (1), add_impure (1), set_online
%                                  ('inapplicable'), also performs EEGLAB-aware input data set
%                                  checks
%
%                    * definining a function that imports or constructs data sets from scratch:
%                      'import'  : functions that construct a data set from scratch based on
%                                  arguments that are not data sets
%
%                    * defining a function that performs purely symbolic manipulations of
%                      expressions:
%                      'symbolic': all inputs are held unevaluated, associated (impure) values are
%                                  ignored (thus no fingerprinting and no memoization); attributes:
%                                  delayed (1), hold ('all'), fingerprints unchecked & passed
%                                  through, memoize (0)
%                                  note: symbolic functions cannot be used in online code.
%
%                    * defining a function that performs immediate (non-lazy) processing of data
%                      structures, and which shall have their results annotated with an executable
%                      (repeatable) description of the processing performed:
%                      'annotate': delayed (0), add_impure (1)
%
%
%   Attributes... : optional name-value pairs to refine attributes of the function being defined:
%
%                     'delayed'   : whether the function has delayed-evaluation semantics, i.e. it
%                                   doesn't evaluate unless passed to exp_eval (default: 1)
%
%                     'hold'      : []/'all'/'first'/'rest'; whether to leave all/the first/the
%                                   remaining arguments or no arguments of the function unevaluated
%                                   before executing the function body (useful for higher-order
%                                   functions), (default: [])
%
%                     'symbolic'  : whether the function being defined should be evaluated even if it
%                                   receives unevaluated expressions as arguments (ignoring those
%                                   arguments that fall under hold) - this is for symbol-manipulating
%                                   functions that can process unevaluated expressions (0/1, default: 0)
%
%                     'add_impure': integrate the original expression into the function value (when
%                                   evaluated), creating an 'impure expression', 0/1
%
%                     'set_online': Relevant only when the settings attribute is 'filter', 'editing'
%                                   or 'offline', may also be specified in exp_endfun. Introduces an
%                                   'online expression' into the function value (when evaluated),
%                                   which is a record of how the resulting data (usually a signal)
%                                   should be calculated online. Usually, this is identical to how
%                                   the filter is applied offline.
%                                   * 'passthrough': the function is skipped in the online
%                                                    processing
%                                   * 'reproduce': the function is called online with the same arguments
%                                                  as in offline processing
%                                   * 'inapplicable': the function gives an error when used in
%                                                     online processing
%                                   * {arg1, arg2, arg3, ...}: the given arguments will be used as
%                                                              function arguments during online use
%                                   * expression: the given expression will be used during online
%                                                 reproduction of the filter
%
%                     'append_online' : {arg1, arg2, arg3, ...}: the given arguments will be
%                                       appended to the function arguments during online use
%                                       (default: [])
%
%                     'argsteps', 'presteps', 'poststeps' : advanced attributes, see also "Advanced
%                                                           Functionality" section.
%
%                     'delayed_online' : whether this function has delayed-evaluation semantics
%                                        during online processing (default: false)
%
% Out:
%   Evaluate     : whether the body of the function being defined shall be evaluated or skipped
%
%
% Advanced functionality:
%   The allowed settings are customizable/extensible via the function exp_settings, which maps a
%   setting string onto a cell array of default attributes (name-value pairs).
%
%   During evaluation, each input argument is sent though a customizable sequence of functions,
%   which can be specified via the 'argsteps' attribute (as a cell array). Further, the attriutes
%   'presteps' and 'poststeps' are cell arrays of functions that should be invoked before and after
%   the function body runs, respectively; they operate on a context structure, which is defined
%   further below. These built-in attributes are the primary mechanism by which custom exp_settings
%   add functionality to exp_beginfun/exp_endfun (e.g., annotating a dataset with an expression or
%   checking fingerprints).
%
%   The function's behavior can be controlled via dynamically scoped variables, as set by exp_block,
%   exp_set and exp_setdelayed, or correspondent attributes in exp_beginfun.
%     These variables can be assigned values, such as 0/1, but may be arbitrary expressions
%      * fingerprint_check: toggle per-argument fingerprint checking (default: true) - may use
%        @expression to refer to the argument expression, to implement expression-dependent checking
%      * fingerprint_create: toggle fingerprint creation (default: true) - may use @expression to
%        refer to the unevaluated complete expression, to implement expression-dependent checking
%      * memoize: selectively enable memoization; this is a cell array of the form
%        {location,expression,location,expression, ...}. location can be either 'disk' or 'memory',
%        expression should evaluate to true/false, may use @expression to refer to the unevaluated
%        complete expression, to implement expression-dependent memoization.
%
% See also:
%   exp_endfun, exp_eval, exp_settings
%
% Examples:
%   % a typical function declaration using the expression system
%   function result = myfunc(myarg,myotherarg,yetanotherarg)
%   if ~exp_beginfun('somesetting'), return; end
%   ...
%   exp_endfun;
%
%
%   % define a filter function using the expression system
%   % (note that, for filters to work with all BCILAB facilities, the varargin arguments should also
%   %  be processed via arg_define)
%   function signal = flt_myfilter(varargin)
%   if ~exp_beginfun('filter'), return; end
%   ...
%   exp_endfun;
%
%
%   % as before, but define that the filter, when applied during online processing, should be passed
%   % a particular list of additional arguments
%   function signal = flt_myfilter(varargin)
%   if ~exp_beginfun('filter'), return; end
%   ...
%   exp_endfun('add_online',{'myinternalarg',10, 'my_reserved_argument','test'})
%
%
%   % like before, but this time replace the entire set of arguments; the assumption here is that
%   % the variable called signal contains the input signal, and that it can be passed to flt_myfilter
%   % by name
%   function signal = flt_myfilter(varargin)
%   if ~exp_beginfun('filter'), return; end
%   ...
%   exp_endfun('set_online',{'signal',signal, 'myparam',10,'myotherparam',200})
%
%
%   % like before, but this time define that the filter, when applied online, should go through an
%   % entirely different function (here assumed to be defined using exp_beginfun itself) with specific
%   % arguments
%   function signal = flt_myfilter(varargin)
%   if ~exp_beginfun('filter'), return; end
%   ...
%   exp_endfun('set_online',flt_myfilter_online('signal',signal, 'specialparam1',10))
%
%
%   % define a dataset editing function
%   function signal = set_myoperation(varargin)
%   if ~exp_beginfun('editing'), return; end
%   ...
%   exp_endfun;
%
%
%   % define a function that operates on expressions symbolically
%   function result = my_operation(varargin)
%   if ~exp_beginfun('symbolic'), return; end
%   ...
%   exp_endfun;
%
%
%   % define a signal processing function that can only be applied in offline processing
%   function signal = set_myoperation(varargin)
%   if ~exp_beginfun('offline'), return; end
%   ...
%   exp_endfun;
%
%
%                                Christian Kothe, Swartz Center for Computational Neuroscience, UCSD
%                                2010-04-14

% --- obtain & interpret calling context ---

% obtain execution context
ctx.stack = dbstack('-completenames');

% debug print
dp(['in ' ctx.stack(2).name]);

% early check of whether the expression system is disabled
if  ~strcmp(setting,'symbolic') && hlp_resolve('disable_expressions',false,ctx) && ~isequal(varargin,{'delayed_online',true})
    % and exit exp_beginfun, if so
    assignin('caller','exp_internal_context',struct('enabled',0, 'exec_ctx',ctx));
    % the function should evaluate
    should_evaluate = 1;
    return;
end

% if we're being called by exp_eval, the function should be immediately evaluated - otherwise it may
% evaluate in a delayed fashion (i.e., only the expression itself is returned)
should_evaluate = any(strncmp('exp_eval',{ctx.stack(3:min(4,length(ctx.stack))).name},8));

% identify the calling function (i.e., the one using exp_beginfun)
callname = ctx.stack(2).name;
callfile = ctx.stack(2).file;
callfunc = str2func(callname);
[codehash,inargs,outargs] = hlp_fileinfo(callname,callfile);

% --- read the caller's function arguments (and unfold varargin, if any) ---
try
    parts = {};
    for i=1:length(inargs)
        val = evalin('caller',inargs{i});
        if ~strcmp(inargs{i},'varargin')
            val = {val}; end
        parts = [parts val]; %#ok<AGROW>
    end
catch, end

% --- parse the attribute settings ---
if isempty(varargin) && strcmp(setting,'symbolic')
    % fast path: symbolic operations with no settings
    if ~should_evaluate
        % -> delayed evaluation: compose & return a new expression from transformed
        % arguments
        assignin('caller',outargs{1},struct('head',{callfunc},'parts',{parts}));
        % also assign [] to all other output arguments
        for k=2:length(outargs)
            assignin('caller',outargs{k},[]); end
    else
        % -> immediate evaluation: opt out of exp_endfun() processing (nothing to do for
        % symbolic functions)
        assignin('caller','exp_internal_context',struct('enabled',{0}, 'exec_ctx',{ctx}));
    end
    return;
else
    % regular path: every other case

    % define general attributes
    attribs = {'delayed',1,'add_impure',1,'hold',[],'symbolic',0,'set_online',[],'append_online',[],'delayed_online',false,'argsteps',[],'presteps',[],'poststeps',[],'fingerprint_check',[],'fingerprint_create',[],'memoize',[]};

    % set defaults according to the chosen setting
    switch setting
        case 'symbolic'
            defaults = {'add_impure',0, 'symbolic',1, 'hold','all', 'fingerprint_check',0, 'fingerprint_create','passthrough', 'memoize',0};
        case 'default'
            defaults = {'add_impure',0, 'argsteps',{@utl_check_fingerprint}};
        case 'annotate'
            defaults = {'delayed',0};
        otherwise
            defaults = exp_settings(setting);
    end

    % concatenate general attribs, per-setting defaults and user overrides into one list of
    % name-value pairs (NVPs)
    if ~iscellstr(varargin(1:2:end))
        error('The given attributes must be name-value pairs.'); end
    nvps = [attribs defaults varargin];
    % retain only the last assignment for each name
    [s,indices] = sort(nvps(1:2:end));
    indices(strcmp(s((1:end-1)'),s((2:end)'))) = [];
    % and turn them into an options struct
    opts = cell2struct(nvps(2*indices),nvps(2*indices-1),2);
end

% if delayed is false, it enforces immediate evaluation
if ~opts.delayed
    should_evaluate = true; end

if should_evaluate

    % apply per-argument transformation steps (and pass the full options)
    expression = struct('head',{callfunc},'parts',{parts},'codehash',{codehash});
    for step = opts.argsteps
        for p=1:length(parts)
            parts{p} = step{1}(parts{p},opts,ctx,expression); end
    end

    % --- check if the result of evaluating this expression is already in the cache ---

    % purify the expression
    expression = utl_purify_expression(struct('head',{callfunc},'parts',{parts},'codehash',{codehash}));

    % quick check if memoization is disabled
    if ~isequal(opts.memoize,0)
        % do a cache lookup
        [action,result] = utl_memoize_lookup(expression,opts.memoize,ctx);
        % act according to outcome
        if strcmp(action,'return')
            % assign the retrieved data to the output argument(s)
            for k=1:length(result)
                assignin('caller',outargs{k},result{k}); end
            % and make sure that the function does not evaluate further
            should_evaluate = 0;
            return;
        elseif strcmp(action,'memoize')
            % remember to memoize the expression later, under the given id(s)
            memoize_id = result;
        else
            % do not memoize the expression
            memoize_id = [];
        end
    else
        memoize_id = [];
    end

    % --- not in cache: evaluate arguments  ---

    % evaluate the function's arguments (recursively), if they are not held
    if ~isempty(parts) && ~any(strcmp(opts.hold,{'all','first'}));
        parts{1} = exp_eval(parts{1},inf);
        if ~opts.symbolic && ~is_evaluated(parts{1})
            % one of the inputs could not be evaluated (e.g. using an undefined symbol): return the
            % expression
            assignin('caller',outargs{1},expression);
            % also assign [] to all other output arguments
            for k=2:length(outargs)
                assignin('caller',outargs{k},[]); end
            should_evaluate = false;
            return;
        end
    end
    if ~any(strcmp(opts.hold,{'all','rest'}))
        for i=2:length(parts)
            parts{i} = exp_eval(parts{i},inf); %#ok<AGROW>
            if ~opts.symbolic && ~is_evaluated(parts{i})
                % one of the inputs could not be evaluated (e.g. using an undefined symbol): return
                % the expression
                assignin('caller',outargs{1},expression);
                % also assign [] to all other output arguments
                for k=2:length(outargs)
                    assignin('caller',outargs{k},[]); end
                should_evaluate = false;
                return;
            end
        end
    end

    % --- prepare function execution ---

    % begin measuring the evaluation time of the expression ...
    eval_time = tic;
    % ... and measure the size of the inputs (both used to determine the best cache location)
    stats = whos('parts');
    input_size = stats.bytes;

    % create a context struct to keep track of the inputs, outputs, and attributes of the function
    % throughout its evaluation
    context = struct( ...
        'inargs',{inargs}, 'outargs',{outargs}, ...                     % input/output argument names
        'ws_input_pre',{struct()}, 'ws_output_post',{struct()}, ...     % input/output function workspace before/after evaluation of the function's body
        'expression_preeval',{expression}, 'expression_posteval', struct('head',{callfunc},'parts',{parts}), ... % expression arguments (original and evaluated)
        'opts',{opts}, 'eval_time',{eval_time}, 'input_size',{input_size}, 'memoize_id',{memoize_id}, 'exec_ctx',{ctx}, 'enabled',{1});    % misc fields

    % re-derive the function's input & output workspace from the now evaluated parts (folding
    % varargs again)
    for i=1:length(parts)
        if strcmp(inargs{i},'varargin')
            context.ws_input_pre.(inargs{i}) = parts(i:end);
            break;
        else
            context.ws_input_pre.(inargs{i}) = parts{i};
        end
    end

    % execute the prestep(s) on the function workspace
    for step = opts.presteps'
        context = step{1}(context); end

    % store the context object...
    assignin('caller','exp_internal_context',context);

    % assign all input for the function body's execution
    for fn=fieldnames(context.ws_input_pre)'
        assignin('caller',fn{1},context.ws_input_pre.(fn{1})); end

    % (from here on, the function body will execute)
else

    % --- delayed evaluation: return the expression ---

    expression = struct('head',{callfunc},'parts',{parts},'codehash',{codehash});
    assignin('caller',outargs{1},expression);
    % also assign [] to all other output arguments
    for k=2:length(outargs)
        assignin('caller',outargs{k},[]); end

    % (the function body will return immediately, as the should_evaluate output is false)
end

repo: REST