text:
function vis_trajgui(trajStruct,arg)

% VIS_TRAJGUI subfuntion for SOM_TRAJECTORY
%
% This function is the actual GUI called by SOM_TRAJECTORY
% function.
%
% See also SOM_TRAJECTORY.

% Contributed code to SOM Toolbox 2.0, February 11th, 2000 by Juha Parhankangas
% Copyright (c) by Juha Parhankangas.
% http://www.cis.hut.fi/projects/somtoolbox/

% Version 2.0beta juha 180699

if nargin == 1

  sM_h=trajStruct.figure;

  if size(trajStruct.bmus,1) ~= 1 && size(trajStruct.bmus,2) ~= 1
    fuzzy_traj(trajStruct,[]);
    return;
  end

  if size(trajStruct.bmus,1) == 1 || size(trajStruct.bmus,2) == 1

    udata.bmus = trajStruct.bmus;
    udata.a_h=[findobj(get(sM_h,'Children'),'Tag','Uplane');...
	       findobj(get(sM_h,'Children'),'Tag','Cplane')];
    udata.sM_h=trajStruct.figure;
    udata.traj=[];
    data1 = trajStruct.primary_data;
    if ~isempty(trajStruct.primary_names)
      names=trajStruct.primary_names;
    else
      for i=1:size(data1,2)
	names{i,1}=sprintf('Var%d',i);
      end
    end

    udata.lattice=trajStruct.lattice;
    form = 0.7*vis_patch(udata.lattice);
    udata.msize = trajStruct.msize;

    %%%%%%%%%%%%%%%%%%%%%%%%
    %
    % forming a patch object, which is placed above every component plane
    %

    l = size(form,1);

    nx = repmat(form(:,1),1,prod(udata.msize));
    ny = repmat(form(:,2),1,prod(udata.msize));

    x=reshape(repmat(1:udata.msize(2),l*udata.msize(1),1),l,prod(udata.msize));
    y=repmat(repmat(1:udata.msize(1),l,1),1,udata.msize(2));

    if strcmp(udata.lattice,'hexa')
      t = find(~rem(y(1,:),2));
      x(:,t)=x(:,t)+.5;
    end
    x=x+nx;
    y=y+ny;

    colors=reshape(ones(prod(udata.msize),1)*[NaN NaN NaN],...
		   [1 prod(udata.msize) 3]);

    set(0,'CurrentFigure',udata.sM_h);

    %%%%%%%%%%%%%%%%%%%%%%
    %
    % drawing patch
    %
    % caxis -commands keep the colormap of the original patch unchanged.
    %

    for i=1:length(udata.a_h)
      udata.real_patch(i)=get(udata.a_h(i),'Children');
      set(udata.real_patch(i),'ButtonDownFcn',...
			'vis_trajgui([],''click'')');
      subplot(udata.a_h(i));
      v=caxis;
      udata.tmp_patch(i)=patch(x,y,colors,'EdgeColor','none',...
			       'ButtonDownFcn',...
			       'vis_trajgui([],''click'')',...
			       'Tag','TmpPatch');
      caxis(v);
    end

    %%%%%%%%%%%%%%%%%%%%

    udata.length_of_traj=length(trajStruct.size);
    udata.size=trajStruct.size;
    udata.color=trajStruct.color;
    udata.poly.x=[];
    udata.poly.y=[];
    udata.poly.h=[];
    udata.new_marks=[];
    udata.all_marks=[];
    udata.d_mark2=[];
    udata.fig1 = figure;
    set(udata.fig1,'KeyPressFcn','vis_trajgui([],''key'')',...
		   'Name','Primary Data');

    %%%%%%%%%%%%%%%%%%%%
    %
    % making the 'Tools' -menu
    %

    udata.m_i=uimenu(udata.fig1,'Label','Trajectoy T&ools');
    udata.m_i(2)=uimenu(udata.m_i,'Label','&Remove Trajectory',...
			'Callback',...
			'vis_trajgui([],''remove_traj'')');
    udata.m_i(3)=uimenu(udata.m_i(1),'Label','&Dye Nodes',...
			'Callback',...
			'vis_trajgui([],''dye_gui'')');
    udata.m_i(4)=uimenu(udata.m_i(1),'Label','&Clear Markers',...
			'Callback',...
			'vis_trajgui([],''clear'')');
    udata.m_i(5)=uimenu(udata.m_i(1),'Label','&Save',...
			'Callback',...
			'vis_trajgui([],''save'')');
    udata.m_i(6)=uimenu(udata.m_i(1),'Label','&Load',...
			'Callback',...
			'vis_trajgui([],''load'')');

    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    %
    % drawing data components to the figure ....
    %
    %

    if nargin < 5 || isempty(comps) || (ischar(comps) && strcmp(comps,'all'))
      comps = 1:size(data1,2);
    end

    x=1:size(data1,1);

    for i=1:length(comps)
      subplot(length(comps),1,i);
      udata.h(i)=gca;
      udata.d_mark(i).h=[];

      udata.d(i)=plot(x,data1(:,comps(i)),...
		      'ButtonDownFcn',...
		      'vis_trajgui([],''line_down'')');
      set(gca,'XLim',[1 size(data1,1)],...
	      'XTick',[],...
	      'ButtonDownFcn','vis_trajgui([],''line_down'')'); %,...
	      %'YLim',[min(data1(:,comps(i))) max(data1(:,comps(i)))]);

      ylabel(names{comps(i)});
      hold on;
      ymin=get(udata.h(i),'YLim');
      pos=mean(get(udata.h(i),'XLim'));
      udata.l(i) = line([pos pos],[ymin(1) ymin(2)],...
			'Color','red',...
			'ButtonDownFcn',...
			'vis_trajgui([],''down'')');
    end
    udata.text1=[];

    udata.fig2=[];
    udata.h2=[];
    udata.l2=[];
    udata.text2=[];

    %%%%%%%%%%%%%%%%%%%%%%%%%%%
    %
    % ... and to the figure 2.
    %

    if ~isempty(trajStruct.secondary_data)
      data2=trajStruct.secondary_data;
      if isempty(trajStruct.secondary_names)
	for i=1:size(data1,2)
	  names2{i,1}=sprintf('Var%d',i);
	end
      else
	names2=trajStruct.secondary_names;
      end

      udata.fig2 = figure;
      set(udata.fig2,'Name','Secondary Data');
      set(udata.fig2,'KeyPressFcn',...
		     'vis_trajgui([],''key'')');
      for i=1:size(data2,2)
	subplot(size(data2,2),1,i);
	udata.h2(i) = gca;
	udata.d_mark2(i).h=[];
	udata.d2(i) = plot(x,data2(:,i),...
			   'ButtonDownFcn',...
			   'vis_trajgui([],''line_down'')');
	set(gca,'XLim',[1 size(data1,1)],'XTick',[],'ButtonDownFcn',...
		'vis_trajgui([],[],[],[],[],[],''line_down'')');
	ylabel(names2{i});
	hold on;
	ymin = get(udata.h2(i),'YLim');
	pos = mean(get(udata.h2(i),'XLim'));
	udata.l2(i) = line([pos pos],ymin,'Color','red',...
			   'ButtonDownFcn','vis_trajgui([],''down'')');
      end
    end

    %%%%%%%%%%%%%%%%%%%%%%%%%%

    set(udata.fig1,'UserData',udata);
    if ~isempty(udata.fig2);
      tmp.fig1=udata.fig1;
      set(udata.fig2,'UserData',tmp);
    end
    tmp=get(udata.sM_h,'UserData');
    tmp.fig1=udata.fig1;
    set(udata.sM_h,'UserData',tmp);
    set_numbers(round(pos));
    return;
  end

end

%%%%%%%%%%%%%%%%%
%
% if figures have been drawn, the only function calls that may exist
% are the ones that change the state of the application.
%
udata=get(gcf,'UserData');
udata=get(udata.fig1,'UserData');

switch arg
 case 'fuzzy'
  fuzzy_traj(tS,[]);
  return;
 case 'move_fuzzy'
  fuzzy_traj([],'move');
  return;
 case 'remove_traj'
  remove_traj;
  return;
 case 'line_down'
  line_bdf('down');
  return;
 case 'line_drag'
  line_bdf('drag');
  return;
 case 'line_up'
  line_bdf('up');
  return;
 case 'dye_gui';
  color_gui(udata.fig1);
  return;
 case {'dye','cyan','magenta','yellow','red','green','blue','white','grey'}
  dye_nodes(arg);
  return;
 case 'clear'
  clear_markers;
  return;
 case 'key'
  key_bdf;
  return;
 case 'click'
  click;
  return;
 case 'save'
  save_data;
  return;
 case 'load'
  load_data;
  return;
end

%%%%%%%%%%
%
% lines in the data figure(s) are dragged ...
%

udata=get(gcf,'UserData');
udata=get(udata.fig1,'UserData');

lims=get(gca,'XLim');
x = getfield(get(gca,'CurrentPoint'),{1});  % the location of the line

if x < lims(1)
  x=lims(1);
elseif x > lims(2)
  x=lims(2);
end

old = gcf;

switch arg
 case 'down',...

  % mouse button is pressed down above the line

  set(gcf,'WindowButtonMotionFcn','vis_trajgui([],''drag'')');
  set(gcf,'WindowButtonUpFcn','vis_trajgui([],''up'')');
  set(udata.l,'EraseMode','xor');
  delete(udata.text1);
  if ~isempty(udata.l2);
    set(udata.l2,'EraseMode','xor');
    delete(udata.text2);
  end
  set(gcf,'Pointer','crosshair');

 case 'drag'
  % change the location of the lines

  set(0,'CurrentFigure',udata.fig1);
  set(udata.l,'XData',[x x]);
  if ~isempty(udata.fig2)
    set(0,'CurrentFigure',udata.fig2);
    set(udata.l2,'XData',[x x]);
  end
  draw_traj(round(x));
  set(0,'CurrentFigure',old);
 case 'up'

  % draw trajectory and set figure to the normal state.

  set(udata.l,'EraseMode','normal');
  set(gcf,'Pointer','arrow','WindowButtonMotionFcn','',...
	  'WindowButtonUpFcn','');
  draw_traj(round(x));
  set_numbers(round(x));
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function draw_traj(point)

udata=get(gcf,'UserData');
udata=get(udata.fig1,'UserData');
color=udata.color;
eMode='normal';
if ischar(udata.color) && strcmp(udata.color,'xor')
  eMode='xor';
  color='black';
end
ind=udata.bmus(point);
[i j] = ind2sub(udata.msize,ind);
if ~mod(i,2)
  j=j+0.5;
end
old = gcf;
set(0,'CurrentFigure',udata.sM_h);
hold on;
if isempty(udata.traj) || length(udata.traj.h) ~= length(udata.a_h)

  % trajectory does not exist

  for i=1:length(udata.a_h)
    subplot(udata.a_h(i));
    hold on;
    new.h = plot(j,i,'Color',color,'EraseMode',eMode,'LineStyle','none');
    udata.traj.h(i)=new;
    udata.traj.j=j;
    udata.traj.i=i;
  end
else
  if length(udata.traj.j) == udata.length_of_traj
    % if the length of trajectory == ...,
    udata.traj.j(1) = [];         % the first (the oldest) coordinate pair
    udata.traj.i(1) = [];         % is removed.
  end
  udata.traj.j=[udata.traj.j;j]; % the new point is added to the
  udata.traj.i=[udata.traj.i;i]; % end of coordinate vectors (i and j)
  for i=1:length(udata.a_h)
    subplot(udata.a_h(i));            % remove the existing trajectory
    delete(udata.traj.h(i).h);          % and plot the new one.
    for j=1:length(udata.traj.j)
      udata.traj.h(i).h(j)=plot(udata.traj.j(j),udata.traj.i(j),...
				'Color',color,...
				'EraseMode',eMode,'Marker','o','LineWidth',2,...
				'MarkerSize',udata.size(udata.length_of_traj-j+1),...
				'LineStyle','none');
    end
  end
end
set(0,'CurrentFigure',udata.fig1);
set(udata.fig1,'UserData',udata);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function set_numbers(x)

% This function writes the numbers beside of the pointer lines

udata=get(gcf,'UserData');
udata=get(udata.fig1,'UserData');

