% Usage:
%
%  terrain = LoadTGTerrain('filename.ter');
%
% Returns terrain as NxN matrix, where N is a typical
% terrain size (129, 257, 513, etc).
%
% Note: Terrain will appear right-side up when using,
% e.g., imshow(), but for PROCESSING purposes, you'll
% probably want to flipud() the resulting matrix to
% get it into Terragen's "upside-down" coordinate system.
% And, of course, don't forget TG's coordinates go from
% [0,0]->[N-1,N-1] instead of Matlab's [1,1]->[N,N].
%
% Sean O'Malley
% http://www.cs.uh.edu/~somalley
% Updated: June 14, 2004

function [ter] = LoadTGTerrain(filename)
    % Open file
    tgfile = fopen(filename, 'rb');
    
    % Dummy result in case of error
    ter = 0;
    
    % Find identifier
    id = fread(tgfile, 16, '*char')';
    
    if (~strcmp(id, 'TERRAGENTERRAIN '))
        fprintf('Error: Not a Terragen terrain file.\n');
        return
    end
    
    % Terrain dimensions
    Width = -1;
    Height = -1;
    
    % Go
    done = 0;
    
    while (~done)
        id = fread(tgfile, 4, '*char')';
        
        switch id
            case 'SIZE'
                val = fread(tgfile, 1, 'short') + 1;
                Width = val;
                Height = val;
                val = fread(tgfile, 1, 'short');
            case 'XPTS'
                Width = fread(tgfile, 1, 'short');
                val = fread(tgfile, 1, 'short');
            case 'YPTS'
                Height = fread(tgfile, 1, 'short');
                val = fread(tgfile, 1, 'short');
            case 'ALTW'
                done = 1;
        end
        
        if (feof(tgfile))
            done = 2;
        end
    end
    
    if (done == 2)
        fprintf('Error: Premature EOF.\n');
        fclose(tgfile);
        return
    end
    
    if (Width == -1 || Height == -1)
        fprintf('Error: One or more terrain dimensions undefined.\n');
        fclose(tgfile);
        return
    end

    % Read last 2 variables
    heightScale = fread(tgfile, 1, 'short');
    baseHeight = fread(tgfile, 1, 'short');
    
    % Read terrain data
    ter = fread(tgfile, [Height Width], 'short')';
    ter = flipud(ter);
    
    % Done with file
    fclose(tgfile);
    
    % Scale/translate
    ter = baseHeight + ((heightScale / 65536) * ter);
return