function [trajmat, outfile] = gdftraj2mat(infile, outfile)
%
% function [trajmat, outfile] = gdftraj2mat(infile, outfile)
%
%    GDFTRAJ2MAT - Convert GPT ASCII data file into Matlab format
%
%    GDFTRAJ2MAT is used to convert GPT trajectory data, from a GDF ASCII
%    file, into a format that Matlab can read.  The input file must be a
%    specific 'trajectory' file created with the GPT command 'gdftrans'
%    must be in ASCII format.
%
%    gdftraj2mat(infile, outfile) - reads in GPT file INFILE, converts the
%    data into Matlab format, and saves it in the file OUTFILE.  If the
%    file OUTFILE already exists, and contains partially written particle
%    data, GDFTRAJ2MAT carries on from where it left off.  This works
%    around Matlabs memory limits.
%
%    [trajmat] = gdftraj2mat(infile) - reads in GPT file INFILE, converts
%    the data into Matlab format, and outputs the data as the structure
%    array TRAJMAT.
%
%    [trajmat, outfile] = gdftraj2mat(infile, outfile) - as above, but
%    echoes the output file OUTFILE.

    if nargin < 1
        error('Too few input variables: syntax is [trajmat, outfile] = gdftraj2mat(infile, outfile)') ;
    elseif nargin > 2
        error('Too many input variables: syntax is [trajmat, outfile] = gdftraj2mat(infile, outfile)') ;
    elseif nargout > 2
        error('Too many output variables: syntax is [trajmat, outfile] = gdftraj2mat(infile, outfile)') ;
%    elseif ( nargin < 2 ) || ( nargout < 1 )
%        error('Incorrect syntax: must specify either output variable or outfile (see help)') ;
    elseif ( nargout >= 2 ) && ( nargin < 2 )
        error('Must specify output file to echo it back: [trajmat, outfile] = gdftraj2mat(infile, outfile)') ;
    elseif ~ischar(infile)
        error('Input filename must be a string') ;
    elseif ( nargin >= 2 ) && ~ischar(outfile)
        error('Output filename must be a string') ;
    elseif ( exist(infile,'file') ~= 2 )
        error('Input file cannot be found')
    end

    partnum = [] ;
    ID = [] ;
    colheads = [] ;
    goodline = 0 ;
    loopnum = 0 ;
    trajmat = [] ;

    fid = fopen(infile) ;
    ID = textscan(fid, '%s %n',1) ;
    colheads = textscan(fid, '%s %s %s %s %s %s %s %s %s %s %s %s %s %s %s %s %s %s %s %s',1) ;

    for i = 1:length(colheads)

        if ~isempty(char(colheads{i}))
            colhead{i} = colheads{i} ;
        else
            break
        end

    end

    numcols = length(colhead) ;
    hstr = [] ;
    nstr = [] ;

    for i = 1:numcols

        hstr = [hstr '%s '] ;
        nstr = [nstr '%n '] ;

    end

    hstr(end) = [] ;
    nstr(end) = [] ;

    fclose(fid) ;

    fid = fopen(infile) ;

    disp(['fid = ' num2str(fid)]) ;

    while goodline == 0

        [ID, filepos] = textscan(fid, '%s %n',1) ;
        idsize = size(ID{2}) ;
        if idsize(1) == 0
            goodline = 1 ;
            break
        else
            goodline = 0 ;
        end
        loopnum = loopnum + 1 ;
        partnum(loopnum) = ID{2} ;
        idnum = ID{2} ;
        [headings, filepos] = textscan(fid, hstr, 1) ;
        [dat, filepos] = textscan(fid, nstr) ;
        datsize = size(dat) ;

        if ( rem(loopnum,1000) == 0 )
            disp(['Reading particle ' num2str(loopnum) ' ... ']) ;
        end

%        disp(['ID ' num2str(idnum)]) ;

        for j = 1:datsize(2)

            eval(['trajmat(' num2str(idnum) ').' char(colhead{j}) ' = dat{j} ;']) ;

        end

        eval(['trajmat(' num2str(idnum) ').xp = atan2(trajmat(' num2str(idnum) ...
            ').Bx,trajmat(' num2str(idnum) ').Bz) ;']) ;
        eval(['trajmat(' num2str(idnum) ').yp = atan2(trajmat(' num2str(idnum) ...
            ').By,trajmat(' num2str(idnum) ').Bz) ;']) ;

    end

    fclose(fid) ;

    if nargin > 1

        save(outfile,'trajmat') ;

    end

    return