function Film = bht_simulation(varargin)
% BHT_SIMULATION
%
% Return movie frames, captured at each time step, of the articulated
% bodies' motion.
%
% Syntax
%
% Film = BHT_SIMULATION('DataFilename',filename)
%
% Film = BHT_SIMULATION(...,options)
%
% Description
%
% This function does the following:
%
% o it opens a figure in which are ploted, at each time step, the
% fluid's boundaries and the articulated bodies at the position
% computed by the function bht_traject_compute with data of the
% DAT-File filename. The colors are those specified in the
% DAT-File.
% o it captures a movie frame.
%
% The function can be ran only after the articulated bodies' motion has
% been successfully computed by the function BHT_TRAJECT_COMPUTE, for
% it requires the MAT-File filename_results. You can specify the
% following options:
%
% * 'axis'
% - Expected value: array [xmin xmax ymin ymax].
% - Default value: Axis limits chosen by MATLAB at the initial
% simulation time.
% - Description: Set the limits for the x- and y-axis.
% * 'CenteredOnFish'
% - Expected value: integer.
% - Default value:None.
% - Description: Center the figure on the specified fish at each
% time step.
% * 'Grid'
% - Expected value: 'on' or 'off'
% - Default value: 'on'
% - Description: Display grid.
% * 'DrawCenterOfMass'
% - Expected value: LineSpec (line specification string).
% - Default value: None.
% - Description: Draw the trajectories of the articulated
% bodies' center of mass using properties' line of LineSpec.
% * 'DrawFishAxis'
% - Expected value: string 'axis_scale axis_color'
% (a number and a character corresponding to a MATLAB color).
% - Default value: None.
% - Description: Draw the moving frames attached to the articulated
% bodies with scale axis_scale and color axis_color.
% * 'TimeRange'
% - Expected value: array [T1 T2].
% - Default value: Time range of the DAT-File.
% - Description: Time range for the animation. If the time range is
% not within the time range of the DAT-File, the animation is
% empty.
% * 'DisplayTime'
% - Expected value: 'on' or 'off'
% - Default value: 'off'
% - Description: Display a clock in the top-right corner.
%
% * 'DrawFishVelocity'
% - Expected value: string 'velocity_scale velocity_color'
% (a number and a character corresponding to a MATLAB color).
% - Default value: None.
% - Description: Draw bodies' velocities vectors at the centers of
% mass of the bodies with scale velocity_scale and color
% velocity_color.
%
% The command line:
%
% MOVIE(Film,options)
%
% allows to replay the movie. Use the built-in MATLAB function
% movie2avi to turn MATLAB movies into avi movies.
%
% Examples
%
% The command line:
%
% Film =
% bht_simulation('DataFilename','two_fishes','DrawCenterOfMass','--b','
% axis',[-10 5 -3 5],'TimeRange',[5 7],
% DrawFishAxis,{4,'-r','LineWidth',2},'DisplayTime','on',
% 'CenteredonFish',2,'Grid','off');
%
% returns movie frames, captured at each time step, of the motion of the
% articulated bodies described in the DAT-File two_fishes. The
% trajectories of the articulated bodies' center of mass are ploted
% with LineSpec = '--b' (dashed blue line). The x-axis limits are -10
% and 5 and the y-axis limits are -3 and 5, centered on fish 2. The
% simulation goes from time = 5 to time = 7. The moving frames attached
% to the bodies are displayed with scale = 4, color = r (red) and
% LineWidth = 2. Eventually, a clock is displayed in the top-right
% corner of the pictures.
%
% See also BHT_TRAJECT_COMPUTE
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% %
% BIOHYDRODYNAMICS MATLAB TOOLBOX %
% %
% A. MUNNIER and B. PINCON %
% %
% alexandre.munnier@iecn.u-nancy.fr bruno.pincon@iecn.u-nancy.fr %
% http://www.iecn.u-nancy.fr/~munnier http://www.iecn.u-nancy.fr/~pincon %
% %
% %
% INSTITUT ELIE CARTAN, NANCY 1 %
% http://www.iecn.u-nancy.fr/ %
% %
% INRIA Lorraine, Projet CORIDA %
% http://www.iecn.u-nancy.fr/~corida/ %
% %
% %
% %
% %
% August 15th 2008 %
% %
% GNU GPL v3 licence %
% %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%##########################################################################
% Preallocating output variable
Film = struct('cdata',{},'colormap',{});
%==========================================================================
% reading options
%-----------------------
% if there is no options
% ==================================
% default values for options
noaxis = 1;
time_range = [];
draw_param = 0;
draw_axis = 0;
display_time = 0;
gridon = 1;
centered_on_fish = 0;
draw_velocity = 0;
no_filename = 1;
%==========================================================================
% reading options
%================
nb_varargin = numel(varargin);
kk = 1;
while kk <= nb_varargin
if ~ischar(varargin{kk})
error(['Unknown option ',num2str(varargin{kk})]) %#ok<ST2NM>
end
%------------------------------------------------------------------
switch lower(varargin{kk})
%--------------------------------------------------------------
case 'datafilename'
if ischar(varargin{kk+1})
filename = varargin{kk+1};
exten = regexp(filename,'\w*\.(\w*)', 'tokens');
no_filename = 0;
if isempty(exten)
filename = [filename,'.dat']; %#ok<AGROW>
end
else
error('Unexpected expression for option DataFilename.')
end
kk = kk+2;
%--------------------------------------------------------------
case 'axis'
if isvector(varargin{kk+1})
if length(varargin{kk+1}) == 4
axis_param0 = varargin{kk+1};
noaxis = 0;
else
error('Unexpected value for option Axis. Must be a vector [Xmin Xmax Ymin Ymax].')
end
else
error('Unexpected value for option Axis. Must be a vector [Xmin Xmax Ymin Ymax].')
end
kk = kk+2;
%-------------
case 'drawcenterofmass'
draw_data = varargin{kk+1};
draw_param = 1;
kk = kk+2;
%-------------
case 'drawfishaxis'
draw_axis_param = varargin{kk+1};
draw_axis = 1;
kk = kk+2;
%-------------
case 'drawfishvelocity'
draw_velocity_param = varargin{kk+1};
draw_velocity = 1;
kk = kk+2;
%----------------------------------------------------------
case 'timerange'
if isvector(varargin{kk+1})
if length(varargin{kk+1}) == 2
time_range = varargin{kk+1};
else
error('Unexpected value for option TimeRange. Must be a vector [Tmin Tmax].')
end
else
error('Unexpected value for option TimeRange. Must be a vector [Tmin Tmax].')
end
kk = kk+2;
%----------------------------------------------------------
case 'centeredonfish'
if isfloat(varargin{kk+1})
centered_on_fish = varargin{kk+1};
else
error('Unexpected value for option CenterdOnFish. Must be an integer.')
end
kk = kk+2;
%-------------
case 'grid'
if ischar(varargin{kk+1})
if strcmpi(varargin{kk+1},'off')
gridon = 0;
elseif strcmpi(varargin{kk+1},'on')
gridon = 1;
else
error('Unexpected value for option Grid. Must be ''on'' or ''off''.')
end
else
error('Unexpected value for option Grid. Must be ''on'' or ''off''.')
end
kk = kk+2;
%-------------
case 'displaytime'
if ischar(varargin{kk+1})
if strcmpi(varargin{kk+1},'off')
display_time = 0;
elseif strcmpi(varargin{kk+1},'on')
display_time = 1;
else
error('Unexpected value for option DisplayTime. Must be ''on'' or ''off''.')
end
else
error('Unexpected value for option DisplayTime. Must be ''on'' or ''off''.')
end
kk = kk+2;
otherwise
error(['Unknown option ''',varargin{kk},''''])
end
end
% =========================================================================
% DataFilename is missing
if no_filename
filename = input('Enter data file''s name (hit Return to browse): ','s');
if isempty(filename);
filename = uigetfile('*.dat','Select data file');
end
end
%==========================================================================
base_filename = regexp(filename,'(\w*)(?:\.)?', 'tokens');
base_filename = base_filename{1}{1};
filename1 = [base_filename,'_data'];
filename2 = [base_filename,'_results'];
% need datacompil =========================================================
load(filename1,'gene_data','discret_data','bound_data','ode_data',...
'phys_data','fish');
load(filename2,'cont_filename','cont_parameters','T','p','dp','q','dq','Q','dQ');
% =========================================================================
nb_boundaries = gene_data.nb_boundaries;
nb_solids = gene_data.nb_solids;
nb_fishes = gene_data.nb_fishes;
% =========================================================================
if ~isempty(time_range)
ind_start = find(T<=time_range(1),1,'last'); %#ok<NODEF>
ind_end = find(T>=time_range(2),1,'first');
else
ind_start = 1;
ind_end = length(T);
end
% =========================================================================
bht_simul_fig = figure(1);
num_frames = 1;
%---------------
for compt = ind_start:ind_end
t = T(compt);
% draw fluid's boundaries =============================================
clf;
hold on;
for j=nb_solids+1:nb_boundaries
fill(bound_data(j).X(1,:),bound_data(j).X(2,:),bound_data(j).color)
end
% draw fishes =========================================================
for j = 1:nb_solids
h = [q(3*(j-1)+1,compt) ; q(3*(j-1)+2,compt)];
theta = q(3*(j-1)+3,compt);
R = [cos(theta),-sin(theta);...
sin(theta), cos(theta)];
Y = R*bound_data(j).X + h*ones(1,bound_data(j).nbr_points);
fill(Y(1,:),Y(2,:),bound_data(j).color)
end
% draw trajectories of fishes =========================================
if draw_param
if iscell(draw_data)
for k = 1:nb_fishes
plot(Q(3*(k-1)+1,:),Q(3*(k-1)+2,:),draw_data{:})%#ok
end
else
for k = 1:nb_fishes
plot(Q(3*(k-1)+1,:),Q(3*(k-1)+2,:),draw_data)%#ok
end
end
end
% draw fishes' frames =================================================
if draw_axis
if iscell(draw_axis_param)
for k = 1:nb_fishes
theta = Q(3*k,compt);
R = [cos(theta),-sin(theta);...
sin(theta), cos(theta)];
% e1 = axis_scale*R*[1;0];
% e2 = axis_scale*R*[0;1];
e1 = R*[1;0];
e2 = R*[0;1];
quiver([Q(3*(k-1)+1,compt),Q(3*(k-1)+1,compt)],...
[Q(3*(k-1)+2,compt),Q(3*(k-1)+2,compt)],...
[e1(1),e2(1)],...
[e1(2),e2(2)],draw_axis_param{:});
end
else
for k = 1:nb_fishes
theta = Q(3*k,compt);
R = [cos(theta),-sin(theta);...
sin(theta), cos(theta)];
% e1 = axis_scale*R*[1;0];
% e2 = axis_scale*R*[0;1];
e1 = R*[1;0];
e2 = R*[0;1];
quiver([Q(3*(k-1)+1,compt),Q(3*(k-1)+1,compt)],...
[Q(3*(k-1)+2,compt),Q(3*(k-1)+2,compt)],...
[e1(1),e2(1)],...
[e1(2),e2(2)],draw_axis_param);
end
end
end
% draw fishes' velocity ================================================
if draw_velocity
if iscell(draw_velocity_param)
for k = 1:nb_fishes
velocity = [dQ(3*(k-1)+1,compt),dQ(3*(k-1)+2,compt)];
quiver(Q(3*(k-1)+1,compt),Q(3*(k-1)+2,compt),...
velocity(1),velocity(2),draw_velocity_param{:});
end
else
for k = 1:nb_fishes
velocity = [dQ(3*(k-1)+1,compt),dQ(3*(k-1)+2,compt)];
quiver(Q(3*(k-1)+1,compt),Q(3*(k-1)+2,compt),...
velocity(1),velocity(2),draw_velocity_param);
end
end
end
% center of mass ======================================================
if draw_param
for k = 1:nb_fishes
plot(Q(3*(k-1)+1,compt),Q(3*(k-1)+2,compt),'Xk','linewidth',2)
end
end
% grid on =============================================================
if gridon
set(gca,'Layer','bottom')
grid on
end
% set axis ============================================================
if compt == ind_start && noaxis
axis_param0 = axis;
end
%----------------------------------------------------------------------
if centered_on_fish
axis_param = axis_param0+(Q(3*(centered_on_fish-1)+1,compt)-...
Q(3*(centered_on_fish-1)+1,1))*[1 1 0 0]...
+(Q(3*(centered_on_fish-1)+2,compt)-...
Q(3*(centered_on_fish-1)+2,1))*[0 0 1 1];
else
axis_param = axis_param0;
end
%----------------------------------------------------------------------
axis equal
axis(axis_param)
set(gca,'units','pixels');
%----------------------------------------------------------------------
%title
title(['Simulation from data file: ',base_filename],'interpreter','none')
%----------------------------------------------------------------------
% displaytime
if display_time
text(0.95*(axis_param(2)-axis_param(1))+axis_param(1),...
0.95*(axis_param(4)-axis_param(3))+axis_param(3),['t = ',num2str(t,'%6.2f')],...
'HorizontalAlignment','right','verticalalignment','top',...
'fontname','arial','BackgroundColor',[.7 .9 .7])
end
drawnow
%----------------------------------------------------------------------
% ensure that all of the captured frames have the same dimension
if num_frames == 1
axis_position = get(gca,'position');
end
%----------------------------------------------------------------------
% capturing frame
%++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Film(num_frames) = getframe(bht_simul_fig,axis_position);
%++++++++++++++++++++++++++++++++++++++++++++++++++++++++
num_frames = num_frames +1;
end
