function ErrCode=DesignF(FrameFile,TotIt,MaxDays)% DesignF   Design the frame and store it in FrameFile, based on training data
%           stored in mat-files, DataXnnn.mat, in current directory.
% The function will return after a total of TotIt iterations.
% 
% ErrCode=DesignF(FrameFile,TotIt,MaxDays);%-----------------------------------------------------------------------------------
% arguments:
%  ErrCode   - 0 is returned if the function execute without error
%              1 an 'unexpected' error occurred, lasterr may explain it
%              2 or larger, another error occurred, 
%  FrameFile - the name of the mat-file used to store the frame
%  TotIt     - total number of iterations to do 
%              If more than TotIt iterations already has been done on this frame
%              then no iterations will be done now, the function just return and
%              do not change the frame in FrameFile.
%  MaxDays   - how many days the design process should be allowed to run
%              if MaxDays=0.25 it is 6 hours, 0.01 is 14 minutes and 24 seconds
%-----------------------------------------------------------------------------------

%----------------------------------------------------------------------
% Copyright (c) 2000.  Karl Skretting.  All rights reserved.
% Hogskolen in Stavanger (Stavanger University), Signal Processing Group
% Mail:  karl.skretting@tn.his.no   Homepage:  http://www.ux.his.no/~karlsk/
% 
% HISTORY:  dd.mm.yyyy
% Ver. 1.0  11.12.2000  KS: function made
% Ver. 1.1  04.10.2001  KS: changed how vector selection is used
%   also the new function NormalizeF is used now
% Ver. 1.2  05.11.2001  KS: Overlapping frame for 2D added, (Mdim==2) & (P>1)
% Ver. 1.3  06.11.2001  KS: Separable frame for 2D added
% Ver. 1.4  03.12.2002  KS: moved from ..\Frames to ..\FrameTools%----------------------------------------------------------------------

Mfile='DesignF';
ErrCode=0; 
Display=1;
UpdateFforEachSignal=0;

if (nargin<3)
   disp([Mfile,': wrong number of arguments, see help.']);
   ErrCode=2;
   return
end
StopTime=now+MaxDays;

if ~exist([FrameFile,'.mat'])
   disp([Mfile,': the FrameFile ',FrameFile,'.mat does not exist.']);
   ErrCode=3;
   return
end
try
   % load the following variables from FrameFile: Class, Type, Mdim, F, SizeF
   % G, Ctab, Dtab, Fbest, Savg, Mdat, PreProc, VecSel, InitialF, History, SNRtot
   load(FrameFile);
catch
   disp([Mfile,': the FrameFile ',FrameFile,'.mat could not be read.']);
   ErrCode=3;
   return
end
History0=char(History,[Mfile,': Training of ',FrameFile,' started ',datestr(now,0)]);

% The data-files used for training of the frame should be ok and
% stored in the correct mat-files.
% If not, MakeX may be called at this point
% but for now we assume that it is ok.

% check and display the frame type
if Type=='b'
   disp([Mfile,': ',FrameFile,' is a block oriented frame.']);
elseif Type=='o'
   disp([Mfile,': ',FrameFile,' is an overlapping frame.']);
elseif Type=='g'
   disp([Mfile,': ',FrameFile,' is a general frame.']);
elseif Type=='s'
   disp([Mfile,': ',FrameFile,' is a separable frame.']);
elseif Type=='t'
   disp([Mfile,': ',FrameFile,' is a tree structured filter bank (frame), ',...
         Mfile,' is not ready for this yet.']);
   ErrCode=5;
   return
else
   disp([Mfile,': ',FrameFile,' is type ',Type,' which is not allowed.']);
   ErrCode=5;
   return
end
if (length(VecSel.arg4)==0); VecSel.arg4=0; end;
   
