%	TMTF.m
%	Demonstration of temporal contrast sensitivity function.
%	
%	History
%	LLK	4/??/98	Wrote it.
%	CCC	12/30/98	Animation output, Write out option.
%	CCC	1/12/99	Add Comments, user friendlier interface, cosmetics 

clear
imageOut=1;	%0: Disable image write out; 1: allow image write out.
ntimes=100; %The number of loops of the movie.
refreshRate=75; % Assigned monitor refresh rate. Can be changed 
					 % the users for best results

% making a linear look-up table
% The default power of inverse Gamma function is set to 0.4
x=(1:256)';
p=0.4;
base=23;
% Inver gamma function
invLut=base+(255-base)*(x.^p)/(256^p);
% making lookup table
map=invLut(1:1:256)*[1 1 1]/255;
[lvl dummy]=size(map);

% number of frames per period
% temp. frequency = (frame/s)/(frame/period)
prd=[60 40 30 24 20 15 12 10 8 6 5 4 3 2];
np=length(prd);

% range of contrasts
nc=14;
cmin=.005;
c=exp((log(cmin)-.0001):(-log(cmin)/(nc-1)):0)'*ones(1,np);

% matrix of periods
prd=ones(nc,1)*prd;

% random phase matrix
phase=rand(nc,np);

% initializing movie
nfr=120;
set(figure(1),'Position',[400 300 250 170])
m=moviein(nfr);
colormap(map);

% making image
for n=1:nfr
   dispWord=['Prepare frame ' num2str(n) '/120'];
   disp(dispWord);
   mat=floor(lvl/2+(lvl/2)*c.*cos(2*pi*(phase+n./prd)));
   % Image write out
   if (imageOut==1)
   	if n <10
   		writeFilename=['tmtf00',num2str(n),'.bmp'];
   	elseif n<100
   		writeFilename=['tmtf0',num2str(n),'.bmp'];
		else
   		writeFilename=['tmtf',num2str(n),'.bmp']; 
  	 	end
   	imwrite(uint8(round(kron(mat,ones(1,1)))),map,writeFilename);
   end   
   % End Image write out
   image(mat)	
	m(:,n)=getframe;
end
disp('Show movie');
disp('Type ctrl-C to stop')
% play movie colormap(map);
movie(m,ntimes,refreshRate)