% drumplots - return functions to do plots for one array of drum tracking data
% Execute each function in its own figure to see the results.
function PX=drumplots(BDT,varargin)

	% columns: time, type:{1,2}, error, beat pos, prob, KL, post entro, prior entro, beat period
	opts=prefs('bar',4,'with_metre',0,'marker_size',2,varargin{:});
	K=find(BDT(:,2)==1); % kick
	S=find(BDT(:,2)==2); % snare
	TM=mod(BDT(:,4)+0.5,opts.bar)-0.5;
	Surp=max(0,-log(BDT(:,5)));
	Info=BDT(:,6);
	HPrior=BDT(:,8);
	HPost=BDT(:,7);
	plotargs={'MarkerSize',opts.marker_size};

	PX={	@()scat_with_metre([HPrior,Surp], 'surprise/nats', 'prior entropy/nats'), ...
			@()scat_with_metre([HPost,Surp],'posterior entropy/nats','surprise/nats'), ...
			@()scat_with_metre([Surp,Info], 'surprise/nats','information gain/nats'), ...
			@()scat_with_metre([Info,HPost],'information gain/nats','posterior entropy/nats'), ...
			@()scat_with_metre([HPrior,Info],'prior entropy/nats','information gain/nats'), ...
			@()scat_with_metre([HPrior,HPost],'posterior entropy/nats', 'prior entropy/nats'), ...
			@()scat_by_type([TM,Info], {'metrical position','information gain/nats'}), ...
			@()scat_by_type([TM,Surp], {'metrical position','surprise/nats'}), ...
			@()scat_by_type([TM,noisy(BDT(:,2),0.02)], {'metrical position','drum type'}), ...
			@metre_entropies, ...
			@()metre_hist, ...
			@()total_info_dist, ...
			@()mean_info_dist, ...
			@()scat3d_by_type([HPrior,HPost,Info],{'prior entro','posterior entro','info gain'}), ...
	%	@()scat_by_type([HPrior,HPost,Info], { 'prior entropy/nats', 'posterior entropy/nats', 'information gain/nats'});
		};


	function y=noisy(x,k), y=x+k*randn(size(x)); end

	function h=metre_hist
		DTime=round(8*TM);
		T1=min(DTime);
		T2=max(DTime);
		ID=accumhist(DTime-T1+1,1,T2-T1+1);
		colormap jet;
		bar((T1:T2)/8,ID);
		xlabel('metrical position');
		title('total event count');
		h=gca;
	end

	function h=total_info_dist
		DTime=round(8*TM);
		T1=min(DTime);
		T2=max(DTime);
		ID=accumhist(DTime-T1+1,Info,T2-T1+1);
		colormap jet;
		bar((T1:T2)/8,ID);
		xlabel('metrical position');
		title('total information gain');
		h=gca;
	end

	function h=mean_info_dist
		DTime=round(8*TM);
		T1=min(DTime);
		T2=max(DTime);
		ID=accumhist(DTime-T1+1,Info,T2-T1+1);
		CC=accumhist(DTime-T1+1,1,T2-T1+1);
		colormap jet;
		bar((T1:T2)/8,ID./CC);
		xlabel('metrical position');
		title('mean information per event');
		h=gca;
	end
	function h=metre_entropies
		newplot; hold on;
		scat1([TM,BDT(:,8)],'b+');
		scat1([TM,BDT(:,7)],'ro');
		hold off;
		legend({'prior entropy','posterior entropy'});
		xlabel('metrical position');
		ylabel('entropy/nats');
		rotate3d off;
		h=gca;
	end

	function h=scat_with_metre(X,lab1,lab2)
		if opts.with_metre
			h=scat_by_type([X,TM],{lab1,lab2,'metrical position'});
			camproj('orthographic');
		else
			h=scat_by_type(X,{lab1,lab2});
		end
	end

	function h=scat3d_by_type(X,labels)
		h=scat_by_type(X,labels);
		camproj('perspective');
		axis vis3d;
	end

	function h=scat_by_type(X,labels)
		newplot; hold on;
		scat1(X(K,:),'ro');
		scat1(X(S,:),'b+');
		hold off;
		legend({'kick','snare'});
		axislabels(labels);
		if length(labels)>2, rotate3d on; 
		else rotate3d off; 
		end
		h=gca;
	end

	function scat1(X,marker)
		if size(X,2)<=2
			plot(X(:,1),X(:,2),marker,plotargs{:});
		else 
			plot3(X(:,1),X(:,2),X(:,3),marker,plotargs{:});
		end
		box on;
	end
end