% now do the main iterations
% try
   K=SizeF(Mdim+1);
   N=prod(SizeF(1:Mdim));
   P=prod(SizeF((Mdim+2):(2*Mdim+1)));
   if (Mdim==2) & (P>1) 
      Ic=(prod(2*SizeF((Mdim+2):(2*Mdim+1))-1)+1)/2;  % Ic==P when Mdim==1
      %TotWWT=zeros(K,K,Ic);
      TotWWT=zeros(K*P,K*P);    % larger than what is strictly needed
      TotXWT=zeros(N,K*P);
   else
      Ic=(prod(2*SizeF((Mdim+2):(2*Mdim+1))-1)+1)/2;  % Ic==P when Mdim==1
      TotWWT=zeros(K,K,Ic);
      TotXWT=zeros(N,K,P);
   end
   ItNo=size(SNRtot,2);
   Fnormfact=ones(K,1);
   if ItNo>0
      SNRbest=max(SNRtot(Mdat+1,:));
   else
      SNRbest=0;
   end
   if Display
      disp([Mfile,': Training started ',datestr(now,0),', ',int2str(TotIt-ItNo),' it.']);
      t3=num2str(SNRbest,4);t3=[blanks(8-length(t3)),t3];
      disp(['Best SNR is now ',t3]);
   end
   %
   for i=1:(TotIt-ItNo)
      if (now>StopTime); break; end;
      if ~UpdateFforEachSignal
         if (Mdim==2) & (P>1)
            TotWWT=zeros(K*P,K*P);
            TotXWT=zeros(N,K*P);
         else
            TotWWT=zeros(K,K,Ic);
            TotXWT=zeros(N,K,P);
         end
      end
      ex=0;er=0;      for m=1:Mdat
         if (Mdat>1) | (i==1) | (Type=='s')
            % read the data file
            t1=int2str(m);t1=['0000',t1];t1=t1((length(t1)-2):length(t1));
            DataFile=['DataX',t1]; 
            if ~exist([DataFile,'.mat'])
               disp([Mfile,': the DataFile ',DataFile,'.mat does not exist.']);
               ErrCode=3;
               return
            end
            clear S X W    % make sure that variables are from the DataFile
            load(DataFile);    
            [n,L]=size(X);     % we need the L variable
         end
         disp([Mfile,': use DataFile ',DataFile,'.mat, iteration ',int2str(i),'.']);
         if UpdateFforEachSignal & (i>1)
            TotWWT=TotWWT-WWT;
            TotXWT=TotXWT-XWT;
         end
         % do vector selection
         if Type=='g'
            Fg=BuildFg(F,G);
            Q=size(F,1);
         elseif Type=='s'
            Fv=F{1};
            Fh=F{2};
            Fg=kron(Fh,Fv);
         else
            Fg=F;
         end
         % Fg is now used as the frame, do vector selection
         if strcmp(VecSel.Prog1,'BlockVS')
            disp([Mfile,': vector selection using ',VecSel.Prog1,' is obsolete.']);
            ErrCode=4;
            return
            % koden er i slutten av fila
         elseif strcmp(VecSel.Prog1,'VSblock')
            % first find S (used as third argument in VSblock)
            if (i==1); S=Savg; end;
            % when i>1, S is stored in DataFile (except for 'VSab2' which store W)
            if strcmp(VecSel.arg1,'VSab2') & (i>1)
               W=(Fnormfact(:)*ones(1,L)).*W;   % W was stored instead of S
               S=W;              % and also used as third argument in VSblock
            end
            %if strcmp(VecSel.arg1,'VSab2') & (i==1)
            %   [W,S]=GlobalMP(X,Fg,Savg);    % first time GlobalMP
            %else
            %   disp([Mfile,': do vector selection using ',VecSel.Prog1,'.']);
            %                      VSalg           B          el          Bf
            W=VSblock(X,Fg,S,VecSel.arg1,VecSel.arg2,VecSel.arg3,VecSel.arg4);   
            %end
            S=full(sum(W~=0));
         elseif strcmp(VecSel.Prog1,'VSolap1')
            % first find third argument in VSolap1
            if (i==1)
               W=Savg; 
            else
               W=(Fnormfact*ones(1,L)).*W;   % W was stored instead of S
            end
            %
            W=VSolap1(X,Fg,W,VecSel.arg1,VecSel.arg2);   
            % S=full(sum(W~=0));  % is not used ??
         elseif strcmp(VecSel.Prog1,'VSolap2')
            % first find third argument in VSolap2
            if (i==1)
               W=Savg; 
            else
               W=(Fnormfact*ones(1,L)).*W;   % W was stored instead of S
            end
            %
            W=VSolap2(X,Fg,W,VecSel.arg1,VecSel.arg2);   
            % S=full(sum(W~=0));  % is not used ??
         elseif strcmp(VecSel.Prog1,'GlobalMP')
            disp([Mfile,': vector selection using ',VecSel.Prog1,' is obsolete.']);
            ErrCode=4;
            return
            % [W,S]=GlobalMP(X,Fg,Savg);  
         else
            disp([Mfile,': not ready for vector selection using ',VecSel.Prog1,'.']);
            ErrCode=4;
            return
         end
         % calculate WWT and XWT
         if Type=='s'
            % do not use WWT or XWT, the weights W are stored directly
            WWT=[];XWT=[];
         elseif (P==1)               % for all values of Mdim !!
            WWT=full(W*W');   
            XWT=full(X*W');   
         else 
            if Mdim==1      
               temp=W;
               for p=1:P
                  WWT(:,:,p)=full(W*temp');   
                  XWT(:,:,p)=full(X*temp');   
                  temp=[temp(:,L),temp(:,1:(L-1))];
               end
            elseif (Mdim==2) & (P==4) & (L==4096)
               % WWT (PKxPK) should be Wstar (PKxL) * Wstar' (LxPK)
               % XWT (NxPK) should be X (NxL) * Wstar' (LxPK)    
               % W (KxL) 
               L1=sqrt(L);L2=sqrt(L);  % ==64
               Wstar=sparse(P*K,L);
               for p1=0:1
                  for p2=0:1
                     p=p1+p2*2;
                     Wt=reshape(full(W),K,L1,L2);    % KxL1xL2
                     if p1
                        Wt1=Wt;            
                        Wt(:,1,:)=Wt1(:,L1,:);
                        Wt(:,2:L1,:)=Wt1(:,1:(L1-1),:);
                     end
                     if p2
                        Wt1=Wt;            
                        Wt(:,:,1)=Wt1(:,:,L2);
                        Wt(:,:,2:L2)=Wt1(:,:,1:(L2-1));
                     end
                     Wt=reshape(Wt,K,L);
                     Wstar((1:K)+p*K,:)=sparse(Wt);
                  end
               end
               clear Wt Wt1
               WWT=full(Wstar*Wstar');
               XWT=full(X*Wstar');
               clear Wstar
            else
               disp([Mfile,': overlapping frame and 2D or 3D signal, not ready yet.']);
               ErrCode=4;
               return
            end
         end
         %
         if strcmp(VecSel.arg1,'VSab2') | strcmp(VecSel.Prog1,'VSolap1') ...
               | strcmp(VecSel.Prog1,'VSolap2') | Type=='s'
            %disp([Mfile,': save ',DataFile,', here W is saved.']);
            save(DataFile,'Name','X_DC','X_LP','X','SizeX','ss2','W','WWT','XWT');
         else
            %disp([Mfile,': save ',DataFile,', here S is saved.']);
            save(DataFile,'Name','X_DC','X_LP','X','SizeX','ss2','S','WWT','XWT');
         end
         %
         ex=ex+ss2;    
         if P==1
            R=X-full(Fg*W);
         elseif Mdim==1
            temp=W;R=X;
            for p=1:P
               R=R-full(Fg(:,:,p)*temp);   
               temp=[temp(:,L),temp(:,1:(L-1))];
            end
         elseif (Mdim==2) & (P==4) & (L==4096)
            % this is the special 2D overlapping case handled
            R=X;
            L1=sqrt(L);L2=sqrt(L);  % ==64
            for p1=0:1
               for p2=0:1
                  p=p1+p2*2;
                  Wt=reshape(full(W),K,L1,L2);  % 3D is not fully handled as sparse
                  if p1
                     Wt1=Wt; 
                     Wt(:,1,:)=Wt1(:,L1,:);
                     Wt(:,2:L1,:)=Wt1(:,1:(L1-1),:);
                  end
                  if p2
                     Wt1=Wt; 
                     Wt(:,:,1)=Wt1(:,:,L2);
                     Wt(:,:,2:L2)=Wt1(:,:,1:(L2-1));
                  end
                  Wt=reshape(Wt,K,L);
                  R=R-Fg(:,:,p+1)*Wt;