function [MiCsMat1,MiCsMat2,XCorrMat,IndVec1,IndVec2]=fMriXCorrStat(FuncMat,Base,NShotPerCycle,MinLag,MaxLag)

%   fMriXCorrStat.m
%
%   This program is designed to calculate xcorr of fMRI signal 

%   This software is protected by german copyright and international treaties.             
%   Copyright 2004 Qun Zhao, Markus Junghöfer & Peter Peyk. All Rights Reserved.                     
%                                                                                          
%   THIS SOFTWARE AND DOCUMENTATION IS PROVIDED "AS IS," AND COPYRIGHT HOLDERS MAKE        
%   NO REPRESENTATIONS OR WARRANTIES, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO,    
%   WARRANTIES OF MERCHANTABILITY OR FITNESS FOR ANY PARTICULAR PURPOSE OR THAT THE USE    
%   OF THE SOFTWARE OR DOCUMENTATION WILL NOT INFRINGE ANY THIRD PARTY PATENTS, COPYRIGHTS,
%   TRADEMARKS OR OTHER RIGHTS. COPYRIGHT HOLDERS WILL NOT BE LIABLE FOR ANY DIRECT,       
%   INDIRECT, SPECIAL OR CONSEQUENTIAL DAMAGES ARISING OUT OF ANY USE OF THE SOFTWARE      
%   OR DOCUMENTATION.

 if nargin<5; MaxLag=[]; end
 if nargin<4; MinLag=[]; end
 if nargin<3; NShotPerCycle=[]; end
 if nargin<2; Base=[]; end
 if nargin<1; FuncMat=[]; end

 if isempty(Base); return; end
 if isempty(FuncMat); return; end
 
[XSizeFunc,YSizeFunc,NSlice,NShot]=size(FuncMat);
FigTitle='XCorr Shot Order';
[IndVec1,IndVec2,Num,LastShot]=CalcInd12Vec(NShot,NShotPerCycle,Base,FigTitle);

IndVec=-1.*ones(NShot,1);
IndVec(IndVec2)=1;
StdIndVec=std(IndVec);
AbsMax=0;
MiCsMat1=zeros(XSizeFunc,YSizeFunc,NSlice);
MiCsMat2=MiCsMat1;
Message=char(zeros(9,100));
hProgressMsgbox=msgbox(Message,'Progress Info:','help'); 
hChildProgressMsgbox=get(hProgressMsgbox,'Children');
set(hChildProgressMsgbox(length(hChildProgressMsgbox)),'Visible','off');
AbsMaxLag=max(abs([MinLag MaxLag]));
NXCorrInd=2.*AbsMaxLag+1;
NXCorrInd2=AbsMaxLag+1;
XCorrMat=zeros(XSizeFunc,YSizeFunc,NSlice,NXCorrInd);
MinInd=AbsMaxLag+1+MinLag;
MaxInd=AbsMaxLag+1+MaxLag;
for SliceInd=1:NSlice
	fprintf('Slice %g of %g ...\n',SliceInd,NSlice);
	Message=char('Calculate Cross-Correlation !');
	Message=char(Message,'');
	Message=char(Message,['Base:  ',int2str(Base)]);
	Message=char(Message,['Min. Lag:  ',int2str(MinLag)]);
	Message=char(Message,['Max. Lag:  ',int2str(MaxLag)]);
	Message=char(Message,'');
	Message=char(Message,['Slice ' int2str(SliceInd) ' of ' int2str(NSlice) ' ...']);
	set(hChildProgressMsgbox(2),'String',Message);
	pause(.001);
	for i=1:XSizeFunc
		for j=1:YSizeFunc
			X1=squeeze(FuncMat(i,j,SliceInd,:));
			StdX1=std(X1);
% 			X1=X1-mean(X1);
% 			Tmp=xcorr(X1,IndVec,MaxLag,'unbiased');
			if StdX1>0
				Tmp=xcov(X1,IndVec,MaxLag,'unbiased');	
				Tmp=Tmp./(StdX1.*StdIndVec);
				XCorrMat(i,j,SliceInd,:)=Tmp;
				[MiCsMat1(i,j,SliceInd),MiCsMat2(i,j,SliceInd)]=max(Tmp(MinInd:MaxInd));
			end
		end
	end
end
MiCsMat2=MiCsMat2-NXCorrInd2.*ones(size(MiCsMat2));
close(hProgressMsgbox);
return;