xlim = get(gca,'XLim');
ylim = get(gca,'YLim');
p = ylim(1) + 0.9*(ylim(2)-ylim(1));

old = gcf;
set(0,'CurrentFigure',udata.fig1);

for i=1:length(udata.h)
  subplot(udata.h(i));

  % check if the text is placed to the left side of the line...

  if abs(x-xlim(1)) > (abs(x-xlim(2)))
    udata.text1(i)=text(x-1,p,sprintf('%d ->',x),...
                        'VerticalAlignment','top',...
                        'HorizontalAlignment','right',...
                        'FontWeight','demi');
  else

    %  or to the right side.

    udata.text1(i)=text(x+1,p,sprintf('<- %d',x),...
			'VerticalAlignment','top',...
			'FontWeight','demi');
  end
end

if ~isempty(udata.fig2)
  set(0,'CurrentFigure',udata.fig2);

  for i=1:length(udata.h2)
    subplot(udata.h2(i));

    if abs(x-xlim(1)) > (abs(x-xlim(2)))
      udata.text2(i)=text(x-1,p,sprintf('%d ->',x),...
			  'VerticalAlignment','top',...
			  'HorizontalAlignment','right',...
			  'FontWeight','demi');
    else
      udata.text2(i)=text(x+1,p,sprintf('<- %d',x),...
			  'VerticalAlignment','top',...
			  'FontWeight','demi');
    end
  end
end

set(0,'CurrentFigure',old);
set(udata.fig1,'UserData',udata);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function remove_traj()

% delete trajectory -object from every component plane.

udata=get(gcf,'UserData');
udata=get(udata.fig1,'UserData');

if isempty(udata.traj)
  return;
end

for i=1:length(udata.traj.h)
  delete(udata.traj.h(i).h);
end

udata.traj=[];
set(udata.fig1,'UserData',udata);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function line_bdf(arg)

% this function takes care of action when region is selected in the
% data figure.

udata=get(gcf,'UserData');
udata=get(udata.fig1,'UserData');
xlim=get(gca,'XLim');

if ~(any(strcmp('THIS',fieldnames(udata))))
  p = getfield(get(gca,'CurrentPoint'),{1});
else
  p = getfield(get(udata.THIS,'CurrentPoint'),{1});
end

if p < xlim(1)
  p = xlim(1);
elseif p > xlim(2)
  p = xlim(2);
end

