function [tile_grid, N_land_ind, land_ind_i, land_ind_j, compr_grid_data ] ...
    = read_grid_data(fname, N_out_fields )
%  
% reichle, 20 Jul 2010
% reichle, 15 Mar 2011 - adapted to read flexible number of N_out_fields
%
% ------------------------------------------------------------------

if ~exist('N_out_fields','var')
  N_out_fields = 44;
end

int_precision   = 'int32';      % precision of fortran tag
float_precision = 'float32';    % precision of data in input file

disp(['reading from ', fname])

ifp = fopen( fname, 'r', 'b' );

% read first record (grid definitions)

fortran_tag   = fread( ifp,  1,    int_precision   ) ;
tmp_data_str  = fread( ifp, 40,    'uint8=>char'   )';
tmp_data_int  = fread( ifp, [1 7], int_precision   ) ;
tmp_data_real = fread( ifp, [1 6], float_precision ) ;
fortran_tag   = fread( ifp, 1,     int_precision   ) ;

tile_grid.gridtype = tmp_data_str; 
                     
tile_grid.ind_base = tmp_data_int(1); 
                     
tile_grid.i_dir    = tmp_data_int(2);  
tile_grid.j_dir    = tmp_data_int(3); 
                     
tile_grid.N_lon    = tmp_data_int(4); 
tile_grid.N_lat    = tmp_data_int(5); 
                     
tile_grid.i_offg   = tmp_data_int(6); 
tile_grid.j_offg   = tmp_data_int(7); 
                     
tile_grid.ll_lon   = tmp_data_real(1);  
tile_grid.ll_lat   = tmp_data_real(2);  

tile_grid.ur_lon   = tmp_data_real(3);  
tile_grid.ur_lat   = tmp_data_real(4);  
                     
tile_grid.dlon     = tmp_data_real(5);  
tile_grid.dlat     = tmp_data_real(6);  

% read second record (N_land_ind)

fortran_tag   = fread( ifp, 1, int_precision );
N_land_ind    = fread( ifp, 1, int_precision );
fortran_tag   = fread( ifp, 1, int_precision );

N_grid = tile_grid.N_lon*tile_grid.N_lat;

% read third and fourth records (land_ind_i/j)

fortran_tag   = fread( ifp, 1, int_precision );
land_ind_i    = fread( ifp, [1 N_land_ind], int_precision );
fortran_tag   = fread( ifp, 1, int_precision );

fortran_tag   = fread( ifp, 1, int_precision );
land_ind_j    = fread( ifp, [1 N_land_ind], int_precision );
fortran_tag   = fread( ifp, 1, int_precision );

% read remaining records

compr_grid_data = NaN*ones(N_out_fields,N_land_ind);

for i=1:N_out_fields
  
  fortran_tag = fread( ifp, 1, int_precision );
  tmp_data    = fread( ifp, [1 N_land_ind], float_precision );
  fortran_tag = fread( ifp, 1, int_precision );
  
  compr_grid_data(i,:) = tmp_data;
  
end

fclose(ifp);
  
% ======= EOF ==========================================================