function warray(fid,var,unit,ftfmt,txt,norec)
%function warray(fid,var,ftfmt,txt [,norec]) --- write out a 1D or 2D MT3D matrix
%general function to write out nD matries. See UD2REL in Modflow manual
% No limit set to number of columns
% ftfmt is a fortran format (10E3.2) (4I1) etc  10G6.1 also works
% TO 070703 081228 090907


% Copyright 2009 Theo Olsthoorn, TU-Delft and Waternet, without any warranty
% under free software foundation GNU license version 3 or later

[mlfmt,N,ftfmt]=ftfmt2ml(ftfmt); if mlfmt(end)=='i', mlfmt(end)='d'; end

if nargin<5 || isempty(txt), txt=''; end

% no control record wanted
if exist('norec','var') && ~isempty(norec)
        foldprint(fid,var,ftfmt,norec,txt);
        return;
end

cfmt='%10d%10d%-20s%10d     %s\n';  % control record format except of %10g yields string of more than 10 chars

% variables are not constant, use CNSTNT=0 and print the actual values
if ~all(var(:)==var(1))
        IREAD=unit; CONSTNT=0; IPRN=3;
        fprintf(fid,cfmt,IREAD,CONSTNT,ftfmt,IPRN,txt);
        foldprint(fid,var,ftfmt);
        return;
end


%% remainder: control record wanted, and all variables are the same

IREAD=0; CNSTNT=var(1); IPRN=-1;

if mlfmt(end)=='d',
   fprintf(fid,cfmt,IREAD,CNSTNT,ftfmt,IPRN,txt);  % no float use as was
else
   cfmt='%10d%10s%-20s%10d     %s\n';   % now uses %10s instead of %10g
   s=sprintf('%g',CNSTNT);  % see what would be printed
   if length(s)>10          % if s requires >10 chars
      i=findstr('e',s);     % find end of field before exponent
      if ~isempty(i)
          s=[s(1:(10-(length(s)-i+1))) s(i:end)];  % truncate and assemble to field of 10 char long
      else
          s=s(1:10);
      end
   end
    fprintf(fid,cfmt,IREAD,s,ftfmt,IPRN,txt);       % print iwth float as string of 10 char wide
end

function foldprint(fid,var,ftfmt,norec,txt)
% foldprint(fid,var,ftfmt [,norec [, txt])
% prints according to fortran format ftfmt
% if norec is given (arbitrary) var is printed as a list else it is printed
% as an array
% if text is present, the text is printed at the end of the list.
% row by row with a newline behind every row.

[mlfmt NW]=ftfmt2ml(ftfmt);

if exist('norec','var')
    % print the variable as list according to the format with a text at the end
    for i=1:numel(var)
        fprintf(fid,mlfmt,var(i));
        if rem(i,NW)==0 && i~=numel(var)
            fprintf(fid,'\n');
        end
    end
    fprintf(fid,'     %s\n',txt);

else
    % print layers, but row by row in chunks of NW
    [Nrow Ncol]=size(var);
    
    for i=1:Nrow
        k1=1; k2=min(k1+NW-1,Ncol);
        while k1<=Ncol
            fprintf(fid,mlfmt,var(i,k1:k2));
            if k2<Ncol
                fprintf(fid,'\n');
                k1=k1+NW;
                k2=min(k2+NW,Ncol);
            else
                break;
            end
        end
        fprintf(fid,'\n');
    end
end