switch arg
 case 'down'

  % the mouse button is pressed down, the pointer lines are drawn.
  % and the state of the figure is changed.

  udata.THIS=gca;
  set(gcf,'WindowButtonMotionFcn',...
	  'vis_trajgui([],''line_drag'')',...
	  'WindowButtonUpFcn','vis_trajgui([],''line_up'')');
  udata.start_p=p;

  old = gcf;
  set(0,'CurrentFigure',udata.fig1);

  for i=1:length(udata.h)
    subplot(udata.h(i));
    udata.t_line.h(i)=line([p p],get(gca,'YLim'),'Color','red');
    udata.t_line.h2(i)=line([p p],get(gca,'YLim'),'Color','red',...
			    'EraseMode','xor');
  end
  if ~isempty(udata.h2)
    set(0,'CurrentFigure',udata.fig2);
    for i=1:length(udata.h2)
      subplot(udata.h2(i));
      udata.t_line2.h(i)=line([p p],get(gca,'YLim'),'Color','red');
      udata.t_line2.h2(i)=line([p p],get(gca,'YLim'),'Color','red',...
			       'EraseMode','xor');
    end
  end

 case 'drag'

  % change the position of the pointer lines

  old = gcf;
  set(0,'CurrentFigure',udata.fig1);
  set(udata.t_line.h2,'XData',[p p]);
  if ~isempty(udata.fig2)
    set(0,'CurrentFigure',udata.fig2);
    set(udata.t_line2.h2,'XData',[p p]);
  end
  set(0,'CurrentFigure',old);
 case 'up'

  % sort the 'points' -vector and draw the markers to the data and nodes

  points=sort([round(udata.start_p) round(p)]);
  draw_markers(points(1):points(2));
  udata=get(udata.fig1,'UserData');
  udata.new_marks=unique([udata.new_marks ;(points(1):points(2))']);
  delete([udata.t_line.h2 udata.t_line.h]);
  if ~isempty(udata.fig2)
    delete([udata.t_line2.h2 udata.t_line2.h]);
  end
  set(get(udata.THIS,'Parent'),'WindowButtonMotionFcn','',...
		    'WindowButtonUpFcn','');
  udata=rmfield(udata,'THIS');
end

set(udata.fig1,'UserData',udata);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function draw_markers(x)

plot2data(x);
plot2plane(x);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function plot2data(x)

% plot black markers to the data figure(s)

udata=get(gcf,'UserData');
udata=get(udata.fig1,'UserData');

old = gcf;

set(0,'CurrentFigure',udata.fig1);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% if there already exist points in the positions that are members
% of the set x, then the old points are removed from data figures...

for i=1:length(udata.d_mark(1).h)
  tmp1 = get(udata.d_mark(1).h(i),'XData');
  tmp2 = setdiff(tmp1,x);
  if length(tmp1) ~= length(tmp2)
    inds=[];
    for j=1:length(tmp2);
      inds=[inds find(tmp2(j)==tmp1)];
    end
    for j=1:length(udata.d_mark)
      ydata=getfield(get(udata.d_mark(j).h(i),'YData'),{inds});
      set(udata.d_mark(j).h(i),'XData',tmp2,'YData',ydata);
    end
    if ~isempty(udata.fig2)
      for j=1:length(udata.d_mark2)
        ydata=getfield(get(udata.d_mark2(j).h(i),'YData'),{inds});
        set(udata.d_mark2(j).h(i),'XData',tmp2,'YData',ydata);
      end
    end
  end
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% ... and the new ones are plotted.

for i=1:length(udata.h)
  subplot(udata.h(i));
  h=plot(x,getfield(get(udata.d(i),'YData'),{x}),'oblack',...
	 'ButtonDownFcn',...
	 'vis_trajgui([],''line_down'')');
  udata.d_mark(i).h=[udata.d_mark(i).h;h];
end

if ~isempty(udata.h2)
  set(0,'CurrentFigure',udata.fig2);

  for i=1:length(udata.h2)
    subplot(udata.h2(i));
    h=plot(x,getfield(get(udata.d2(i),'YData'),{x}),'oblack',...
	   'ButtonDownFcn',...
	   'vis_trajgui([],''line_down'')');
    udata.d_mark2(i).h=[udata.d_mark2(i).h;h];
  end
end

set(0,'CurrentFigure',old);
set(udata.fig1,'UserData',udata);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function plot2plane(x)

% sets markers to the component planes.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% actually new markers are never plotted, but the color of the patch
% lying above the original component plane patch is changed black in
% the right positions.

udata=get(gcf,'UserData');
udata=get(udata.fig1,'UserData');

udata.new_marks=unique([udata.new_marks ;x']);

for i=1:length(udata.a_h)
  col=get(udata.tmp_patch(i),'FaceVertexCData');
  if length(size(col)) == 3
    col = reshape(col,[size(col,1) 3]);
  end
  for j=1:length(udata.new_marks)
    col(udata.bmus(udata.new_marks(j)),:)=[0 0 0];
  end
  set(udata.tmp_patch(i),'FaceVertexCData',col);
end

set(udata.fig1,'UserData',udata);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function color_gui(fig1)

% construct the graphical user interface for changing the color of the
% black (marked) nodes.

udata=get(fig1,'UserData');

a = figure('Color',[0.8 0.8 0.8], ...
	   'Name','Colors', ...
	   'PaperType','a4letter', ...
	   'Position',[518 456 120 311], ...
	   'Tag','Fig1');

udata.c_struct.fig=a;

b = uicontrol('Parent',a, ...
	      'Units','normalized', ...
	      'BackgroundColor',[0.701961 0.701961 0.701961], ...
	      'Position',[0.0700415 0.28956 0.830492 0.594566], ...
	      'Style','frame', ...
	      'Tag','Frame1');
b = uicontrol('Parent',a, ...
	      'Units','normalized', ...
	      'BackgroundColor',[0.8 0.8 0.8], ...
	      'Position',[0.100059 0.301143 0.770456 0.571399], ...
	      'Style','frame', ...
	      'Tag','Frame2');

b = uicontrol('Parent',a, ...
	      'Units','normalized', ...
	      'Callback','vis_trajgui([],''cyan'')', ...
	      'Position',[0.130077 0.795326 0.170101 0.0617729], ...
	      'Style','radiobutton', ...
	      'Tag','cyan', ...
	      'Value',1);

b = uicontrol('Parent',a, ...
	      'Units','normalized', ...
	      'Callback','vis_trajgui([],''magenta'')', ...
	      'Position',[0.130077 0.733553 0.170101 0.057912], ...
	      'Style','radiobutton', ...
	      'Tag','magenta');

b = uicontrol('Parent',a, ...
	      'Units','normalized', ...
	      'Callback','vis_trajgui([],''yellow'')', ...
	      'Position',[0.130077 0.664059 0.170101 0.0617729], ...
	      'Style','radiobutton', ...
	      'Tag','yellow');

b = uicontrol('Parent',a, ...
	      'Units','normalized', ...
	      'Callback','vis_trajgui([],''red'')', ...
	      'Position',[0.130077 0.590703 0.170101 0.0617729], ...
	      'Style','radiobutton', ...
	      'Tag','red');

b = uicontrol('Parent',a, ...
	      'Units','normalized', ...
	      'Callback','vis_trajgui([],''green'')', ...
	      'Position',[0.130077 0.525068 0.170101 0.057912], ...
	      'Style','radiobutton', ...
	      'Tag','green');

b = uicontrol('Parent',a, ...
	      'Units','normalized', ...
	      'Callback','vis_trajgui([],''blue'')', ...
	      'Position',[0.130077 0.455575 0.170101 0.0617729], ...
	      'Style','radiobutton', ...
	      'Tag','blue');

b = uicontrol('Parent',a, ...
	      'Units','normalized', ...
	      'Callback','vis_trajgui([],''white'')', ...
	      'Position',[0.130077 0.38608 0.170101 0.0617729], ...
	      'Style','radiobutton', ...
	      'Tag','white');

b = uicontrol('Parent',a, ...
	      'Units','normalized', ...
	      'Callback','vis_trajgui([],''grey'')', ...
	      'Position',[0.130077 0.320447 0.170101 0.057912], ...
	      'Style','radiobutton', ...
	      'Tag','grey');

b = uicontrol('Parent',a, ...
	      'Units','normalized', ...
	      'BackgroundColor',[0.8 0.8 0.8], ...
	      'FontWeight','demi', ...
	      'HorizontalAlignment','left', ...
	      'Position',[0.32019 0.795326 0.470278 0.0501905], ...
	      'String','Cyan', ...
	      'Style','text', ...
	      'Tag','StaticText1');
b = uicontrol('Parent',a, ...
	      'Units','normalized', ...
	      'BackgroundColor',[0.8 0.8 0.8], ...
	      'FontWeight','demi', ...
	      'HorizontalAlignment','left', ...
	      'Position',[0.32019 0.733553 0.520308 0.0463296], ...
	      'String','Magenta', ...
	      'Style','text', ...
	      'Tag','StaticText2');
b = uicontrol('Parent',a, ...
	      'Units','normalized', ...
	      'BackgroundColor',[0.8 0.8 0.8], ...
	      'FontWeight','demi', ...
	      'HorizontalAlignment','left', ...
	      'Position',[0.32019 0.664059 0.470278 0.0501905], ...
	      'String','Yellow', ...
	      'Style','text', ...
	      'Tag','StaticText3');
b = uicontrol('Parent',a, ...
	      'Units','normalized', ...
	      'BackgroundColor',[0.8 0.8 0.8], ...
	      'FontWeight','demi', ...
	      'HorizontalAlignment','left', ...
	      'Position',[0.32019 0.590703 0.470278 0.0501905], ...
	      'String','Red', ...
	      'Style','text', ...
	      'Tag','StaticText4');
b = uicontrol('Parent',a, ...
	      'Units','normalized', ...
	      'BackgroundColor',[0.8 0.8 0.8], ...
	      'FontWeight','demi', ...
	      'HorizontalAlignment','left', ...
	      'Position',[0.32019 0.525068 0.470278 0.0463296], ...
	      'String','Green', ...
	      'Style','text', ...
	      'Tag','StaticText5');
b = uicontrol('Parent',a, ...
	      'Units','normalized', ...
	      'BackgroundColor',[0.8 0.8 0.8], ...
	      'FontWeight','demi', ...
	      'HorizontalAlignment','left', ...
	      'Position',[0.32019 0.455575 0.470278 0.0463296], ...
	      'String','Blue', ...
	      'Style','text', ...
	      'Tag','StaticText6');
b = uicontrol('Parent',a, ...
	      'Units','normalized', ...
	      'BackgroundColor',[0.8 0.8 0.8], ...
	      'FontWeight','demi', ...
	      'HorizontalAlignment','left', ...
	      'Position',[0.32019 0.38608 0.470278 0.0501905], ...
	      'String','White', ...
	      'Style','text', ...
	      'Tag','StaticText7');
b = uicontrol('Parent',a, ...
	      'Units','normalized', ...
	      'BackgroundColor',[0.8 0.8 0.8], ...
	      'FontWeight','demi', ...
	      'HorizontalAlignment','left', ...
	      'Position',[0.32019 0.320447 0.470278 0.0463296], ...
	      'String','Grey', ...
	      'Style','text', ...
	      'Tag','StaticText8');
b = uicontrol('Parent',a, ...
	      'Units','normalized', ...
	      'Position',[0.0700415 0.146711 0.830492 0.135128], ...
	      'Style','frame', ...
	      'Tag','Frame3');
b = uicontrol('Parent',a, ...
	      'Units','normalized', ...
	      'BackgroundColor',[0.8 0.8 0.8], ...
	      'Position',[0.100059 0.158293 0.770456 0.111963], ...
	      'Style','frame', ...
	      'Tag','Frame4');
b = uicontrol('Parent',a, ...
	      'Units','normalized', ...
	      'BackgroundColor',[0.8 0.8 0.8], ...
	      'FontWeight','demi', ...
	      'HorizontalAlignment','left', ...
	      'Position',[0.130077 0.177597 0.270833 0.0617729], ...
	      'String','RGB', ...
	      'Style','text', ...
	      'Tag','StaticText9');
b = uicontrol('Parent',a, ...
	      'Units','normalized', ...
	      'BackgroundColor',[1 1 1], ...
	      'Position',[0.410243 0.173736 0.420249 0.0810768], ...
	      'Style','edit', ...
	      'Tag','EditText1');

udata.c_struct.RGB=b;

b = uicontrol('Parent',a, ...
	      'Units','normalized', ...
	      'Callback','vis_trajgui([],''dye'')', ...
	      'FontWeight','demi', ...
	      'Position',[0.0700415 0.0270256 0.360214 0.0772162], ...
	      'String','OK', ...
	      'Tag','Pushbutton1');
b = uicontrol('Parent',a, ...
	      'Units','normalized', ...
	      'Callback','close gcf', ...
	      'FontWeight','demi', ...
	      'Position',[0.54032 0.0270256 0.360214 0.0772162], ...
	      'String','Close', ...
	      'Tag','Pushbutton2');

udata.c_struct.color=[0 1 1];

tmp.fig1=fig1;
set(a,'UserData',tmp);
set(udata.fig1,'UserData',udata);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function dye_nodes(arg)

% takes care of the action, when radiobuttons are pressed
% (or the RGB value is set) in the color_gui -figure.
% It also handles the starting of dying nodes and plots.

udata=get(gcf,'UserData');
udata=get(udata.fig1,'UserData');

switch arg
 case {'cyan','magenta','yellow','red','green','blue','white','grey'}
  h=findobj(get(gcf,'Children'),'Style','radiobutton');
  set(h,'Value',0);
  set(gcbo,'Value',1);
end

switch arg
 case 'cyan'
  RGB = [0 1 1];
 case 'magenta'
  RGB = [1 0 1];
 case 'yellow'
  RGB = [1 1 0];
 case 'red'
  RGB = [1 0 0];
 case 'green'
  RGB = [0 1 0];
 case 'blue'
  RGB = [0 0 1];
 case 'white'
  RGB = [1 1 1];
 case 'grey'
  RGB = [0.4 0.4 0.4];
 case 'dye'

  RGB = get(udata.c_struct.RGB,'String');
  if isempty(RGB)
    dye;
    return;
  else
    str1='The value of RGB must be vector containing three scalars';
    str2='between 0 and 1.';
    color = str2num(RGB);
    set(udata.c_struct.RGB,'String','');
    if isempty(color)
      close gcf;
      udata=rmfield(udata,'c_struct');
      set(udata.fig1,'UserData',udata);
      errordlg([{str1};{str2}]);
      return;
    end
    if ~all([1 3] == size(color)) && ~all([3 1] == size(color))
      close gcf;
      errordlg([{str1};{str2}]);
      udata=rmfield(udata,'c_struct',udata);
      set(udata.fig1,'UserData',udata);
      return;
    end
    if ~isempty(cat(2,find(color>1),find(color<0)))
      close gcf
      errordlg([{str1};{str2}]);
      udata=rmfield(udata,'c_struct',udata);
      set(udata.fig1,'UserData',udata);
      return;
    end
    udata.c_struct.color=color;
    set(udata.fig1,'UserData',udata);
    dye;
    return;
  end
end

udata.c_struct.color=RGB;
set(udata.fig1,'UserData',udata);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function dye()

% dyes black markers in the component planes and in the data figures

udata=get(gcf,'UserData');
udata=get(udata.fig1,'UserData');

inds=unique([udata.all_marks ; udata.new_marks]);

for i=1:length(udata.d_mark);
  for j=1:length(udata.d_mark(i).h)
    if all(get(udata.d_mark(i).h(j),'Color') == [0 0 0])
      set(udata.d_mark(i).h(j),'Color',udata.c_struct.color);
    end
  end
end

if ~isempty(udata.fig2);
  for i=1:length(udata.d_mark2)
    for j=1:length(udata.d_mark2(i).h)
      if all(get(udata.d_mark2(i).h(j),'Color') == [0 0 0])
        set(udata.d_mark2(i).h(j),'Color',udata.c_struct.color);
      end
    end
  end
end

for i=1:length(udata.a_h)
  col=get(udata.tmp_patch(i),'FaceVertexCData');
  for j=1:length(udata.new_marks)
    col(udata.bmus(udata.new_marks(j)),:)=udata.c_struct.color;
  end
  set(udata.tmp_patch(i),'FaceVertexCData',col);
end

udata.all_marks=unique([udata.all_marks;udata.new_marks]);
udata.new_marks=[];
close gcf;
udata=rmfield(udata,'c_struct');
set(udata.fig1,'UserData',udata);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function clear_markers()

% removes markers from the componentplanes and the data figure(s).

udata=get(gcf,'UserData');
udata=get(udata.fig1,'UserData');

for i=1:length(udata.d_mark)
  delete(udata.d_mark(i).h);
  udata.d_mark(i).h=[];
end

for i=1:length(udata.d_mark2)
  delete(udata.d_mark2(i).h);
  udata.d_mark2(i).h=[];
end

col=NaN*get(udata.tmp_patch(1),'FaceVertexCData');
col=reshape(col,[size(col,1) 3]);

for i=1:length(udata.tmp_patch)
  set(udata.tmp_patch(i),'FaceVertexCData',col);
end

udata.new_marks=[];
udata.all_marks=[];

if any(strcmp('c_struct',fieldnames(udata)))
  udata=rmfield(udata,'c_struct');
end

set(udata.fig1,'UserData',udata);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function key_bdf

% moves trajectory and pointer lines, when either of
% the keys '>' or '<' is pressed.

udata=get(gcf,'UserData');
udata=get(udata.fig1,'UserData');

key=get(gcbo,'CurrentCharacter');

% The easiest way to get a new coordinates is to get them from the texts...
% The texts are either '<- x' or 'x ->'

x=get(udata.text1(1),'String');
x=str2num(x(4:length(x)));

if isempty(x)
  x=get(udata.text1(1),'String');
  x=str2num(x(1:length(x)-3));
end

switch(key)
 case '<'
  if x ~= 1
    x= x-1;
  end
 case '>'
  if x ~= getfield(get(get(udata.text1(1),'Parent'),'XLim'),{2})
    x = x+1;
  end
 otherwise
  return;
end

set(udata.l,'XData',[x x]);
if ~isempty(udata.fig2)
  set(udata.l2,'XData',[x x]);
end

delete(udata.text1);
delete(udata.text2);

set_numbers(x);
draw_traj(x);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function click()

switch get(gcf,'SelectionType')
 case 'open'
  return;
 case {'normal','alt'}
  draw_poly;
 case 'extend'
  click_node;
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function click_node()

% takes care of the action, when the middle mouse button is
% pressed (mouse pointer is above some component plane).

udata=get(gcf,'UserData');
udata=get(udata.fig1,'UserData');
new_marks=[];

old=gcf;

NEW=0;
AGAIN = 0;

coords=get(gca,'CurrentPoint');
row=round(coords(1,2));

if strcmp(udata.lattice,'hexa') && ~mod(row,2)
  col = round(coords(1,1) - 0.5);
else
  col = round(coords(1,1));
end

ind = sub2ind(udata.msize,row,col);
new_marks=find(udata.bmus==ind);

if strcmp(get(gcbo,'Tag'),'TmpPatch');

  % if the callback is made via temporary patch object, node is marked
  % (node is black) => the mark is to be removed

  node_color = getfield(get(gcbo,'FaceVertexCData'),{ind,[1:3]});
  AGAIN = 1;
end

for i=1:length(udata.tmp_patch)
  color = get(udata.tmp_patch(i),'FaceVertexCData');
  if length(size(color)) ~= 2
    color = reshape(color,[size(color,1) 3]);
  end
  if all(isnan(color(ind,:)))
    NEW=1;
    color(ind,:)=[0 0 0];
  else
    color(ind,:)=[NaN NaN NaN];
  end
  set(udata.tmp_patch(i),'FaceVertexCData',color);
end

set(0,'CurrentFigure',udata.fig1);

for j=1:length(udata.h)
  subplot(udata.h(j));
  if NEW
    y=getfield(get(udata.d(j),'YData'),{new_marks});
    udata.d_mark(j).h=[udata.d_mark(j).h;plot(new_marks,y,'Color',[0 0 0],...
                                              'LineStyle','none',...
                                              'Marker','o')];
  end
end

if ~isempty(udata.fig2)
  set(0,'CurrentFigure',udata.fig2);
  for j=1:length(udata.h2);
    subplot(udata.h2(j));
    if NEW
      y=getfield(get(udata.d2(j),'YData'),{new_marks});
      udata.d_mark2(j).h=[udata.d_mark2(j).h;plot(new_marks,y,...
						  'LineStyle','none',...
						  'Color','black',...
						  'Marker','o')];
    end
  end
end

if NEW
  udata.new_marks=[udata.new_marks; new_marks];
end

if AGAIN

  % find marks from the data that map to the clicked node. if the color
  % of the mark(s) is the same as the node's color, remove mark(s), else
  % let mark be unchanged.

  for i=1:length(udata.d_mark(1).h)
    if all(node_color==get(udata.d_mark(1).h(i),'Color'))
      tmp1 = get(udata.d_mark(1).h(i),'XData');
      tmp2 = setdiff(tmp1,new_marks);
      if length(tmp1) ~= length(tmp2)
        inds=[];
        for j=1:length(tmp2);
          inds=[inds find(tmp2(j)==tmp1)];
        end
        for j=1:length(udata.d_mark)
          ydata=getfield(get(udata.d_mark(j).h(i),'YData'),{inds});
          set(udata.d_mark(j).h(i),'XData',tmp2,'YData',ydata);
        end
        if ~isempty(udata.fig2)
          for j=1:length(udata.d_mark2)
            ydata=getfield(get(udata.d_mark2(j).h(i),'YData'),{inds});
            set(udata.d_mark2(j).h(i),'XData',tmp2,'YData',ydata);
          end
        end
      end
    end
  end
  udata.new_marks=setdiff(udata.new_marks, new_marks);
  udata.all_marks=setdiff(udata.all_marks,new_marks);
end

set(udata.fig1,'UserData',udata);
set(0,'CurrentFigure',old);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function draw_poly()

udata=get(gcf,'UserData');
udata=get(udata.fig1,'UserData');

if isempty(udata.poly.x)
  if strcmp(get(gcf,'SelectionType'),'alt')
    return;
  end
  udata.poly.THIS = gca;
end

% 'THIS' indicates what was the axes where the polygon was meant to
% drawn. It is not possible to add points, that lie in another axes, to the
% polygon.

if gca ~= udata.poly.THIS
  return;
end

coords(1,1) = getfield(get(gca,'CurrentPoint'),{3});
coords(1,2) = getfield(get(gca,'CurrentPoint'),{1});

udata.poly.x=cat(1,udata.poly.x,coords(2));
udata.poly.y=cat(1,udata.poly.y,coords(1));

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% remove old 'polygon' from axis

delete(udata.poly.h);

switch get(gcf,'SelectionType')
 case 'normal'

  % add point to the 'polygon' and draw it

  for i=1:length(udata.a_h)
    subplot(udata.a_h(i));
    hold on;
    udata.poly.h(i) = plot(udata.poly.x,udata.poly.y,'black',...
			   'EraseMode','xor',...
			   'ButtonDownFcn',...
			   'vis_trajgui([],''click'')',...
			   'LineWidth',2);
  end
 case 'alt'

  % The polygon is ready.

  udata.poly.x=cat(1,udata.poly.x,udata.poly.x(1));
  udata.poly.y=cat(1,udata.poly.y,udata.poly.y(1));

  for i=1:length(udata.a_h)
    subplot(udata.a_h(i));
    udata.poly.h(i) = plot(udata.poly.x,udata.poly.y,'black',...
			   'EraseMode','xor',...
			   'ButtonDownFcn',...
			   'vis_trajgui([],''click'')',...
			   'LineWidth',2);
  end

  tmp=sort(repmat((1:udata.msize(1))',udata.msize(2),1));
  tmp(:,2)=repmat((1:udata.msize(2))',udata.msize(1),1);
  tmp2=tmp;
  if strcmp(udata.lattice,'hexa');
    t=find(~rem(tmp(:,1),2));
    tmp(t,2)=tmp(t,2)+0.5;
  end

  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  % find the nodes that lie inside polygon and change coordinates to
  % linear indices.

  in = find(inpolygon(tmp(:,2),tmp(:,1),udata.poly.x,udata.poly.y));
  in = sub2ind(udata.msize,tmp2(in,1),tmp2(in,2));

  colors=get(udata.tmp_patch(1),'FaceVertexCData');
  colors=reshape(colors,[size(colors,1) 3]);
  tmp=ones(length(in),1);

  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  % set the color of the nodes just selected, black.

  colors(in,:)=tmp*[0 0 0];
  set(udata.tmp_patch,'FaceVertexCData',colors);

  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  % find the points mapping to the nodes from data

  inds = [];
  for i=1:length(in)
    inds=[inds;find(in(i) == udata.bmus)];
  end

  %%%%%%%%%%%%%%%%%%%
  % plot marks to data

  set(udata.fig1,'UserData',udata);
  plot2data(inds);
  udata=get(udata.fig1,'UserData');
  udata.new_marks=union(udata.new_marks,inds);
  delete(udata.poly.h);
  udata.poly.h=[];
  udata.poly.x=[];
  udata.poly.y=[];
  udata.poly.THIS=[];
end

set(udata.fig1,'UserData',udata);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function save_data()

udata=get(gcf,'UserData');
udata=get(udata.fig1,'UserData');

data.points=[];
data.nodes=[];
k=1;

for i=1:length(udata.d_mark(1).h)
  data.points(i).inds=get(udata.d_mark(1).h(i),'XData');
  data.points(i).color=get(udata.d_mark(1).h(i),'Color');
end

color=get(udata.tmp_patch(1),'FaceVertexCData');
color=reshape(color,[size(color,1) 3]);

for i=1:size(color,1)
  if all(~isnan(color(i,:)))
    tmp.ind=i;
    tmp.color=color(i,:);
    data.nodes(k)=tmp;
    k=k+1;
  end
end

answer=inputdlg('Enter the name of the output variable:','',1);

if isempty(answer) || isempty(answer{1})
  msgbox('Output is not set to workspace.');
  return;
else
  assignin('base',answer{1},data);
  sprintf(1,'Struct is set to the workspace as ''%s''.',answer{1});
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function load_data()

answer = inputdlg('Enter the name of the struct to be loaded:','',1);

if isempty(answer) || isempty(answer{1})
  msgbox('Data is not loaded.');
  return;
end

data=evalin('base',answer{1});

if ~isstruct(data)
  errordlg('Input variable must be a struct.');
  return;
end

tmp1 = fieldnames(data);
tmp2 = {'nodes','points'};

for i=1:length(tmp1)
  for j=1:length(tmp2);
    if ~any(strcmp(tmp2{j},tmp1))
      errordlg('Wrong type of struct.');
      return;
    end
  end
end

if ~isempty(data.points)
  tmp1=fieldnames(data.points(1));
end
if ~isempty(data.nodes)
  tmp2=fieldnames(data.nodes(1));
end

for i=1:length(tmp1)
  if ~any(strcmp(tmp1{i},{'inds','color'}))
    errordlg('Wrong type of struct.');
    return;
  end
end

for i=1:length(tmp2)
  if ~any(strcmp(tmp2{i},{'ind','color'}))
    errordlg('Wrong type of struct.');
    return;
  end
end

udata=get(gcf,'UserData');
udata=get(udata.fig1,'UserData');

clear_markers;
remove_traj;

old = gcf;

for i=1:length(data.points)
  for j=1:length(udata.h);
    set(0,'CurrentFigure',udata.fig1);
    subplot(udata.h(j));
    ydata=getfield(get(udata.d(j),'YData'),{data.points(i).inds});
    udata.d_mark(j).h=[udata.d_mark(j).h;...
                       plot(data.points(i).inds,ydata,...
			    'Color',data.points(i).color,...
			    'LineStyle','none',...
			    'Marker','o',...
			    'ButtonDownFcn',...
			    'vis_trajgui([],''line_down'')')];
    if all(data.points(i).color == [0 0 0])
      udata.new_marks=unique([udata.new_marks; (data.points(i).inds)']);
    else
      udata.all_marks=unique([udata.all_marks; (data.points(i).inds)']);
    end
  end
  if ~isempty(udata.fig2)
    set(0,'CurrentFigure',udata.fig2);
    for j=1:length(udata.h2)
      subplot(udata.h2(j));
      ydata=getfield(get(udata.d2(j),'YData'),{data.points(i).inds});
      udata.d_mark2(j).h=[udata.d_mark2(j).h;...
                          plot(data.points(i).inds,ydata,...
			       'Color',data.points(i).color,...
			       'LineStyle','none',...
			       'Marker','o',...
			       'ButtonDownFcn',...
			       'vis_trajgui([],''line_down'')')];
    end
  end
end

set(0,'CurrentFigure',udata.sM_h);
color=get(udata.tmp_patch(1),'FaceVertexCData');
color=reshape(color,[size(color,1) 3]);
for i=1:length(data.nodes)
  color(data.nodes(i).ind,:)=data.nodes(i).color;
end
for i=1:length(udata.tmp_patch);
  set(udata.tmp_patch(i),'FaceVertexCData',color);
end

set(0,'CurrentFigure',old);
set(udata.fig1,'UserData',udata);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function fuzzy_traj(trajStruct,arg)
%function fuzzy_traj(sM_h,sM,sD,interval,arg)

if isempty(arg)

  if strcmp(trajStruct.lattice,'hexa')
    udata.lattice='hexa';
    udata.form=vis_patch('hexa');
  else
    data.lattice='rect';
    udata.form=vis_patch('rect');
  end

  % interval=[1 size(trajStruct.primary_data,1)];

  l=size(udata.form,1);
  dim = size(trajStruct.primary_data,2);
  udata.a_h=[findobj(get(trajStruct.figure,'Children'),'Tag','Uplane');...
             findobj(get(trajStruct.figure,'Children'),'Tag','Cplane')];

  udata.sM_h=trajStruct.figure;
  udata.msize=trajStruct.msize;

  %%%%%%%%%%%%%%%
  %
  % constructing patch that is drawn above every plane in map
  %

  nx = repmat(udata.form(:,1),1,prod(udata.msize));
  ny = repmat(udata.form(:,2),1,prod(udata.msize));

  x_c=reshape(repmat(1:udata.msize(2),l*udata.msize(1),1),l,prod(udata.msize));
  y_c=repmat(repmat(1:udata.msize(1),l,1),1,udata.msize(2));

  if strcmp(udata.lattice,'hexa')
    t = find(~rem(y_c(1,:),2));
    x_c(:,t)=x_c(:,t)+.5;
  end

  x_c=x_c+nx;
  y_c=y_c+ny;

  udata.orig_c=ones(prod(udata.msize),1)*[NaN NaN NaN];
  colors=reshape(udata.orig_c,[1 size(udata.orig_c,1) 3]);
  set(0,'CurrentFigure',trajStruct.figure);

  %%%%%%%%%
  % drawing

  for i=1:length(udata.a_h);
    subplot(udata.a_h(i));
    v=caxis;
    udata.patch_h(i) =patch(x_c,y_c,colors,'EdgeColor','none');
    caxis(v);
  end

  udata.orig_x=get(udata.patch_h(1),'XData');
  udata.orig_y=get(udata.patch_h(1),'YData');

  %  if interval(1) < 1 | interval(2) > size(trajStruct.primary_data,1)
  %    error('Invalid argument ''interval''.');
  %  end

  x=1:size(trajStruct.primary_data,1);
  udata.fig1=figure;
  set(udata.fig1,'KeyPressFcn',...
                 'vis_trajgui([],''move_fuzzy'')');
  for i=1:size(trajStruct.primary_data,2)
    subplot(size(trajStruct.primary_data,2),1,i);
    udata.h(i)=gca;
    set(udata.h(1),'XTick',[]);
    udata.d(i)=plot(x,trajStruct.primary_data(:,i));
    l_x=1;
    lims(1) = round(min(trajStruct.primary_data(:,i)));
    lims(2) = round(max(trajStruct.primary_data(:,i)));
    udata.l(i) = line([l_x l_x],lims,'Color','red','EraseMode','xor');
  end

  udata.l_x = l_x;
  udata.interval=[1 size(trajStruct.bmus,2)];

  tmp.fig1=udata.fig1;
  %  [K,P] = estimate_kernels(sM,sD);
  %  udata.K=K;
  %  udata.P=P;
  %  udata.codebook=sM.codebook;
  %  udata.data=sD.data;
  udata.bmus=trajStruct.bmus;
  set(udata.fig1,'UserData',udata);
  set(trajStruct.figure,'UserData',tmp);
  draw_fuzzy(l_x);
  return;
end

move_fuzzy;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function draw_fuzzy(x)

udata=get(gcf,'UserData');
udata=get(udata.fig1,'UserData');

inds = find(udata.bmus(:,x));
[row col] = ind2sub(udata.msize,inds);
if strcmp(udata.lattice,'hexa')
  t=find(~mod(row,2));
  col(t)=col(t)+0.5;
end

color=udata.orig_c;
color=reshape(color,[size(color,1) 3]);
xdata=udata.orig_x;
ydata=udata.orig_y;
tmp= ones(size(xdata(:,1),1),1)*udata.bmus(inds,x)';
color(inds,:) = ones(length(inds),1)*[0 0 0];
xdata(:,inds) = udata.form(:,1)*ones(1,length(inds)).*tmp+ones(6,1)*col';
ydata(:,inds) = udata.form(:,2)*ones(1,length(inds)).*tmp+ones(6,1)*row';

set(udata.patch_h,'FaceVertexCData',color,'XData',xdata,'YData',ydata);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function move_fuzzy

% moves pointer lines and draws fuzzy response.

udata=get(gcf,'UserData');
udata=get(udata.fig1,'UserData');

switch get(gcf,'CurrentCharacter');
 case {'<','>'}
  key = get(gcf,'CurrentCharacter');
  if key == '>'
    if udata.l_x + 1 > udata.interval(2)
      return;
    end
    l_x = udata.l_x + 1;
  else
    if udata.l_x - 1 < udata.interval(1)
      return;
    end
    l_x = udata.l_x - 1;
  end
  draw_fuzzy(l_x);
  set(udata.l,'XData',[l_x l_x]);
  udata.l_x=l_x;
  set(udata.fig1,'UserData',udata);
 otherwise
  return;
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

author: ilarinieminen

repo: SOM-Toolbox


text:
function h = p13_h ( n, x )

%*****************************************************************************80
%
%% P13_H evaluates the Hessian for problem 13.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    18 October 2011
%
%  Author:
%
%    John Burkardt
%
%  Parameters:
%
%    Input, integer N, the number of variables.
%
%    Input, real X(N), the values of the variables.
%
%    Output, real H(N,N), the N by N Hessian matrix.
%
  h = zeros ( n, n );

  s1 = sum ( cos ( x(1:n) ) );

  for j = 1 : n
    h(j,j) = sin ( x(j) );
  end

  s2 = 0.0;

  for j = 1 : n

    th = cos ( x(j) );
    t = ( n + j ) - h(j,j) - s1 - j * th;
    s2 = s2 + t;
    for k = 1 : j - 1
      h(j,k) = 2.0 * ( sin ( x(k) ) * ( ( n + j + k ) * h(j,j) - th ) - ...
        h(j,j) * cos ( x(k) ) );
    end

    h(j,j) =  ( j * ( j + 2 ) + n ) * h(j,j) * h(j,j) + th * ...
      ( th -  ( 2 * j + 2 ) * h(j,j) ) + t * ( j * th + h(j,j) );

  end

  for j = 1 : n
    h(j,j) = 2.0 * ( h(j,j) + cos ( x(j) ) * s2 );
  end

  for i = 1 : n
    h(i,i+1:n) = h(i+1:n,i);
  end

  return
end

author: johannesgerer

repo: jburkardt-m


text:
function test_ft_plot_topo

% MEM 1gb
% WALLTIME 00:10:00
% DEPENDENCY ft_plot_topo

x = randn(30,1);
y = randn(30,1);

% make something round between 0 and 1
val = -sqrt(x.^2 + y.^2);
val = val-min(val);
val = val./max(val);
isolines = 0:.1:1;

% 1 default style (surfiso) (without isolines is equal to style surf)
figure
ft_plot_topo(x, y, val);
hold on; plot(x, y, 'k*')
axis tight

% 2 style surfiso shading interp
figure
ft_plot_topo(x, y, val, 'style', 'surfiso', 'shading', 'interp');
hold on; plot(x, y, 'k*')
axis tight

% 3 default style (surfiso) with isolines
figure
ft_plot_topo(x, y, val, 'isolines', isolines);
hold on; plot(x, y, 'k*')
axis tight

% 4 style isofill with isolines
figure
ft_plot_topo(x, y, val, 'style', 'isofill', 'isolines', isolines);
hold on; plot(x, y, 'k*')
axis tight

%%% MASKING
% 5 style surf with binary mask
figure
ft_plot_topo(x, y, val, 'datmask', val>0.5, 'style', 'surf', 'clim', [0 1]);
hold on; plot(x, y, 'k*')
axis tight

% 6 style imsat with binary mask
figure
ft_plot_topo(x, y, val, 'datmask', val>0.5, 'style', 'imsat', 'clim', [0 1]);
hold on; plot(x, y, 'k*')
axis tight

% 7 style surf with continuous mask  (part of image will always be masked out because val is always 0 in one spot)
figure
ft_plot_topo(x, y, val, 'datmask', val, 'style', 'surf', 'clim', [0 1]);
hold on; plot(x, y, 'k*')
axis tight

% 8 style imsat with continuous mask  (part of image will always be masked out because val is always 0 in one spot)
figure
ft_plot_topo(x, y, val, 'datmask', val, 'style', 'imsat', 'clim', [0 1]);
hold on; plot(x, y, 'k*')
axis tight

% 9 style imsatiso with continuous mask and isolines (part of image will always be masked out because val is always 0 in one spot)
figure
ft_plot_topo(x, y, val, 'datmask', val, 'style', 'imsatiso', 'clim', [0 1], 'isolines', isolines);
hold on; plot(x, y, 'k*')
axis tight

author: fieldtrip

repo: fieldtrip

save_path: github-repos/MATLAB/fieldtrip-fieldtrip


text:
function res = bf_inverse_nutmeg(BF, S)
% Interface to NUTMEG inverse methods
% http://www.nitrc.org/plugins/mwiki/index.php/nutmeg:MainPage
% Copyright (C) 2015 Wellcome Trust Centre for Neuroimaging

% $Id: bf_inverse_nutmeg.m 7703 2019-11-22 12:06:29Z guillaume $

%--------------------------------------------------------------------------
if nargin == 0
    method = cfg_menu;
    method.tag = 'method';
    method.name = 'Method';
    method.labels = {
        'sLORETA'
        'swLORETA'
        'dSPM'};
    method.values = method.labels;
    method.help = {'Select one of NUTMEG methods'};

    snr = cfg_entry;
    snr.tag = 'snr';
    snr.name = 'SNR';
    snr.strtype = 'e';
    snr.num = [0 0];
    snr.val = {[]};
    snr.help = {'The assumed ratio of variances of signal and noise,',...
        'used for setting the regularisation parameter.'};

    regularisation = cfg_entry;
    regularisation.tag = 'regularisation';
    regularisation.name = 'Regularisation parameter';
    regularisation.strtype = 'e';
    regularisation.num = [0 0];
    regularisation.val = {[]};
    regularisation.help = {'Optional regularization parameter.'};

    nutmeg      = cfg_branch;
    nutmeg.tag  = 'nutmeg';
    nutmeg.name = 'NUTMEG methods';
    nutmeg.val  = {method, snr, regularisation};
    nutmeg.help = {'Methods ported from NUTMEG toolbox',...
        'See http://www.nitrc.org/plugins/mwiki/index.php/nutmeg:MainPage'};

    res = nutmeg;

    return
elseif nargin < 2
    error('Two input arguments are required');
end

C     =  BF.features.(S.modality).C;
U     =  BF.features.(S.modality).U;

L = S.L;
W = cell(size(L));
nvert = numel(W);

data.Ryy = C;
flags    = [];
if ~isempty(S.snr)
    flags.snr = S.snr;
end
if ~isempty(S.regularisation)
    flags.gamma = S.regularisation;
end

LL = [];
spm('Pointer', 'Watch');drawnow;
spm_progress_bar('Init', nvert, ['Preparing ' S.modality ' leadfields']); drawnow;
if nvert > 100, Ibar = floor(linspace(1, nvert,100));
else Ibar = 1:nvert; end

for i = 1:nvert
    lf = U'*L{i};

    LL = cat(3, LL, lf);

    if ismember(i, Ibar)
        spm_progress_bar('Set', i); drawnow;
    end
end

switch S.method
    case 'sLORETA'
        w = nut_sLORETA(LL,data,flags);
    case 'swLORETA'
        w = nut_swLORETA(LL,data,flags);
    case  'dSPM'
        w = nut_dSPM(LL,data,flags);
end

spm_progress_bar('Init', nvert, ['Preparing ' S.modality ' filters']); drawnow;
if nvert > 100, Ibar = floor(linspace(1, nvert,100));
else Ibar = 1:nvert; end

for i = 1:nvert
    W{i} = spm_squeeze(w(:, :, i), 3)';

    if ismember(i, Ibar)
        spm_progress_bar('Set', i); drawnow;
    end
end

spm_progress_bar('Clear');

res.W = W;

author: spm


text:
function [imageVol, lastfile] = GE_readVolume(startDir, passnum, volSize, depth, im_offset)
%
%GE_readVolume
%
% [imageVol, lastfile] = GE_readVolume(startDir, passnum, [nX nY nZ], depth, im_offset)
%
% reads the volume for passnum from the series which is stored
% starting in startDir and returns the name of the last file read
%
% Souheil J. Inati
% Dartmouth College
% May 2000
% souheil.inati@dartmouth.edu
%

% initialize some variables
nX = volSize(1);
nY = volSize(2);
nZ = volSize(3);
sliceSize = nX*nY;
imageVol = zeros(nX, nY, nZ);
[path_stem, path_start] = fileparts(startDir);

for i = 1:nZ
        % Make the filename
        filenum = (passnum-1)*nZ + i;   % file no.
        filep = ceil(filenum/999) - 1;
        filen = filenum - 999*filep;
        path_num = str2num(path_start) + 20*filep;
        path_now = sprintf('00%d',path_num);
        disp(['DEBUD: ',path_now]);
        path_now = path_now(length(path_now)-2:length(path_now));
        path = fullfile(path_stem, path_now);
	stub = sprintf('00%d',filen);
	stub = stub(length(stub)-2:length(stub));
        imageFile = fullfile(path,['I.' stub]);

        % Open the file
	[fid,message] = fopen(imageFile,'r','b');
	if (fid == -1)
          fprintf('Cannot Open %s.\n',imageFile);
          passnum
	  break
	end

	% Skip to the data
	fseek(fid,im_offset,-1);
	% Read the slice
	buffer = fread(fid,sliceSize,sprintf('int%d',depth));
	% append the slice to the imageSet
	imageVol(:,:,i) = reshape(buffer, nX, nY);

        % Close the file
	status = fclose(fid);

        % Set the lastfile
        lastfile = imageFile;
end

return

author: vistalab

repo: vistasoft

save_path: github-repos/MATLAB/vistalab-vistasoft


text:
function compare_implementations()

n_trials = 50;
n_obstacles = 5;
matlab_times = nan(1, n_trials);
cpp_times = nan(1, n_trials);

iter = 1;
while iter <= n_trials
  dim = 1 + ceil(iter / 20);
  A_bounds = [eye(dim); -eye(dim)];
  b_bounds = [ones(dim, 1); zeros(dim, 1)];

  obstacle_pts = iris.test.random_obstacles(dim, n_obstacles, zeros(dim,1), ones(dim,1), 0.05);
  obstacles = mat2cell(obstacle_pts, size(obstacle_pts, 1), size(obstacle_pts, 2), ones(1, size(obstacle_pts, 3)));
  start = 0.5 * ones(dim,1);
  options = struct('require_containment', false,...
                   'error_on_infeasible_start', false,...
                   'termination_threshold', 2e-2,...
                   'iter_limit', 100);

  try
    t0 = tic();
    [A_c, b_c, C_c, d_c] = iris.inflate_regionmex(obstacles, A_bounds, b_bounds, start, options);
    cpp_times(iter) = toc(t0);

    t0 = tic();
    [A_m, b_m, C_m, d_m] = iris.inflate_region_fallback(obstacle_pts, A_bounds, b_bounds, start, options);
    matlab_times(iter) = toc(t0);

    iter = iter + 1
  catch
    continue
  end

  % assert(iris.util.equal_up_to_permutations(A_c, A_m, 1e-2));
  % assert(iris.util.equal_up_to_permutations(b_c, b_m, 1e-2));
  assert(all(all(abs(C_c - C_m) < 1e-2)));
  assert(all(all(abs(d_c - d_m) < 1e-2)));

end

fprintf(1, 'c++ mean: %f std: %f\n', mean(cpp_times), std(cpp_times));
fprintf(1, 'matlab mean: %f std: %f\n', mean(matlab_times), std(matlab_times));

author: rdeits

repo: iris-distro

save_path: github-repos/MATLAB/rdeits-iris-distro


text:
function [bankAng,angOfAttack,angOfSideslip] = computeInertialAeroAnglesFromBodyAxes(ut, rVect, vVect, bodyInfo, bodyX, bodyY, bodyZ)
    %Source: http://www.dept.aoe.vt.edu/~cdhall/courses/aoe5204/AircraftMotion.pdf

    [R_wind_2_inert, ~, ~, ~] = computeWindFrame(rVect,vVect);
    Rtotal = horzcat(bodyX, bodyY, bodyZ);

    angles = rotm2eulARH(R_wind_2_inert' * Rtotal, 'zyx');

    angles = real(angles);

    bankAng = angles(3);
	angOfAttack = angles(2);
	angOfSideslip = angles(1);
end

author: Arrowstar

repo: ksptot

save_path: github-repos/MATLAB/Arrowstar-ksptot


text:
%convert images to double and grayscale
I1 = impyramid( imInit('car1.png'), 'reduce' );
I2 = impyramid( imInit('car2.png'), 'reduce' );

%%

[flowHor flowVer] = pyramidFlow(I1, I2, 5, 3, 3);  %pyramidFlow( I1, I2, winSize, ITER_NO, PYRE_NO )

%show the output flows
imS(flowHor,1,[-10 10])
imS(flowVer,2,[-10 10])

%%
%show the warped image to check the quality of registration
Iw = imWarp(flowHor, flowVer, I2);
imS(I1,10)
imS(Iw,11)
imS(I2,12)

author: Sable

repo: mcbench-benchmarks


text:
function varargout = ubcread(varargin)
% SIFTREAD Read Lowe's SIFT implementation data files
%   [F,D] = VL_UBCREAD(FILE) reads the frames F and the descriptors D
%   from FILE in UBC (Lowe's original implementation of SIFT) format
%   and returns F and D as defined by VL_SIFT().
%
%   VL_UBCREAD(FILE, 'FORMAT', 'OXFORD') assumes the format used by
%   Oxford VGG implementations .
%
%   See also: VL_SIFT(), VL_HELP().
[varargout{1:nargout}] = vl_ubcread(varargin{:});

author: yihui-he

repo: panorama

save_path: github-repos/MATLAB/yihui-he-panorama


text:
%System model for geodetic navigation frame mechanization for both phi and psi parametrization.

%This script can be used for both quaternion and DCM attitude parametrizations.

%Position errors are represented as delta_r^n (position error with respect to
%earth defined in navigation frame) in meters rather than
%Delta-Llh (llh perturbation model see:sys_llh_phipsi) or Delta-Theta (Cen perturbation: to be added later)

%Delta_r^n, delta_llh and Delta_theta are assumed to be convertible to each
%other. That is why this system model can be used for any of these position
%mechanizations.

%There is no difference between the psi form of this script and the "sys_wander_psi". This is because
%theta_z does not appear anywhere in the error propagation model for the
%psi model. (Ccn*dv^c is assumed to be equal to dv^c in position error model
%and theta*g is independent of theta_z). Therefore, you can use
%sys_wander_psi instead of this script. (Also note that sys_wander_psi is a
%nonsingular implementation, whereas this script is singular at the pole as geoparam is called with arg 2)

%Returns F and Q of discrete time model. The IMU error model must be
%specified here. (Note that in the previous implementations I had seperated
%nav and imu models in different m-files. I no longer follow that approach)

%mode==1 --> PHI model
%mode==2 --> PSI model

%imutype is directly transferred to imu_err_defs_v000

function [STM Qd]=sys_metric_phipsi(Cen, h, Vn, att, acc, gyro, dt, mode, imutype, modelNo)

%%%Navigation error model
%1-3:Position error in meters
%4-6:vel errors
%7-9:Attitude errors in psi form

%N[da:dw]':Effect of IMU errors on nav states
if(nargin<10)
    modelNo=3; % random walk bias and scale factor
end
Anav=zeros(9); %pos, vel, attitude
Nnav=zeros(9,6);

[Fc, wen_n, wie_n, g]=geoparam_v001(2, Cen(:,3), h, Vn);

if (size(att,1)==1 || size(att,2)==1) %att is qbn
    Cbn=quat2dcm_v000(att);
else %att is dcm
    Cbn=att;
end

%%%Navigation System model
%%Common terms for both phi and psi
%position errors
Anav(1:3,1:3)=-skew(wen_n);
Anav(1:3,4:6)=eye(3);

%Velocity Errors
Anav(6,3)=2*g*Fc(1,1);  %effect of height on gravity (note that, drz is the negative of dh. Therefore, this is positive)
Anav(4:6,4:6)=-skew(wen_n+2*wie_n);
Anav(4:6,7:9)=skew(Cbn*acc);

Nnav(4:6,1:3)=Cbn; %%effect of accelerometer errors

%Attitude errors
Anav(7:9,7:9)=-skew(wen_n+ wie_n);
Nnav(7:9,4:6)=-Cbn; %Effect of gyroscope errors

%PHI model
if (mode==1)
    sL=-Cen(3,3);
    cL=Cen(3,1);    %Note that I assume wander=0 for Cen. (This is not a non-singular implementation)
    %%position errors
    Anav(1:3,1:2)=Anav(1:3,1:2)-skew(Vn)*[0 Fc(1,1);-Fc(2,2) 0; 0 Fc(3,1)];

    %%Velocity errors
    %Effect of 2*d_wie_n
    Anav(4:6,1)=Anav(4:6,1)+2*skew(Vn)*[wie_n(3);0;-wie_n(1)]*Fc(2,2);
    %Effect of d_w_en_n
    Anav(4:6,3:5)=Anav(4:6,3:5)+skew(Vn)*[Vn(2)*Fc(1,1)^2 0 Fc(1,1);-Vn(1)*Fc(2,2)^2 -Fc(2,2) 0;-Vn(2)*Fc(1,1)^2*sL/cL 0 Fc(3,1)];
    Anav(4:6,1)=Anav(4:6,1)+skew(Vn)*[0;0;-(Vn(2)*Fc(1,1)/cL/cL)]*Fc(2,2);

    %%Attitude errors
    %effect of dw_ie_n
    Anav(7:9,1)=Anav(7:9,1)+[wie_n(3);0;-wie_n(1)]*Fc(2,2);

    %effect of dw_en_n
    Anav(7:9,3:5)=Anav(7:9,3:5)+[-Vn(2)*Fc(1,1)^2 0 Fc(1,1);Vn(1)*Fc(2,2)^2 -Fc(2,2) 0;Vn(2)*Fc(1,1)^2*sL/cL 0 Fc(3,1)];
    Anav(9,1)=Anav(9,1)-(Vn(2)*Fc(1,1)/cL/cL)*Fc(2,2);
end

%PSI model
if (mode==2)
    %Velocity errors
    Anav(4:5,1:2)=Anav(4:5,1:2)+g*[-Fc(2,2) 0;0 -Fc(1,1)];    %effect of horizontal position error on gravity
end

%%%%Imu error model parameters
[Aimu_d, Qimu_d, Cimu, Rimu]=imu_err_model_v001(acc, gyro, dt, imutype, modelNo);

%%Convert continuous time Anav into discrete time
%I used van loan's method here to perform conversion so that I can obtain a
%correct Qnav_d. If you use a 1st order taylor series approximation, do not
%forget to correct the cross-correlation effects caused by arw/vrw. (In,
%the final implementation, it would be better to manually compute each
%element rather than using this method)
mx_a=dt*[-Anav,Nnav*Rimu*Nnav';zeros(9),Anav'];
mx_b = expm(mx_a);
Anav_d = mx_b(10:18,10:18)';
Qnav_d = Anav_d*mx_b(1:9,10:18);

%%Combine everything
nst_imu=size(Aimu_d);
STM=zeros(9+nst_imu);
Qd=zeros(9+nst_imu);

STM(1:9,1:9)=Anav_d;
STM(1:9,10:end)=Nnav*Cimu*dt;
STM(10:end,10:end)=Aimu_d;

Qd(1:9,1:9)=Qnav_d;
Qd(10:end,10:end)=Qimu_d;
Qd(1:9,10:end)=Nnav*Cimu*Qimu_d*dt/2;   %not necessary, but let's keep this
Qd(10:end,1:9)=Qd(1:9,10:end)';
end

repo: ekfmonoslam


text:
function newState = mrvPanelToggle(panel, state);
%
% newState = mrvPanelToggle(panel, [state]);
%
% Toggles the visibility of an mrvPanel, rescaling the parent
% figure and other objects in the figure, such that the other
% objects remain the same size and the figure grows/shrinks.
%
% panel is a handle to the mrvPanel.
%
% state can be:
%   handle to a uimenu: toggle the 'checked' status of the
%       uimenu and, depending on whether the menu has become
%       checked or unchecked, set the panel visibility
%       appropriately.
%   'on' or 'off': force the respective visibility state.
%   omitted: toggles based on the 'visible' property of the
%       panel.
%
% Returns the new visibility state (1 for on, 0 for off) of
% the panel.
%
% ras, 07/06/05.
% ras, 07/15/05. Now deals w/ panels that reside in their
% own figures.
if nargin<1, help(mfilename); error('Not enough args.'); end

% parse the state argument; figure out new state
if ~exist('state','var') | isempty(state)
    % get from panel's visibility property
    currState = get(panel,'Visible');
    if isequal(currState,'on'),
        state = 'off';
    else
        state = 'on';
    end
    newState = (isequal(state,'on'));

elseif ischar(state)
    if ~ismember(lower(state),{'on' 'off'})
        help(mfilename);
        error('Invalid state argument.')
    end
    state = lower(state);
    newState = (isequal(state,'on'));

elseif ishandle(state)
    % shuffle the variables such that
    % state is a string 'on' or 'off',
    % and hmenu points to the handle
    hmenu = state;
    newState = umtoggle(hmenu);
    if newState==1, state = 'on'; else, state = 'off'; end

else
    error('Invalid state argument.')

end

% first, if the panel resides in its own figure,
% just hide/show the figure
if isequal(get(panel, 'Position'), [0 0 1 1])
    set( get(panel, 'Parent'), 'Visible', state );
    return
end

% if the panel is already in the specified state, return quietly
if isequal( get(panel, 'Visible'), state )
    return
end

% show/hide panel
set(panel,'Visible',state);

% the location of the panel w.r.t. other objects is
% kept as the panel's user data
loc = lower(get(panel,'UserData'));

% get the parent object
parent = get(panel,'Parent');

% get all existing child objects of the parent object,
% so we can adjust their relative sizes
objs = findobj('Parent', parent, 'Type', 'uicontrol');
objs = [objs; findobj('Parent', parent, 'Type', 'axes')];
objs = [objs; findobj('Parent', parent, 'Type', 'uipanel')];
objs = setdiff(objs, panel); % don't alter panel any further

% get the dimension along which things will be scaled
switch loc
    case {'above', 'below'}, dim = 4;
    case {'left', 'right'}, dim = 3;
end
corner = dim-2;

%%%%%figure out scaling sizes
figUnits = get(parent, 'Units');
panelUnits = get(panel, 'Units');
set(parent,'Units', 'normalized'); % note these norm units are diff't --
set(panel,'Units', 'normalized');  % parent is rel to screen, panel is
figPos = get(parent, 'Position');   % rel to parent. Need to be careful
figSize = figPos(dim);              % about math.
panelPos = get(panel, 'Position');
panelSize = panelPos(dim);

%%%%%if showing panel, re-initialize its position
if isequal(state, 'on')
    switch loc
        case 'above', pos = [0 1-panelSize 1 panelSize];
        case 'below', pos = [0 0 1 panelSize];
        case 'left', pos = [0 0 panelSize 1];
        case 'right', pos = [1-panelSize 0 panelSize 1];
    end
    set(panel,'Position',pos);
end

if isequal(state, 'on')
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    % show panel: make figure larger, shrink other objs %
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    % extent to grow figure such that the panel size
    % is the right proportion:
    growFactor = 1 / (1/panelSize - 1);

    % resize the parent object -- make larger
    newSize = figSize * (1+growFactor); % new figure size
    newPos = figPos;
    newPos(dim) = newSize;
    if ismember(loc,{'below' 'left'})
        % move the corner of the parent back to accommodate panel
        newPos(corner) = figPos(corner) - (newSize - figPos(dim)); % figPos(dim)*(1-growFactor);
    end
    set(parent,'Position',newPos);

    % existing objects in the parent will need to be shrunk
    % along the relevant dimension by an amount corresponding
    % to the growth of the parent (e.g., they'll stay the same
    % absolute size, but be a smaller proportion of the parent):
    for i = 1:length(objs)
        exUnits = get(objs(i), 'Units');
        set(objs(i), 'Units', 'normalized');
        objPos = get(objs(i), 'Position');

        % shrink object along relevant dim
        objPos(dim) = objPos(dim) * (1-panelSize);
        objPos(dim) = min(1,objPos(dim)); % off-figure check

        % to keep existing panels flush with one another,
        % even after shrinking, we need to nudge panels back:
        if objPos(corner) > 0
            objPos(corner) = objPos(corner) * (1-panelSize);
        end

        % for these locations, need to move the corners
        % of existing objects forward as well:
        if ismember(lower(loc), {'left' 'below'})
            objPos(corner) = objPos(corner) + panelSize;
        end

        objPos(corner) = max(0,objPos(corner)); % off-figure check
        set(objs(i), 'Position', objPos);
        set(objs(i), 'Units', exUnits); % restore previous unit convention
    end

else
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    % hide panel: make figure smaller, grow other objs  %
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

    % existing objects in the parent will need to be grown
    % along the relevant dimension by an amount corresponding
    % to the size decrease of the parent:
    for i = 1:length(objs)
        exUnits = get(objs(i),'Units');
        set(objs(i),'Units','normalized');
        objPos = get(objs(i),'Position');

        % grow object along relevant dim
        objPos(dim) = objPos(dim) / (1-panelSize);
        objPos(dim) = min(1,objPos(dim)); % off-figure check

        % for these locations, need to move the corners
        % of existing objects backward as well
        if objPos(corner) > panelPos(corner)
            objPos(corner) = (objPos(corner)-panelSize);
        end

        % to keep existing panels flush with one another,
        % even after growing, we need to nudge panels forward:
        if objPos(corner) > 0
            objPos(corner) = objPos(corner) / (1-panelSize);
        end

        objPos(corner) = max(0,objPos(corner)); % off-figure check
        set(objs(i),'Position',objPos);
        set(objs(i),'Units',exUnits); % restore previous unit convention
    end

    % resize the parent object -- make smaller
    newSize = figSize * (1-panelSize);
    newPos = figPos;
    newPos(dim) = newSize;
    if ismember(loc,{'below' 'left'})
        newPos(corner) = newPos(corner) + figSize*panelSize;
    end
    set(parent,'Position',newPos);
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% A quick fix for a complicated problem:
% if many panels are added in many locations, (e.g. some above,
% some to the right), and they are toggled in a random order,
% some panels may be made too large. To help correct for this,
% correct panels which are larger than the parent object.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for i = 1:length(objs)
    exUnits = get(objs(i),'Units');
    set(objs(i),'Units','normalized');
    objPos = get(objs(i),'Position');
    if objPos(3)+objPos(1) > 1
        objPos(3) = 1-objPos(1);
        set(objs(i),'Position',objPos);
    end
    if objPos(4)+objPos(2) > 1
        objPos(4) = 1-objPos(2);
        set(objs(i),'Position',objPos);
    end
    set(objs(i),'Units',exUnits); % restore previous unit convention
end

% restore previous units
set(parent,'Units',figUnits);
set(panel,'Units',panelUnits);

figOffscreenCheck(parent);

return


text:
% This script evaluates the percentage of space time coevered by
%alarms
%
re = [];

% Stefan Wiemer    4/95

report_this_filefun(mfilename('fullpath'));

global abo;
abo = abo2;

for tre2 = min(abo(:,4))+1.5:0.2:max(abo(:,4)-0.1)
    tre2
    abo = abo2;
    abo(:,5) = abo(:,5)* par1/365 + a(1,3);
    l = abo(:,4) >= tre2;
    abo = abo(l,:);
    l = abo(:,3) < tresh;
    abo = abo(l,:);
    size(abo)
    hold on

    j = 0;
    tmp = abo;

    while length(abo) > 1
        j = j+1;
        global iala
        [k,m] = findnei(1);
        po = [k];
        for i = 1:length(k)
            [k2,m2]  = findnei(k(i));
            po = [po ; k2];
            po = sort(po);
            po2 = [0;  po(1:length(po)-1)] ;
            l = find(po-po2 > 0) ;
            po = [po(l) ] ;
        end
        do = ['an' num2str(j) ' = abo(po3,:);'];
        disp([num2str(j) '  Anomalie groups  found'])
        eval(do)
        abo(po3,:) =[];
    end   % while j

    re = [re ; tre2 j ];
end   % for tre2

figure

matdraw
axis off

uicontrol('Units','normal',...
    'Position',[.0 .65 .08 .06],'String','Save ',...
     'Callback',{@calSave9, re(:,1), re(:,2)})

rect = [0.20,  0.10, 0.70, 0.60];
axes('position',rect)
hold on
pl = plot(re(:,1),re(:,2),'r');
set(pl,'LineWidth',1.5)
pl = plot(re(:,1),re(:,2),'ob');
set(pl,'LineWidth',1.5,'MarkerSize',10)

set(gca,'visible','on','FontSize',ZmapGlobal.Data.fontsz.m,'FontWeight','bold',...
    'FontWeight','bold','LineWidth',1.5,...
    'Box','on')
grid

ylabel('Number of Alarm Groups')
xlabel('Zalarm ')
watchoff

author: CelsoReyes


text:
% DERIVEST demo script

% This script file is designed to be used in cell mode
% from the matlab editor, or best of all, use the publish
% to HTML feature from the matlab editor. Older versions
% of matlab can copy and paste entire blocks of code into
% the Matlab command window.

% DERIVEST is property/value is driven for its arguments.
% Properties can be shortened to the

%% derivative of exp(x), at x == 0
[deriv,err] = derivest(@(x) exp(x),0)

%% DERIVEST can also use an inline function
[deriv,err] = derivest(inline('exp(x)'),0)

%% Higher order derivatives (second derivative)
% Truth: 0
[deriv,err] = derivest(@(x) sin(x),pi,'deriv',2)

%% Higher order derivatives (third derivative)
% Truth: 1
[deriv,err] = derivest(@(x) cos(x),pi/2,'der',3)

%% Higher order derivatives (up to the fourth derivative)
% Truth: sqrt(2)/2 = 0.707106781186548
[deriv,err] = derivest(@(x) sin(x),pi/4,'d',4)

%% Evaluate the indicated (default = first) derivative at multiple points
[deriv,err] = derivest(@(x) sin(x),linspace(0,2*pi,13))

%% Specify the step size (default stepsize = 0.1)
deriv = derivest(@(x) polyval(1:5,x),1,'deriv',4,'FixedStep',1)

%% Provide other parameters via an anonymous function
% At a minimizer of a function, its derivative should be
% essentially zero. So, first, find a local minima of a
% first kind bessel function of order nu.
nu = 0;
fun = @(t) besselj(nu,t);
fplot(fun,[0,10])
x0 = fminbnd(fun,0,10,optimset('TolX',1.e-15))
hold on
plot(x0,fun(x0),'ro')
hold off

deriv = derivest(fun,x0,'d',1)

%% The second derivative should be positive at a minimizer.
deriv = derivest(fun,x0,'d',2)

%% Compute the numerical gradient vector of a 2-d function
% Note: the gradient at this point should be [4 6]
fun = @(x,y) x.^2 + y.^2;
xy = [2 3];
gradvec = [derivest(@(x) fun(x,xy(2)),xy(1),'d',1), ...
           derivest(@(y) fun(xy(1),y),xy(2),'d',1)]

%% Compute the numerical Laplacian function of a 2-d function
% Note: The Laplacian of this function should be everywhere == 4
fun = @(x,y) x.^2 + y.^2;
xy = [2 3];
lapval = derivest(@(x) fun(x,xy(2)),xy(1),'d',2) + ...
           derivest(@(y) fun(xy(1),y),xy(2),'d',2)

%% Compute the derivative of a function using a central difference scheme
% Sometimes you may not want your function to be evaluated
% above or below a given point. A 'central' difference scheme will
% look in both directions equally.
[deriv,err] = derivest(@(x) sinh(x),0,'Style','central')

%% Compute the derivative of a function using a forward difference scheme
% But a forward scheme will only look above x0.
[deriv,err] = derivest(@(x) sinh(x),0,'Style','forward')

%% Compute the derivative of a function using a backward difference scheme
% And a backward scheme will only look below x0.
[deriv,err] = derivest(@(x) sinh(x),0,'Style','backward')

%% Although a central rule may put some samples in the wrong places, it may still succeed
[d,e,del]=derivest(@(x) log(x),.001,'style','central')

%% But forcing the use of a one-sided rule may be smart anyway
[d,e,del]=derivest(@(x) log(x),.001,'style','forward')

%% Control the behavior of DERIVEST - forward 2nd order method, with only 1 Romberg term
% Compute the first derivative, also return the final stepsize chosen
[deriv,err,fdelta] = derivest(@(x) tan(x),pi,'deriv',1,'Style','for','MethodOrder',2,'RombergTerms',1)

%% Functions should be vectorized for speed, but its not always easy to do.
[deriv,err] = derivest(@(x) x.^2,0:5,'deriv',1)
[deriv,err] = derivest(@(x) x^2,0:5,'deriv',1,'vectorized','no')


text:
%NONNEWTONIAN Non-Newtonian flow physics mode.

% Copyright 2013-2022 Precise Simulation, Ltd.

author: precise-simulation

repo: featool-multiphysics


text:
classdef GaussianMixtureD
%%GAUSSIANMIXTURED Function to handle scalar or multivariate real Gaussian
%                  mixture distributions.
%Implemented methods are: mean, cov, PDF, PDFS, PDFI, rand, randS.
%
%(UNCLASSIFIED) DISTRIBUTION STATEMENT A. Approved for public release.

methods(Static)
function meanVal=mean(w,mu)
%%MEAN Obtain the mean of the Gaussiax mixture distribution.
%
%INPUTS: w The 1XN or NX1 set of weights of the distribution components
%          such that all w>=0 and sum(w)=1.
%       mu The xDimXN set of mean values of the distribution
%          components.
%
%OUTPUTS: meanVal The xDimX1 mean of the Gaussian mixture distribution.
%
%It is simple to show that the mean of the distribution is the weighted
%sum of the means of the components.
%
%EXAMPLE:
%In this example,
% w=[0.4,0.6];
% mu=[1,-0.5;
%    -1, 1];
% Sigma=zeros(2,2,2);
% Sigma(:,:,1)=[4/9,  14/45;
%               14/45,4/9];
% Sigma(:,:,2)=[4/9, 0;
%                 0, 4/9];
% N=10000;
% analyticMean=GaussianMixtureD.mean(w,mu);
% sampMean=mean(GaussianMixtureD.rand(N,w,mu,Sigma),2);
% RelErr=abs((sampMean-analyticMean)./analyticMean)
%In this example, the relative error will typically be 3% or better.
%Increasing the number of samples will typically decrease the relative
%error, indicating that the analytic and sample means are converging.
%
%July 2019 David F. Crouse, Naval Research Laboratory, Washington D.C.

    meanVal=sum(bsxfun(@times,w(:).',mu),2);
end

function [covMat,meanVal]=cov(w,mu,Sigma)
%%COV Obtain the covariance matrix of the Gaussian mixture distribution
%     (the variance if it is scalar).
%
%INPUTS: w The 1XN or NX1 set of weights of the distribution components
%          such that all w>=0 and sum(w)=1.
%       mu The xDimXN set of mean values of the distribution
%          components.
%    Sigma The xDimXxDimXN set of covariance matrices of the
%          distribution components.
%
%OUTPUTS: covMat The xDimXxDim covariance matrix of the Gaussian
%                mixture.
%        meanVal The xDimX1 mean value of the distribution.
%
%The function calcMixtureMoments is called to find the covariance
%matrix.
%
%EXAMPLE:
% w=[0.4,0.6];
% mu=[1,-0.5;
%    -1, 1];
% Sigma=zeros(2,2,2);
% Sigma(:,:,1)=[4/9,  14/45;
%               14/45,4/9];
% Sigma(:,:,2)=[4/9, 0;
%                 0, 4/9];
% N=10000;
% analyticCov=GaussianMixtureD.cov(w,mu,Sigma);
% [~,sampleCov]=calcMixtureMoments(GaussianMixtureD.rand(N,w,mu,Sigma));
% RelErr=abs((sampleCov-analyticCov)./analyticCov)
%The relative error of each element will tend to be less than 1% and
%will decrease as the number of samples increases
%
%July 2019 David F. Crouse, Naval Research Laboratory, Washington D.C.

    [meanVal,covMat]=calcMixtureMoments(mu,w,Sigma);
end

function vals=PDF(z,w,mu,Sigma)
%%PDF Evaluate the PDF of a scalar or multivariate Gaussian
%     distribution at specified points.
%
%INPUTS: z The xDimXnumPoints set of points at which the PDF should be
%          evaluated.
%        w The 1XN or NX1 set of weights of the distribution components
%          such that all w>=0 and sum(w)=1.
%       mu The xDimXN set of mean values of the distribution
%          components.
%    Sigma The xDimXxDimXN set of covariance matrices of the
%          distribution components.
%
%OUTPUTS: vals The 1XN set of PDF values of the Gaussian mixture PDF
%              evaluated at the points in z.
%
%The PDF is the weighted sum of the PDFs of the individual components.
%
%EXAMPLE:
%Here, we validate the PDF by generating random samples and comparing
%the PDF plot with a histogram of the random samples.
% w=[0.4,0.6];
% mu=[1,-0.5];
% sigma=zeros(1,1,2);
% sigma(:,:,1)=4/9;
% sigma(:,:,2)=5/9;
% N=5000;
% figure(1)
% clf
% histogram(GaussianMixtureD.randS(N,w,mu,sigma),'Normalization','pdf','BinLimits',[-3,3])
% hold on
% numPoints=1000;
% x=linspace(-3,3,numPoints);
% vals=GaussianMixtureD.PDF(x,w,mu,sigma.^2);
% plot(x,vals,'linewidth',2)
% axis([-3,3,0,0.45])
% h1=xlabel('x');
% h2=ylabel('PDF(x)');
% set(gca,'FontSize',14,'FontWeight','bold','FontName','Times')
% set(h1,'FontSize',14,'FontWeight','bold','FontName','Times')
% set(h2,'FontSize',14,'FontWeight','bold','FontName','Times')
%
%July 2019 David F. Crouse, Naval Research Laboratory, Washington D.C.

    n=length(w);
    numMeas=size(z,2);

    vals=zeros(1,numMeas);
    for k=1:n
        vals=vals+w(k)*GaussianD.PDF(z,mu(:,k),Sigma(:,:,k));
    end
end

function vals=PDFS(z,w,mu,S)
%%PDF Evaluate the PDF of a scalar or multivariate Gaussian
%     distribution at specified points given the weights, means and
%     lower-triangular covariance matrices of the individual
%     components.
%
%INPUTS: z The xDimXnumPoints set of points at which the PDF should be
%          evaluated.
%        w The 1XN or NX1 set of weights of the distribution components
%          such that all w>=0 and sum(w)=1.
%       mu The xDimXN set of mean values of the distribution
%          components.
%        S The xDimXxDimXN set of lower-triangular square-root
%          covariance matrices of the distribution components.
%
%OUTPUTS: vals The 1XN set of PDF values of the Gaussian mixture PDF
%              evaluated at the points in z.
%
%The PDF is the weighted sum of the PDFs of the individual components.
%
%EXAMPLE:
%Here, we validate the PDF by generating random samples and comparing
%the PDF plot with a histogram of the random samples.
% w=[0.4,0.6];
% mu=[1,-0.5];
% sigma=zeros(1,1,2);
% sigma(:,:,1)=4/9;
% sigma(:,:,2)=5/9;
% N=5000;
% figure(1)
% clf
% histogram(GaussianMixtureD.randS(N,w,mu,sigma),'Normalization','pdf','BinLimits',[-3,3])
% hold on
% numPoints=1000;
% x=linspace(-3,3,numPoints);
% vals=GaussianMixtureD.PDFS(x,w,mu,sigma);
% plot(x,vals,'linewidth',2)
% axis([-3,3,0,0.45])
% h1=xlabel('x');
% h2=ylabel('PDF(x)');
% set(gca,'FontSize',14,'FontWeight','bold','FontName','Times')
% set(h1,'FontSize',14,'FontWeight','bold','FontName','Times')
% set(h2,'FontSize',14,'FontWeight','bold','FontName','Times')
%One will see that the histogram matches well with the plot.
%
%July 2019 David F. Crouse, Naval Research Laboratory, Washington D.C.

    n=length(w);
    numMeas=size(z,2);

    vals=zeros(1,numMeas);
    for k=1:n
        vals=vals+w(k)*GaussianD.PDFS(z,mu(:,k),S(:,:,k));
    end
end

function vals=PDFI(z,w,mu,SigmaInv,SigmaInvDet)
%%PDFI Evaluate the PDF of a scalar or multivariate Gaussian
%      distribution at specified points given the weights, means and
%      the inverses of the covariance matrices of the individual
%      components.
%
%INPUTS: z The xDimXnumPoints set of points at which the PDF should be
%          evaluated.
%        w The 1XN or NX1 set of weights of the distribution components
%          such that all w>=0 and sum(w)=1.
%       mu The xDimXN set of mean values of the distribution
%          components.
% SigmaInv The xDimXxDimXN set of inverse of covariance matrices of the
%          distribution components.
% SigmaInvDet Optionally, a length-N set of determinants of the matrices
%          in SigmaInv can be passed so as to speed up the computation. If
%          omitted or an empty matrix is passed, determinants will be taken
%          as needed.
%
%OUTPUTS: vals The 1XN set of PDF values of the Gaussian mixture PDF
%              evaluated at the points in z.
%
%The PDF is the weighted sum of the PDFs of the individual components.
%
%July 2019 David F. Crouse, Naval Research Laboratory, Washington D.C.

    n=length(w);
    numMeas=size(z,2);

    vals=zeros(1,numMeas);
    if(nargin<5||isempty(SigmaInvDet))
        for k=1:n
            vals=vals+w(k)*GaussianD.PDFI(z,mu(:,k),SigmaInv(:,:,k));
        end
    else
        for k=1:n
            vals=vals+w(k)*GaussianD.PDFI(z,mu(:,k),SigmaInv(:,:,k),SigmaInvDet(k));
        end
    end
end

function vals=CDF(z,w,mu,varVals)
%%CDF Evaluate the cumulative distribution function (CDF) of a scalar
%     Gaussian mixture at the specified points.
%
%INPUTS: z A matrix of the point(s) at which the CDF should be
%          evaluated.
%        w The 1XN or NX1 set of weights of the distribution components
%          such that all w>=0 and sum(w)=1.
%       mu The length N set of mean values of the distribution
%          components.
%  varVals The length-N set of variances of the distribution
%          components.
%
%OUTPUTS: val The scalar value(s) of the Gaussian mixture distribution
%             evaluated at the point(s) z.
%
%The CDF of the mixture is the weighted sum of the CDFs of the
%individual Gaussian components.
%
%July 2019 David F. Crouse, Naval Research Laboratory, Washington D.C.

    n=length(w);

    vals=zeros(size(z));
    for k=1:n
        vals=vals+w(k)*GaussianD.CDF(z,mu(k),varVals(k));
    end
end

function vals=rand(N,w,mu,P)
%%RAND Generate multivariate Gaussian mixture samples with given
%      weights, means and lower-triangular square root covariance
%      matrices.
%
%INPUTS: N The number of random variables to generate.
%        w The 1XN or NX1 set of weights of the distribution components
%          such that all w>=0 and sum(w)=1.
%       mu The xDimXN set of mean values of the distribution
%          components.
%        P The xDimXxDimXN set of covariance matrices of the
%          distribution components. If all of them are the same, then a
%          single xDimXxDim matrix can be passed.
%
%OUTPUT: x An xDimXN matrix of random instances of the multivariate
%          Gaussian distribution.
%
%The distribution is sampled by sampling the empirical distribution of
%weights, to determine which Gaussian component to sample, and then
%sampling the chosen Gausian distribution. Lower-triangular square
%roots of the componets must be found for sampling, so if these are
%available, then the randS method will be faster.
%
%July 2019 David F. Crouse, Naval Research Laboratory, Washington D.C.
%(UNCLASSIFIED) DISTRIBUTION STATEMENT A. Approved for public release.

    %The number of components of the mixture.
    n=length(w);
    xDim=size(mu,1);

    if(size(P,3)==1)
        P=repmat(P,[1,1,n]);
    end

    vals=zeros(xDim,N);

    %Determine which components are to be sampled.
    idx=EmpiricalD.rand([N,1],1:n,w);

    %Sample the Gaussian components that were selected.
    S=zeros(xDim,xDim,n);
    componentUsed=false(n,1);
    for curSamp=1:N
        curComp=idx(curSamp);
        %Get the lower-triangular square roots only for the components
        %that are used.
        if(componentUsed(curComp)==false)
            S(:,:,curComp)=chol(P(:,:,curComp),'lower');
            componentUsed(curComp)=true;
        end

        vals(:,curSamp)=mu(:,curComp)+S(:,:,curComp)*randn(xDim,1);
    end
end

function vals=randS(N,w,mu,S)
%%RANDS Generate multivariate Gaussian mixture samples with given
%       weights, means and lower-triangular square root covariance
%       matrices.
%
%INPUTS: N The number of random variables to generate.
%        w The 1XN or NX1 set of weights of the distribution components
%          such that all w>=0 and sum(w)=1.
%       mu The xDimXN set of mean values of the distribution
%          components.
%        S The xDimXxDimXN set of lower-triangular square-root
%          covariance matrices of the distribution components.
%
%OUTPUT: x An xDimXN matrix of random instances of the multivariate
%          Gaussian distribution.
%
%The distribution is sampled by sampling the empirical distribution of
%weights, to determine which Gaussian component to sample, and then
%sampling the chosen Gausian distribution.
%
%July 2019 David F. Crouse, Naval Research Laboratory, Washington D.C.
%(UNCLASSIFIED) DISTRIBUTION STATEMENT A. Approved for public release.

    %The number of components of the mixture.
    n=length(w);
    xDim=size(mu,1);

    if(size(S,3)==1)
        S=repmat(S,[1,1,n]);
    end

    vals=zeros(xDim,N);

    %Determine which components are to be sampled.
    idx=EmpiricalD.rand([N,1],1:n,w);

    %Sample the Gaussian components that were selected.
    for curSamp=1:N
        curComp=idx(curSamp);

        vals(:,curSamp)=mu(:,curComp)+S(:,:,curComp)*randn(xDim,1);
    end
end
end
end

%LICENSE:
%
%The source code is in the public domain and not licensed or under
%copyright. The information and software may be used freely by the public.
%As required by 17 U.S.C. 403, third parties producing copyrighted works
%consisting predominantly of the material produced by U.S. government
%agencies must provide notice with such work(s) identifying the U.S.
%Government material incorporated and stating that such material is not
%subject to copyright protection.
%
%Derived works shall not identify themselves in a manner that implies an
%endorsement by or an affiliation with the Naval Research Laboratory.
%
%RECIPIENT BEARS ALL RISK RELATING TO QUALITY AND PERFORMANCE OF THE
%SOFTWARE AND ANY RELATED MATERIALS, AND AGREES TO INDEMNIFY THE NAVAL
%RESEARCH LABORATORY FOR ALL THIRD-PARTY CLAIMS RESULTING FROM THE ACTIONS
%OF RECIPIENT IN THE USE OF THE SOFTWARE.

author: USNavalResearchLaboratory

repo: TrackerComponentLibrary


text:
function A = warmUpExercise()
%WARMUPEXERCISE Example function in octave
%   A = WARMUPEXERCISE() is an example function that returns the 5x5 identity matrix

A = eye(5);
% ============= YOUR CODE HERE ==============
% Instructions: Return the 5x5 identity matrix
%               In octave, we return values by defining which variables
%               represent the return values (at the top of the file)
%               and then set them accordingly.

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

end

author: AvaisP


text:
function [rwproc] = ranwalk(npoints, p)
% RANWALK generate and plot a trajectory of a random walk which
%    starts at 0, jumps up with probability p and down with
%    probability 1-p
%
% [rwproc] = ranwalk(npoints, p)
%
% Inputs: npoints - length of the trajectory
%         p - probability of the jump upwards
%
% Outputs: rwproc - trajectory of the process

% Authors: R.Gaigalas, I.Kaj
% v1.2 07-Oct-02

  % default parameter values
  if (nargin==0)
    npoints = 100;
    p = 0.5;
  end

  % generate jumps and sum up
  rwproc = [0 cumsum(2.*(rand(1, npoints-1)<p)-1)];

  stairs(rwproc);